Fick’s Laws of Diffusion


I’ve realized that Fick’s Laws for diffusion are very important for a lot of different fields, but aren’t very easily searchable.

Below is a short compilation of my *favorite* problems and their solutions. Just for fun!


Steady State Diffusion through the Walls of a Cylindrical Pipe (From Jeff Hoyt’s MSE 103 Course at UC Berkeley)

Consider hydrogen diffusing through an infinitely long cylindrical pipe of Pd. Let the inside radius of the pipe be r = Ri and the outer surface be r = Ra. Assume the concentration of H is maintained at c = ci on the inside surface and c = ca on the outside surface of the cylinder. Determine the steady state concentration profile c(r).


First we can write out Fick’s second law for a cylinder, knowing that we are in steady state for this problem.

\frac{\partial c}{\partial t} = \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial c}{\partial r} \right )

\frac{\partial c}{\partial t} = 0

0 = \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial c}{\partial r} \right )

Then we note that the solution to this equation is logarithmic.

\int \int \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial c}{\partial r} \right ) = A \ln r + B

Finally, we need to apply the described boundary conditions. These are the outer and inner concentrations of the cylinder.

c(R_i) = c_i

c(R_a) = c_a

We solve for A and B.

A \ln(R_i) + B = c_i

A \ln(R_a) + B = c_a

A \ln(R_a/R_i) = (c_a - c_i)

A = \frac{(c_a - c_i)}{\ln(R_a/R_i)}

B = c_a - \frac{(c_a - c_i)}{\ln(R_a/R_i)} \ln(R_a)

Finally, we can plug in the constants and get the final concentration for steady state.

c(r) = c_a - \frac{(c_a - c_i)}{\ln(R_a/R_i)} \ln(R_a/r)


Steady State Diffusion through Stacked Cylindrical Pipes, a Mass Balance (From Terry Johnson’s BioE 104 Course at UC Berkeley)


Screen Shot 2015-02-20 at 1.19.27 PM

You are interested in steady-state mass transfer (without reaction) through the walls of a cylinder that consists of two membranes (as shown). This is exactly like Example 6.5, only in cylindrical coordinates. For simplicity’s sake, assume that the partition coefficients for each membrane are ~ 1.0.

I’ve also neglected to include the substrate ‘i’s to denote solute i – there is only one solute in the system.

Membrane 1 has a diffusivity D1 from r = R0 to R1 and membrane 2 has a diffusivity D2 from r = R1 to R2.

Given that the concentration profile through a single membrane is C(r) = A•ln(r) + B (Example 6.7), solve for the flux through the membranes in terms of a sum of resistors in series. Your final equation should be in a form that makes it easy to recognize the individual membrane resistances.



If we were dealing with a normal profile through a cylinder, we know that the equation describing the concentration profile is logarithmic.

C(r) = A \ln(r) + B

First we can note that the mass transfer of solute is equal at the membrane interfaces otherwise there would be a buildup of solute.

J(R_1) = J_1(R_1) = J_2(R_1)

Now we can write out the concentration profile for each membrane independently.

C_1(r) = A_1 \ln(r) + B_1

C_2(r) = A_2 \ln(r) + B_2

Then we list out boundary conditions (there are four) and the equations we get from them.

C_1(R_0) = C_0

C_0 = A_1 \ln(R_0) + B_1

C_2(R_2) = C_L

C_L = A_2 \ln(R_2) + B_2

C_1(R_1) = C_2

A_1 \ln(R_1) + B_1 = A_2 \ln(R_2) + B_2

J_1(R_1) = J_2

D_1 A_1 = D_2A_2

Solve for constants, A1, A2, B1, and B2. (this gets a little messy).

A_1 = \frac{(C_0-C_L)}{\ln\left(\frac{R_0}{R_1}\right)-\frac{R_0}{R_1} \ln\left(\frac{R_2}{R_1}\right)}

A_2 = \frac{(C_0-C_L)}{\frac{D_2}{D_1}\ln\left(\frac{R_0}{R_1}\right)-\ln\left(\frac{R_2}{R_1}\right)}

