Kris's Research Notes

October 12, 2010

Diffusion through a Droplet

Filed under: GaAs Simulations — Kris Reyes @ 6:14 pm

Recall we wish to simulate the crystallization of Ga droplets with the presence of As flux. Let us simplify a bit, and suppose the Ga droplet is hemispherical of radius R. Let’s describe the concentration of  As in the droplet by a function u(r, \theta), where we have used spherical coordinates to describe the location inside the hemisphere (so that our domain is r \in [0, R], \theta \in [0, \pi]).

We wish to solve Laplace’s equation (WHY??):

\displaystyle \nabla^2 u =\frac{1}{r}\frac{\partial}{\partial r}\left( r \frac{\partial u}{\partial r}\right) + \frac{1}{r^2}\frac{\partial^2 u}{\partial \theta ^2} = 0

subject to boundary conditions

  • \displaystyle \frac{\partial u}{\partial r}\big|_{r = R} = c, i.e. constant flux on the surface of the droplet;
  • u(r, 0) = u(r,\pi)= 0, i.e. no concentration on the bottom of the droplet.

To solve this, we assume a separable solution of the form:

u(r, \theta) = F(r)G(\theta).

Applying the Laplacian we get

\displaystyle \nabla^2\left[F(r)G(\theta)\right] = \frac{1}{r} \frac{dF}{dr} G + \frac{d^2F}{dr^2}G + \frac{1}{r^2}F\frac{d^2G}{d\theta^2} = 0.

Separating variables, we get:

\displaystyle \frac{r^2}{F}\left(\frac{1}{r}\frac{dF}{dR} + \frac{d^2F}{dr^2}\right) = \lambda = -\frac{1}{G} \frac{d^2G}{d\theta^2},

for some constant \lambda. This gives us our ODEs:

\displaystyle \lambda F - r\frac{dF}{dr} - r^2\frac{d^2F}{dr^2} = 0.

\displaystyle \lambda G + \frac{d^2G}{d\theta^2} = 0.

We consider the second equation first, which is an eigenvalue problem.

If \lambda > 0, we have

G(\theta) = c_1 \cos(\sqrt{\lambda} \theta) + c_2 \sin(\sqrt{\lambda}\theta)

and given the boundary conditions G(0) = G(\pi) = 0, we see that \lambda = n^2 for n = \mathbb Z^+ and the corresponding eigenvector is of the form:

G_n(\theta) = c_n \sin(n \theta).

Now if \lambda < 0, we have a general solution of the form

G(\theta) = b_1\cosh(\sqrt{|\lambda|}\theta) + b_2\sinh(\sqrt{|\lambda|}\theta)

However, only the trivial solution satisfies the boundary conditions.

Lastly if \lambda = 0, we have G(\theta) = a_1 + a_2\theta, but once again only the trivial solution satisfies boundary conditions.

Therefore the only eigenvalues of the above operator are of the form \lambda = n^2 for n \in \mathbb Z^+ with corresponding eigenvectors G_n(\theta) = c_n \sin(n\theta).

Then we have corresponding solutions to the first ODE: F_n(r) = a_n r^p_n, where p satisfies

n^2 r^{p_n} - p_ nr^{p_n} - p_n(p_n-1)r^{p_n} = 0,

which is what is obtained if we plug in this form of F into the above ODE. That is, p_n = \pm n and F(r) = a_nr^{\pm n}. We require F(r) is bounded near r = 0, and so we ignore solutions F(r) = a_nr^{-n}.

Therefore we have a general solution of the form (incorporating all undetermined coefficients)

\displaystyle u(r,\theta) = \sum_{n=1}^\infty A_n r^{n}\sin(n\theta).

We now apply the final boundary condition \frac{\partial u}{\partial r}\big|_{r = R} = c:

\displaystyle \frac{\partial u}{\partial r}\big|_{r = R} = \sum_{n=1}^\infty nA_n R^{n-1} \sin(n\theta) = c.

By Fourier Inversion, we can solve

\displaystyle n A_n R^{n-1} = \frac{c}{\pi} \int_0^\pi \sin(n\theta) d\theta = \frac{c}{n\pi} \left(1-\cos(n\pi)\right)

and so

\displaystyle A_n = \frac{c}{\pi}\frac{(1-\cos(n\pi))}{n^2R^{n-1}}.

Plugging this in to the expression for u(r, \theta) we get:

\displaystyle u(r, \theta) = \frac{c}{\pi} \sum_{n=1}^\infty \frac{1-\cos(n\pi)}{n^2R^{n-1}} r^{n}\sin(n\theta).

We may compute the flux through the bottom of the crystal by calculating:

\displaystyle \frac{\partial u}{\partial y} = \frac{\partial u}{\partial r}\frac{\partial r}{\partial y} + \frac{\partial u}{\partial \theta}\frac{\partial \theta}{\partial y}.

That is:

\displaystyle \frac{\partial u}{\partial y} = \frac{c}{\pi} \sum_{n=1}^\infty \frac{1-\cos(n\pi)}{n R^{n-1}} \left[\sin\theta \sin(n\theta) r^{n-1} + \frac{\cos(n\theta)}{\cos\theta} r^{n-1} - \frac{\cos(n\theta)\sin^2\theta}{\cos\theta}r^{n-1}\right]

and when so \theta = 0, i.e. when y=0, x \geq 0, we get:

\displaystyle \frac{\partial u}{\partial y}\bigg|_{\theta = 0} = \frac{c}{\pi} \sum_{n=1}^\infty \frac{1-\cos(n\pi)}{n R^{n-1}} r^{n-1} = \frac{2c}{\pi}\sum_{k=1}^\infty \frac{1}{2k+1}\left(\frac{r}{R}\right)^{2k}.

I don’t know of a closed form for this, but we may plot truncated sums (for e.g. n \leq 10000) and examine the flux vs r. We may also vary the flux c. Here are the plots of \frac{\partial u}{\partial y}\big|_{\theta = 0} for R = 1, 0 \leq r \leq R and c varying from 0.1 to 1.0.



Note this is what we want: the amount of material being deposited is greater near the edge of the droplet than the middle.


1 Comment »

  1. […] we do not get the distribution of crystal inside the droplet that we expected in this post. That is, we had previously calculated that there should be more atoms diffusing to the edges of […]

    Pingback by Droplet Crystallization – Part 3 « Kris's Research Notes — October 19, 2010 @ 11:32 pm

RSS feed for comments on this post. TrackBack URI

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: