As we discussed in lab, we can exploit the linearity of the
Schrödinger equation to avoid overflow of Psi and Psip in schint()
by renormalizing the solution by multiplying it through by a small
constant each time that the magnitude of the solution exceeds a
critical value. The IEEE range of double precision values is somewhat
beyond 
to 
. To ensure enough room to avoid
overflows, it is safest to keep the magnitude of the solution between
and 
. You may do this by scaling the
solution by every time its absolute value exceeds
. Once again, it is best to collect such special numerical
constants in '#define' statements near to top of your file. You
should thus add the following lines to the top of ``physics.c''
next to the '#define TOL (1.e-12)' statement,
Then, you should modify your main loop in schint(), so that whenever the absolute value of Psi (probably called y[1]) exceeds ``SQRTBIG'', your entire solution is multiplied by ``SMALL''.#define SMALL (1.e-300) #define SQRTSMALL (1.e-150) #define BIG (1.e300) #define SQRTBIG (1.e150)
in your differential equation (i.e., both y[1]
and y[2]), and you must also rescale all of the values
previously generated and stored in Psiout[]. The rescaling of
Psiout[] will require writing another loop and some care to be certain
that all of the appropriate elements of Psiout[] are rescaled, for
both dk>0 and dk<0.