Episode 08.
Radial Diffusion in Spherical Solids - Unsteady State
Summary
This episode covers radial diffusion in spherical solids under unsteady states. First, we define a new variable for rearranging the PDE and then discuss the diffusion of various cases based on the PDE.
1. Introduction
In this episode, we describe radial diffusion in spherical solids under unsteady states. Since we will be dealing with radial diffusion, we will use (Eq 2.9) for discussion. Let’s rewrite (Eq 2.9) in a more familiar form to find a solution. We start by defining a new variable, \({\small\displaystyle u}\), as the product of concentration and radius.
\( {\small\displaystyle\begin{align} \frac{\partial C}{\partial t} \quad&=\quad D\left(\frac{\partial^2 C}{\partial r^2}\,+\,\frac{2}{r}\frac{\partial C}{\partial r}\right) \tag{Eq 2.9} \\[15pt] u \quad&=\quad C\cdot r \end{align} } \)
(Figure 8.1)
Using the newly defined \({\small\displaystyle u}\), we substitute \({\small\displaystyle u/r}\) for concentration \({\small\displaystyle C}\). We can see that \({\small\displaystyle u/r}\) has the following relationships.
\( {\small\displaystyle\begin{align} \frac{\partial(u/r)}{\partial t} \quad&=\quad \frac{1}{r}\frac{\partial u}{\partial t} \\[10pt] \frac{\partial(u/r)}{\partial r} \quad&=\quad \frac{1}{r}\frac{\partial u}{\partial r}\,-\,\frac{u}{r^2} \\[10pt] \frac{\partial^2(u/r)}{\partial r^2} \quad&=\quad \frac{1}{r}\frac{\partial^2 u}{\partial r^2}\,-\,\frac{2}{r^2}\frac{\partial u}{\partial r}\,+\,\frac{2u}{r^3} \tag*{} \end{align} } \)
To convert (Eq 2.9) to PDE in terms of \({\small\displaystyle u}\), put the above equations into (Eq 2.9). Then, we get (Eq 8.1).
\( {\small\displaystyle\begin{align} \frac{\partial u}{\partial t} \quad=\quad D\frac{\partial^2 u}{\partial r^2} \tag{Eq 8.1} \end{align} } \)
2. Various cases of radial diffusion in spherical solids
Now, let’s discuss how to calculate the concentration when the boundary condition is zero, when it is constant, and when it varies with time using (Eq 8.1).
Case 1 : when boundary conditions are 0.
\( {\small\displaystyle\begin{matrix} \dfrac{\partial u}{\partial t} &=& D\dfrac{\partial^2 u}{\partial r^2}&\text{for }0\lt r\lt a \\[10pt] u(0,t) &=& 0&\text{for all }t \\[10pt] u(a,t) &=& 0&\text{for all }t \\[10pt] u(r,0) &=& rf(r)&\text{for }0\lt r\lt a \tag*{} \end{matrix} } \)
As in previous episodes, we can use the separation of variables method to find the solution under these conditions.
\( {\small\displaystyle\begin{align} u \quad&=\quad R(r)\cdot T(t) \\[10pt] &=\quad \sum_{n=1}^{\infty}\,A_n \sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right),\quad A_n \,=\, \frac{2}{a}\int_0^a\,rf(r)\sin\left(\frac{n\pi r}{a}\right)\,dr \tag*{} \end{align} } \)
Hence, according to the definition of \({\small\displaystyle u}\), we can see that the concentration \({\small\displaystyle C}\) becomes (Eq 8.3).
\( {\small\displaystyle\begin{align} u \quad&=\quad \frac{2}{a}\sum_{n=1}^{\infty}\,\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)\int_0^a\,rf(r)\sin\left(\frac{n\pi r}{a}\right)\,dr \tag{Eq 8.2} \\[15pt] C \quad&=\quad \frac{2}{ra}\sum_{n=1}^{\infty}\,\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)\int_0^a\,rf(r)\sin\left(\frac{n\pi r}{a}\right)\,dr \tag{Eq 8.3} \end{align} } \)
Case 2: when boundary condition is constant.
\( {\small\displaystyle\begin{matrix} \dfrac{\partial u}{\partial t} &=& D\dfrac{\partial^2 u}{\partial r^2}&\text{for }0\lt r\lt a \\[10pt] u(0,t) &=& 0&\text{for all }t \\[10pt] u(a,t) &=& aC_0&\text{for all }t \\[10pt] u(r,0) &=& rf(r)&\text{for }0\lt r\lt a \tag*{} \end{matrix} } \)
Let’s assume that \({\small\displaystyle u}\) is the sum of variables \({\small\displaystyle v}\) and \({\small\displaystyle w}\). (Table 8.1) summarizes the PDEs and conditions for each variable.
| Variables | \( {\small\displaystyle u \quad=\quad v\,+\,w } \) | \( {\small\displaystyle v } \) | \( {\small\displaystyle w } \) |
|---|---|---|---|
| PDEs | \( {\small\displaystyle \frac{\partial u}{\partial t} \quad=\quad D\frac{\partial^2 u}{\partial r^2} } \) | \( {\small\displaystyle \frac{\partial v}{\partial t} \quad=\quad D\frac{\partial^2 v}{\partial r^2} } \) | \( {\small\displaystyle \frac{\partial^2 w}{\partial r^2} \quad=\quad 0 } \) |
| BCs | \( {\small\displaystyle\begin{align} u(0,t) \quad&=\quad 0 \\[10pt] u(a,t) \quad&=\quad aC_0 \end{align} } \) | \( {\small\displaystyle\begin{align} v(0,t) \quad&=\quad 0 \\[10pt] v(a,t) \quad&=\quad 0 \end{align} } \) | \( {\small\displaystyle\begin{align} w(0,t) \quad&=\quad 0 \\[10pt] w(a,t) \quad&=\quad aC_0 \end{align} } \) |
| IC | \( {\small\displaystyle u(r,0) \quad=\quad rf(r) } \) | \( {\small\displaystyle v(r,0) \quad=\quad rf(r)\,-\,w } \) | - |
Depending on the boundary conditions, the variable \({\small\displaystyle w}\) is:
\( {\small\displaystyle\begin{align} w \quad=\quad \frac{W}{a}r,\qquad W\,=\,aC_0 \tag*{} \end{align} } \)
For variable \({\small\displaystyle v}\), we can express it using (Eq 8.2) since its boundary conditions are the same as Case 1.
\( {\small\displaystyle\begin{align} v \quad&=\quad \frac{2}{a}\sum_{n=1}^{\infty}\,\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)\int_0^a\,[rf(r)\,-\,w]\sin\left(\frac{n\pi r}{a}\right)\,dr \\[10pt] &=\quad \frac{2}{a}\sum_{n=1}^{\infty}\,\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)\int_0^a\,\left[rf(r)\,-\,\frac{Wr}{a}\right]\sin\left(\frac{n\pi r}{a}\right)\,dr \tag*{} \end{align} } \)
By rearranging the above equation, we obtain (Eq 8.4).
\( {\small\displaystyle\begin{align} v \quad=\quad \frac{2}{\pi}&\sum_{n=1}^{\infty}\,\frac{(-1)^nW}{n}\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right) \\ \tag{Eq 8.4} \\ &\qquad \,+\,\frac{2}{a}\sum_{n=1}^{\infty}\,\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)\int_0^a\,rf(r)\sin\left(\frac{n\pi r}{a}\right)\,dr \end{align} } \)
Since \({\small\displaystyle u}\) is the sum of \({\small\displaystyle v}\) and \({\small\displaystyle w}\), \({\small\displaystyle u}\) becomes:
\( {\small\displaystyle\begin{align} u \quad=\quad \frac{Wr}{a}\,+\,\frac{2}{\pi}&\sum_{n=1}^{\infty}\,\frac{(-1)^nW}{n}\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right) \\ \tag{Eq 8.5} \\ & \,+\,\frac{2}{a}\sum_{n=1}^{\infty}\,\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)\int_0^a\,rf(r)\sin\left(\frac{n\pi r}{a}\right)\,dr \end{align} } \)
Case 3: when boundary condition is variable with time.
\( {\small\displaystyle\begin{matrix} \dfrac{\partial u}{\partial t} &=& D\dfrac{\partial^2 u}{\partial r^2}&\text{for }0\lt r\lt a \\[10pt] u(0,t) &=& 0&\text{for all }t \\[10pt] u(a,t) &=& a\phi(t)&\text{for all }t \\[10pt] u(r,0) &=& rf(r)&\text{for }0\lt r\lt a \tag*{} \end{matrix} } \)
Likewise, assume that \({\small\displaystyle u}\) is the sum of variables \({\small\displaystyle v}\) and \({\small\displaystyle w}\). (Table 8.2) summarizes the PDEs and conditions for each variable.
| Variables | \( {\small\displaystyle u } \) | \( {\small\displaystyle v } \) | \( {\small\displaystyle w } \) |
|---|---|---|---|
| PDEs | \( {\small\displaystyle \frac{\partial u}{\partial t} \quad=\quad D\frac{\partial^2 u}{\partial r^2} } \) | \( {\small\displaystyle \frac{\partial v}{\partial t} \quad=\quad D\frac{\partial^2 v}{\partial r^2} } \) | \( {\small\displaystyle \frac{\partial w}{\partial t} \quad=\quad D\frac{\partial^2 w}{\partial r^2} } \) |
| BCs | \( {\small\displaystyle\begin{align} u(0,t) \quad&=\quad 0 \\[10pt] u(a,t) \quad&=\quad a\phi(t) \end{align} } \) | \( {\small\displaystyle\begin{align} v(0,t) \quad&=\quad 0 \\[10pt] v(a,t) \quad&=\quad 0 \end{align} } \) | \( {\small\displaystyle\begin{align} w(0,t) \quad&=\quad 0 \\[10pt] w(a,t) \quad&=\quad a\phi(t) \end{align} } \) |
| IC | \( {\small\displaystyle u(r,0) \quad=\quad rf(r) } \) | \( {\small\displaystyle v(r,0) \quad=\quad rf(r) } \) | \( {\small\displaystyle w(r,0) \quad=\quad 0 } \) |
By Duhamel’s method, we can represent variable \({\small\displaystyle w}\) as follows.
\( {\small\displaystyle\begin{align} w \quad=\quad \int_0^t\,w_0^{\prime}(t\,-\,\lambda)a\phi(\lambda)\,d\lambda \tag*{} \end{align} } \)
\({\small\displaystyle w_0}\) is the function with the boundary condition of 1 at \({\small\displaystyle r\,=\,a}\) and the initial condition of 0. So, we can express \({\small\displaystyle w_0}\) using (Eq 8.5). And by substituting \({\small\displaystyle w_0}\) into the above equation, we can find the solution to the PDE for \({\small\displaystyle w}\).
\( {\small\displaystyle\begin{align} w_0 \quad&=\quad \frac{r}{a}\,+\,\frac{2}{\pi}\sum_{n=1}^{\infty}\,\frac{(-1)^n}{n}\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right) \\[10pt] \frac{\partial w_0}{\partial t} \quad&=\quad -\frac{2D}{a^2}\sum_{n=1}^{\infty}\,(-1)^n\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)n\pi \\[15pt] w \quad&=\quad -\frac{2D}{a}\sum_{n=1}^{\infty}\,(-1)^n\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right)n\pi\int_0^t\,\exp\left(\frac{Dn^2\pi^2\lambda}{a^2}\right)\phi(\lambda)\,d\lambda \tag{Eq 8.6} \end{align} } \)
Expressing \({\small\displaystyle u}\) using (Eq 8.2) and (Eq 8.6), we get:
\( {\small\displaystyle\begin{align} u \quad=\quad \frac{2}{a}\sum_{n=1}^{\infty}\,&\sin\left(\frac{n\pi r}{a}\right)\exp\left(-\frac{Dn^2\pi^2t}{a^2}\right) \\ \tag{Eq 8.7} \\ &\times\, \left[\int_0^a\,rf(r)\sin\left(\frac{n\pi r}{a}\right)\,dr\,-\,(-1)^nDn\pi\int_0^t\,\exp\left(\frac{Dn^2\pi^2\lambda}{a^2}\right)\phi(\lambda)\,d\lambda\right] \end{align} } \)
Nomenclature
\({\small\displaystyle BC}\) : boundary condition
\({\small\displaystyle C}\) : concentration
\({\small\displaystyle D}\) : diffusion coefficient
\({\small\displaystyle IC}\) : initial condition
\({\small\displaystyle PDE}\) : partial differential equation
\({\small\displaystyle r}\) : radius
\({\small\displaystyle t}\) : time
\({\small\displaystyle u,\,v,\,w}\) : variables
\({\small\displaystyle \lambda}\) : variable