B_1 = C_0 - \left(\frac{(C_0-C_L)}{\ln\left(\frac{R_0}{R_1}\right)-\frac{D_2}{D_1}\ln\left(\frac{R_0}{R_1}\right)}\right)\ln(R_0)

B_2 = C_L - \left(\frac{(C_0-C_L)}{\frac{D_2}{D_1}\ln\left(\frac{R_2}{R_1}\right)-\ln\left(\frac{R_0}{R_1}\right)}\right)\ln(R_2)

C_1 = C_0 - \left(\frac{(C_0-C_L)}{\ln\left(\frac{R_0}{R_1}\right)-\frac{D_1}{D_2}\ln\left(\frac{R_0}{R_1}\right)}\right)\ln\left(\frac{r}{R_0}\right)

C_2 = C_L - \left(\frac{(C_0-C_L)}{\frac{D_2}{D_1}\ln\left(\frac{R_0}{R_1}\right)-\ln\left(\frac{R_2}{R_1}\right)}\right)\ln\left(\frac{r}{R_2}\right)

Use that to solve for Flux (also messy).

J = -D_j \frac{dC_j}{dr}

N_1 = -\frac{D_1}{r} \left(\frac{(C_0-C_L)}{\ln\left(\frac{R_0}{R_1}\right)-\frac{D_1}{D_2}\ln\left(\frac{R_2}{R_1}\right)}\right)

N_2 = -\frac{D_2}{r} \left(\frac{(C_0-C_L)}{\frac{D_2}{D_1}\ln\left(\frac{R_0}{R_1}\right)-\ln\left(\frac{R_2}{R_1}\right)}\right)

Finally we can write out the flux profile.

N_1 = \left( \frac{(C_0-C_L)}{r} \frac{D_1D_2}{D_2\ln\left(\frac{R_0}{R_1}\right)-D_1\ln\left(\frac{R_2}{R_1}\right)}\right)


In this case, the driving force is the change in concentration over distance (first part), and the resistance is the second part with all the diffusivity in it. Looks a lot like Ohms law…


Concentration Dependent Steady State Diffusion through a Plate (From Jeff Hoyt’s MSE 103 Course at UC Berkeley)

A thin plate of a substance A with a thickness L is exposed on one side to a gaseous environment of species B such that the concentration of B at the surface is held constant at c = ci . The other side of the plate is held constant at a concentration of zero. For B diffusing in A the diffusion coefficient is known to vary with composition according to D(c) = D1exp (αc) where D1 and α are constants. Determine the steady state concentration profile in the plate.



This question requires a re-derivation of Fick’s second law.

\frac{\partial c}{\partial t} = \frac{\partial}{\partial x} \left(D(x) \frac{\partial c}{\partial x} \right )

\frac{\partial c}{\partial t} = \frac{\partial}{\partial x} \left(D_1\textup{e}^{\alpha c}\frac{\partial c}{\partial x} \right ) = D_1\frac{\partial}{\partial x} \left(\textup{e}^{\alpha c}\frac{\partial c}{\partial x} \right )

We are still at steady state, so concentration no longer depends on time.

0 = \frac{\partial}{\partial x} \left(D_1\textup{e}^{\alpha c}\frac{\partial c}{\partial x} \right )

Now we integrate once and rearrange.

\textup{e}^{\alpha c}\frac{\partial c}{\partial x} = A

\int \textup{e}^{\alpha c}\partial c = \int A\partial x

\frac{\textup{e}^{\alpha c}}{\alpha} = Ax + B

\ln(\textup{e}^{\alpha c}) = \ln(Ax + B)

\alpha c = \ln(Ax + B)

Now we solve for c with the new system. The form of the solution will be somewhat logarithmic.

c =\frac{\ln(Ax + B)}{\alpha}

Boundary conditions:

c(0) = 0

c(L) = c_i

Then solve for constants.

c(L)=0 = \frac{\ln (AL+B)}{\alpha}

