schint()

int schint(double *Psi, double *Psip, double Psiout[],
           int k1, int k2, double V[], double E,
           double r[], double dr[], int N)

Return value:

• number of ``nodes'' (zero crossings) in solution between r[k1] and r[k2]

Output:

• Psi, Psip (passed by reference): values of $\Psi(r[k2])$ and $\Psi'(r[k2])$
• Psiout[]: Psi(r[k]) along the entire integration

Input:

• Psi, Psip (passed by reference): values of $\Psi(r[k1])$ and $\Psi'(r[k1])$
• k1, k2: indices of initial and final integration points (r[k1] $\rightarrow$ r[k2])
• V[]: potential on corresponding grid points, r[]
• E: electron energy
• r[], dr[], N: standard grid information

As we learned in Lab2, solving for the energies can be accomplished by setting the boundary conditions at one end point of the problem (e.g., $r(x=0)=0$), and investigating the resulting behavior at the other end point (e.g., $r(x=1)=\infty$). As we will see later in this problem set, it is also very useful to count the number of times that the solution crosses zero (number of ``nodes''). Accordingly, write a function of the above form in ``physics.c'' and include its prototype in ``physics.h''.

Note that to write this routine, you should be able to take the appropriate parts from your ``quicksch.c'' program and repackage them as a subroutine. There are a number of details to consider. Among these are

1. the integration range is no longer k=N to k=0, but k=k1 to k=k2;
2. there are two cases k1$>$k2 and k1$<$k2;
3. you are to store in Psiout[k] all the values for $\Psi(r[k])$ as they are generated;
4. you must set up a counter for the number of times that $\Psi(r)$ (i.e., y[1]) changes sign;
5. you must construct T[]=E-V[] on your own.

The last of these considerations is a bit tricky because T[] is not passed in as a variable, and the routine thus must allocate the space for the array it without ``knowing'' how much space is needed until the actual call is made. One solution would be to declare a fixed, large amount of space, but this is wasteful and can cause bugs if the size is exceeded. The best way to handle this kind of situation is ``dynamic memory allocation,'' where the space for T[] is allocated while the programming is running. Fortunately, nrutil.c includes functions for doing just this. To accomplish this, you should declare the variable T (technically, a pointer) in the declarations section of schint() using

double *T;

and, then, allocate the space for T[] as one of the first statements in the program using

T=dvector(0,N);

This latter statement sets up space for an array T[] which holds values T[0], T[1], ..., T[N]. The last, but most important, detail in dynamic memory allocation is to be sure to free up the space. Otherwise, each time the routine is called, it will eat up more program space, and eventually your code will slow down and crash (and possibly even crash the machine...)!

To avoid this, you should clean up the space when you are done with it, by including the following statement at the very end of schint(),

free_dvector(T,0,N);