The evaluation of the integrals to derive the PKN governing equation
The average velocity is:
(1)\[\bar u_x = \frac{\partial p}{\partial x} \frac{2ab}{\mu \pi (a^2+b^2)}
[\int_0^{a} \int_0^{b\sqrt{1-(\frac{z}{a})^2}}
(\frac{z^2}{a^2}+\frac{y^2}{b^2})dydz
- \frac{\pi}{4}ab]\]
And we are looking for the integral:
(2)\[I = \int_0^{a} \int_0^{b\sqrt{1-(\frac{z}{a})^2}}
(\frac{z^2}{a^2}+\frac{y^2}{b^2})dydz\]
The procedure is as follows:
(3)\[\begin{split}\begin{equation}
\begin{split}
I& = \int_0^a \int_0^{b\sqrt{1-(\frac{z}{a})^2}}
(\frac{z^2}{a^2}+\frac{y^2}{b^2})dydz \\
& = \int_0^a ( \frac{z^2}{a^2}y + \frac{y^3}{3b^2} )
\vert_0^{ b\sqrt{ 1-(\frac{z}{a})^2} }dz \\
& = \int_0^a \lbrace \frac{z^2}{a^2} b \sqrt{1-(\frac{z}{a})^2}
+\frac{b}{3}[1-(\frac{z}{a})^2]^{\frac{3}{2}}\rbrace
dz \\
& = \frac{b}{a^3}\int_0^a z^2 \sqrt{a^2-z^2}dz +
\frac{b}{3a^3}\int_0^a(a^2-z^2)^{\frac{3}{2}}dz \\
& = \frac{b}{a^3}I_1 + \frac{b}{3a^3}I_2 \\
& = \frac{b}{a^3}(I_1 + \frac{I_2}{3})
\end{split}
\end{equation}\end{split}\]
And:
(4)\[\begin{split}\begin{equation}
\begin{split}
I_1& = \int_0^a z^2 \sqrt{a^2-z^2}dz \\
& = \frac{z}{8}(2z^2-a^2)\sqrt{a^2-z^2}
+\frac{a^4}{8} \textrm{arcsin} \frac{z}{a} \vert_0^a \\
& = \frac{a^4}{8} \textrm{arcsin}(1) \\
& = \frac{a^4}{8} \times \frac{\pi}{2} \\
& = \frac{\pi a^4}{16}
\end{split}
\end{equation}\end{split}\]
(5)\[\begin{split}\begin{equation}
\begin{split}
I_1& = \int_0^a(a^2-z^2)^{\frac{3}{2}}dz \\
& = \frac{z}{8}(5a^2-2z^2)\sqrt{a^2-z^2}
+\frac{3a^4}{8} \textrm{arcsin} \frac{z}{a} \vert_0^a \\
& = \frac{3a^4}{8} \textrm{arcsin}(1) \\
& = \frac{3a^4}{8} \times \frac{\pi}{2} \\
& = \frac{3\pi a^4}{16}
\end{split}
\end{equation}\end{split}\]
Substitute Eqn (4) and Eqn (5) into Eqn (3):
(6)\[\begin{split}\begin{equation}
\begin{split}
I& = \frac{b}{a^3}(I_1 + \frac{I_2}{3}) \\
& = \frac{b}{a^3}( \frac{\pi a^4}{16} +
\frac{1}{3} \times \frac{3\pi a^4}{16} ) \\
& = \frac{\pi ab}{8}
\end{split}
\end{equation}\end{split}\]
Substitute Eqn (6) into Eqn (1):
(7)\[\begin{split}\begin{equation}
\begin{split}
\bar u_x& = \frac{\partial p}{\partial x} \frac{2ab}{\mu \pi (a^2+b^2)}
(I- \frac{\pi}{4}ab) \\
& = -\frac{\partial p}{\partial x} \frac{a^2b^2}{4\mu (a^2+b^2)}
\end{split}
\end{equation}\end{split}\]
Obtaining \(\bar{T}\)
The governing equation with the boundary condition and initial condition is:
(8)\[\begin{split}\left\{
\begin{align}
& \frac{d^2 \bar{T}}{dr^2} + \frac{1}{r}\frac{d \bar{T}}{dr} - q^2 \bar{T}
= -\frac{f(r)}{\kappa}, \ \ \ \ r>r_w \\
& \bar{T}=\frac{T_0}{s}, \ \ \ \ r=r_w \\
& \bar{T}=0, \ \ \ \ \ \ \ \ t=0 \\
\end{align}
\right.\end{split}\]
Because only \(K_0(qr)\) is involved in the solution,
the solution should be in the form of \(c_1K_0(qr)\) with \(c_1\) an
arbitrary constant.
With the boundary condition, \(K_0(qr)\) is reduced to a constant
\(K_0(qr_w)\) . Thus, a term that equals to \(K_0(qr_w)\) must be
included in the denominator of \(c_1\) to cancel itself so that the value
\(\frac{T_0}{s}\) can be obtained. In this way,
\(c_1\) may be written as \(c_1 = \frac{T_0}{sK_0(qr_w)}\) .
Therefore, the solution is:
(9)\[\bar{T}=\frac{T_0K_0(qr)}{sK_0(qr_w)}\]
Similar procedures can be used for the other two boundary conditions.
Asymptotic expansion of the temperature solution
If \(\bar{f}(\lambda)\) has a branch point at the origin,
a series of equation is available:
(10)\[f(t) = \frac{1}{2\pi i} \int_{\gamma-i\infty}^{\gamma+i\infty} e^{\lambda t}
\bar{f}(\lambda)d\lambda
= \frac{1}{2\pi i} \int_{\infty}^0 e^{\lambda t}
\bar{f}(\lambda)d\lambda\]
(11)\[\frac{1}{2\pi i} \int_{\infty}^0 e^{\lambda t} \lambda^{-\alpha-1} d\lambda
=\frac{t^\alpha}{\Gamma(1+\alpha)}, \ \ \ \ \alpha \in \mathbb{R}\]
(12)\[\frac{1}{2\pi i} \int_{\infty}^0 e^{\lambda t} \lambda^n d\lambda
=0, \ \ \ \ \ n \in \mathbb{N^0}\]
(13)\[\frac{1}{2\pi i} \int_{\infty}^0 e^{\lambda t} \frac{\ln (k\lambda)}
{\lambda} d\lambda = -\ln{ \frac{ct}{m} }\]
(14)\[\frac{1}{2\pi i} \int_{\infty}^0 e^{\lambda t} \ln (k\lambda) d\lambda
= -\frac{1}{t}\]
(15)\[\frac{1}{2\pi i} \int_{\infty}^0 e^{\lambda t} \ln^2 (k\lambda) d\lambda
= -\frac{2}{t} \ln{ \frac{ct}{m} }\]
In addition, the modified nth order Bessel function of second kind
can be expanded as:
(16)\[K_0(z) = -\ln(\frac{zc}{2})I_0(z) + \sum_{p=1}^\infty \sum_{r=1}^p r^{-1}
\frac{ (\frac{z}{2})^{2p} }{p!}\]
and
(17)\[\begin{split}\begin{equation}
\begin{split}
K_n(z)&= (-1)^{n+1} \ln (\frac{zc}{2}) + \frac{1}{2}(-1)^n
\sum_{r=0}^\infty \frac{ (\frac{z}{2})^{n+2r} }{r!(n+r)!}
[ \sum_{p=1}^{n+r}m^{-1} + \sum_{p=1}^r m^{-1}] \\
& \ \ \ \ + \frac{1}{2} \sum_{r=0}^{n-1}(-1)^r( \frac{z}{2} )^{-n+2r}
\frac{ (n-r-1)! }{r!} , \ \ \ \ n > 0 \\
\end{split}
\end{equation}\end{split}\]
where:
- \(\ln c = 0.5722\) is known as Euler’s constant;
- \(m\) is real and positive;
- \(p \in \mathbb{N}\) .
And the temperature solution is:
(18)\[\bar{T} = \frac{q_0''K_0(qr)}{k s q K_1(qr_w)}\]
Apply Eqn (10) on Eqn (18) with writing \(\mu = \sqrt{
\frac{\lambda}{\kappa} }\):
(19)\[T=\frac{q_0''}{2\pi k i} \int_{-\infty}^{0+} \frac{e^{\lambda t}K_0(qr)}
{\lambda \mu K_1(qr_w)}d\lambda\]
Then the application of Eqn (11) through Eqn (17) on Eqn
(19) yields:
(20)\[\begin{split}\begin{equation}
\begin{split}
T&= \frac{q_0''r_w}{2k} \lbrace \ln\frac{4\kappa t}{cr^2} +
\frac{r_w^2}{\kappa t}\ln \frac{4\kappa t}{cr^2} +
\frac{1}{4\kappa t}[r_w^2-r^2-2r_w^2 \ln \frac{r_w}{r}+..
.]\rbrace \\
&= \frac{q_0''r_w}{2k}\ln\frac{4\kappa t}{cr^2} + O(\frac{r_w^2}{\kappa t})
\end{split}
\end{equation}\end{split}\]
Eqn (20) is suitable for large \(t_D\) defined as
\(t_D = \frac{\kappa t}{r_w^2}\) as the higher order terms vanish with
large \(t_D\) . When \(t_D\) becomes smaller, say
\(t_D \le 0.02\) , the asymptotic expansions of the Bessel
functions are needed. User may refer Chapter 13.3 of Ref.1 for details.