B = 1- AL

c_i = \frac{1}{\alpha}\ln(1-AL)

A = \frac{1}{L}\left[1-(\textup{e}^{\alpha c_i}) \right ]

B = 1- \left[1-(\textup{e}^{\alpha c_i} \right ] =(\textup{e}^{\alpha c_i})

Finally plug in those constants for the concentration profile.

c = \frac{1}{\alpha}\ln \left[ \frac{1-(\textup{e}^{\alpha c_i})}{L} x+(\textup{e}^{\alpha c_i}) \right]


Superposition to Solve Time-Dependent Step Concentration Diffusion in 2 Dimensions (From Jeff Hoyt’s MSE 103 Course at UC Berkeley)

Consider one dimensional diffusion in an infinitely long rod. At t = 0 a region given by –L < x < +L has a solute concentration of c = c* and all other points in the rod are solute free. Solve for c(x,t). Notice, the dimension L is sufficiently large such that the thin film, planar source solution can NOT be used. Finally express all integrals in the solution in terms of known functions.



We’re going to use the superposition rule to solve this problem.

We can describe the spreading at any time of solute at any individual point using the initial solute concentration. This is using the planar source solution (known).

c(x, t=0) = c^*\delta(x)

c(x, t) = \frac{c^*}{2\sqrt{\pi Dt}} \textup{e}^{-x^2/4Dt}

We can express this concept over a range of positions using a sum.

c(x,t) = \sum^N_{i=-N} \frac{c(x,t=0)\Delta x}{2\sqrt{\pi Dt}} \textup{e}^{-(x-x_i)^2/4Dt}

To make this description truly continuous, we turn it into an integral.

c(x,t) =\frac{1}{2\sqrt{\pi Dt}} \int^\infty_{-\infty} c(\xi, t=0 )\textup{e}^{-(x-\xi)^2/4Dt}d\xi

Our initial condition in this problem is a concentration step between –L and L. Therefore, a piecewise function is representative.

c(x, t=0) = \left\{\begin{matrix} c^* & -L < x < L\\ 0 & \textup{otherwise } \end{matrix}\right.

So we can break down the integral into three parts representing each part of the system.

c(x,t) =\frac{1}{2\sqrt{\pi Dt}} \int^{-L}_{-\infty} 0\textup{e}^{-(x-\xi)^2/4Dt}d\xi+\frac{1}{2\sqrt{\pi Dt}} \int^L_{-L} c^*\textup{e}^{-(x-\xi)^2/4Dt}d\xi+\frac{1}{2\sqrt{\pi Dt}} \int^\infty_{L} 0\textup{e}^{-(x-\xi)^2/4Dt}d\xi

c(x,t) =\frac{c^*}{2\sqrt{\pi Dt}} \int^L_{-L} \textup{e}^{-(x-\xi)^2/4Dt}d\xi

To make this easier to solve, we’re going to substitute a new variable in. It was also make known functions more recognizable.

y = \frac{\xi-x}{\sqrt{4Dt}}

d\xi = \sqrt{4Dt} \;dy

L \rightarrow \frac{x+L}{2\sqrt{Dt}}

-L \rightarrow \frac{x-L}{2\sqrt{Dt}}

And we substitute:

c(x,t) = \frac{c^*}{\sqrt{\pi}}\int_{\frac{x-L}{2\sqrt{Dt}}}^{\frac{x+L}{2\sqrt{Dt}}} \textup{e}^{-y^2}dy

Now we can convert the form to something a little more familiar (the error function!)

c(x,t) = \frac{c^*}{2} \left[ \textup{erf}\left(\frac{x+L}{2\sqrt{Dt}}\right)-\textup{erf}\left(\frac{x-L}{2\sqrt{Dt}}\right)\right]

For those who just need a reminder, the error function form is:

\textup{erf}(z) = \frac{2}{\sqrt{\pi}}\int_0^z c^{-t^2}dt


Here‘s an excerpt from a book with a good review of other types of problems.