You raise some interesting issues.
With regard to the Boltzmann-Matano technique, the equation I gave in class
is directly applicable to the case of pre-dep diffusion, but it must be
modified to apply to drive-in diffusion. Take a look at the article at
http://doc.tms.org/ezMerchant/prodtms.nsf/ProductLookupItemID/MMTA-9910-2605
/$FILE/MMTA-9910-2605F.pdf?OpenElement
Actually, only the first page of this article is displayed at that URL, but
the essence of the method is described. It says that the point x=0 must
occur at the "Matano interface", which is the plane at which the amount of
impurity diffusing out of the left-hand region is equal to the amount of
impurity diffusing into the right-hand region. For a pre-dep diffusion,
this plane is at the surface of the wafer. But for a drive-in, the Matano
interface is below the wafer surface (within the silicon wafer), and
although it is hard to locate exactly, it eliminates the problem you point
out at the surface of a gaussian distribution, where the derivative of N(x)
goes to zero.
With regard to the digital integration of the diffusion equation with
non-constant D, there are certainly other complications to consider. One
obvious problem you did not address is what happens if the medium in which
the impurity is diffusing is made up of layers having abrupt interfaces, and
different diffusion coefficients. Then, mathematically, the derivative
dD/dx is a delta function at those interfaces. However, there are boundary
conditions at the interface. The impurities, in equilibrium, must satisfy
certain thermodynamic properties, similar to the segregation coefficient,
ko, that we saw at liquid-solid interfaces. For solid-solid interfaces, it
is usually called m (equal to C1/C2, the ratio of the concentrations on the
two sides of the interface). Also, one must guarantee that the flux on the
two sides are equal, so that there is no continuing build-up of
concentration at the interface. Thus, these conditions need to be
incorporated into the simulation of diffusion. (The problem is similar to
what happens when an electric field impinges on an interface between two
materials having a different dielectric constant. The electric flux normal
to the interface is continuous.) If the diffusivities vary slowly with x or
with N, then the extra terms that you suggest are needed would be second
order corrections. You might try to do this type of problem using your
computer model, one with the extra terms and one without, and compare the
results and see which makes sense.
Thank you for your thoughts on this problem.
Prof. Afromowitz