Package your previous Poisson solver program into a routine of the
prototype below in ``physics.c''. Note that the return value is the actual
electrostatic potential, phi(r) 
Phi(r)/r, whereas the input is
Rho(r) 
r 
n(r). Be sure to use ``interpolate'' to find
the values on the odd grid points after you have done the Runge-Kutta
integration, but before dividing Phi[k] by r[k] to get phi[k].
void getphi(double phi[], double Rho[], double r[], double dr[], int N)Output:
- phi[k]: electrostatic potential
on the grid
Input:
- Rho[k]:
for all grid points r[k]
- r[], dr[], N: standard grid information
Finally, because your input is now Rho[k] 4*pi*r[k]*r[k]*rho[k]
instead of just rho[k], you must modify derivs_Poisson() (be sure to
include this in ``physics.c'' before getphi()!) so that it now
contains a code fragment of the form
.
.
dydx[1]=y[2]*dr[k];
dydx[2]=-Rho[k]/r[k]*dr[k];
.
.