Heat diffusion in semi-infinite rod with piecewise I.C.

Question

I am trying to solve a specific heat diffusion problem with a semi-infinite, isolated rod, where the initial temperature distribution is $10^\circ$ C for $0 < x < 1$, and $0$ for $x > 1$. The problem is thus of the form

$$ \frac{\partial u}{\partial t} - \frac{\partial^2u}{\partial x^2} = 0 $$ $$ \frac{\partial u}{\partial x}(0,t) = 0 $$ $$ u(x,0) = \begin{cases} 10 & \mathrm{for} \:\: 0 < x < 1 \\ 0 & \mathrm{for} \:\: x > 1 \end{cases}. $$

One-sided Laplace transform in $t$ and rearranging yields

$$ \frac{\partial^2 \bar{u}}{\partial x^2}(x,s) - s\bar{u}(x,s) = -u(x,0) = \begin{cases} -10 & \mathrm{for} \:\: 0 < x < 1 \\ 0 & \mathrm{for} \:\: x > 1 \end{cases}.$$

I'm not sure how to approach this, but I've tried solving for the first region,

$$ \frac{\partial^2 \bar{u}}{\partial x^2} - s\bar{u} = -10 \Rightarrow $$

$$ \bar{u} = \bar{u}_h + \bar{u}_p = a(s)e^{\sqrt{s}x} + b(s)e^{-\sqrt{s}x} + \frac{10}{s} $$

and then applying the boundary conditions,

$$ u \:\: \mathrm{bounded} \Rightarrow a(s) = 0$$ $$ \frac{\partial u}{\partial x}(0,t) = 0 \Rightarrow \frac{\partial \bar{u}}{\partial x}(0,s) = 0 \Rightarrow -\sqrt{s}b(s) = 0 \Rightarrow b(s) = 0. $$

But this means that

$$ \bar{u}(x,s) = \frac{10}{s} \Rightarrow u(x,t) = 10H(t) \:\:\:\: (0 < x < 1) $$

which does not seem to make much sense since this means that the temperature remains constant in the first region and one would expect it to decrease as $t \rightarrow \infty$. Maybe some kind soul could point out what I'm doing wrong and how to best approach this problem?

(A similar problem with a piecewise I.C. was posted here, but the solution is attempted using the Fourier transform, which doesn't seem to work.)

Answer 1

After some further research, I think I have found a solution. One can extend the problem to all $x$ by reflecting $u$ evenly into negative $x$. The even extension of $u$ is then

$$ u^+(x,t) = \begin{cases} u(x,t) & \mathrm{for} \:\: x>0 \\ u(-x,t) & \mathrm{for} \:\: x < 0 \end{cases}.$$

It can then be shown that

$$ \bigg(\frac{\partial^2u}{\partial x^2}\bigg)^+ = \frac{\partial^2u^+}{\partial x^2} - 2\frac{\partial^2u}{\partial x^2}(0,t)\delta(x). $$

Using this gives a new problem

$$ \frac{\partial u^+}{\partial t} - \frac{\partial^2u^+}{\partial x^2} = 0$$ $$ u^+(x,0) = 10H(x+1) - 10H(x-1). $$

Fourier transform in $x$ yields

$$ \frac{\partial \widehat{u^+}}{\partial t}(\omega,t) + \omega^2\widehat{u^+}(\omega,t) = 0 \Rightarrow $$ $$ \widehat{u^+} = C(\omega)e^{-\omega^2t}. $$

Applying the transformed initial condition yields

$$ \widehat{u^+}(\omega,0) = 10\frac{2\sin{\omega}}{\omega} = C(\omega) $$ $$ \therefore \:\: \widehat{u^+}(\omega,t) = 10\frac{2\sin{\omega}}{\omega} e^{-\omega^2t}.$$

Inverse transform using the convolution theorem (in $x$) now gives the extended solution

$$ u^+(x,t) = 10\big[H(x+1) - H(x-1)\big] \:\: \ast \:\: \frac{1}{\sqrt{4\pi t}}e^{-x^2/4t} = $$ $$\frac{10}{\sqrt{4\pi t}} \int_{-\infty}^\infty \big[H(\tau+1) - H(\tau-1)\big]e^{-(x-\tau)^2/4t} \:d\tau = $$ $$ \frac{10}{\sqrt{4\pi t}} \int_{-1}^1 e^{-(x-\tau)^2/4t} \:d\tau = \Bigg[\mathrm{Sub.} \:\: y=\frac{x-\tau}{\sqrt{4t}}\Bigg] = $$ $$ -\frac{10}{\sqrt{\pi}} \int_{(x+1)/\sqrt{4t}}^{(x-1)/\sqrt{4t}} e^{-y^2} \:dy = -5 \,\mathrm{erf}(y)\Bigg|_{(x+1)/\sqrt{4t}}^{(x-1)/\sqrt{4t}} =$$

$$ 5\Bigg[\mathrm{erf}\bigg(\frac{x+1}{\sqrt{4t}}\bigg) - \mathrm{erf}\bigg(\frac{x-1}{\sqrt{4t}}\bigg) \Bigg].$$

This is also the solution to the original problem ($x>0$). Plotting it for a few different values of $t$ gives a result that looks sensible:

Further comments are still welcome.

Answer 2

The Laplace transform approach also works, but your conclusion $a(s) = 0$ is wrong, since $\mathrm{e}^{\sqrt{s} x}$ is bounded for $x < 1$ anyway. The boundedness of $u$ can only be used to discard the term with the positive exponent in the region where $x > 1$. The general solution then reads $$\bar{u}(x,s) = \begin{cases} a(s) \mathrm{e}^{\sqrt{s} x} + b(s) \mathrm{e}^{-\sqrt{s} x} + \frac{10}{s} &, \, x<1 \\ c(s) \mathrm{e}^{-\sqrt{s} x} & , \, x>1 \end{cases} \, .$$ From the boundary conditions $\partial_1 \bar{u}(0,s) = 0$, $\bar{u}(1^+,s) = \bar{u}(1^-,s)$, $\partial_1 \bar{u}(1^+,s) = \partial_1 \bar{u}(1^-,s)$ we obtain $$ a(s) = b(s) = - \frac{\mathrm{e}^{-\sqrt{s}}}{\mathrm{e}^{\sqrt{s}}-\mathrm{e}^{-\sqrt{s}}} c(s) = - 5\frac{\mathrm{e}^{-\sqrt{s}}}{s} \, ,$$ so \begin{align} \bar{u}(x,s) &= \frac{10}{s} \begin{cases} 1 - \frac{1}{2} \left[\mathrm{e}^{\sqrt{s} (x-1)} + \mathrm{e}^{-\sqrt{s} (x+1)}\right] &, \, x \leq 1 \\ \frac{1}{2} \left[\mathrm{e}^{-\sqrt{s} (x-1)} - \mathrm{e}^{-\sqrt{s} (x+1)}\right] & , \, x>1 \end{cases} \\ &= \frac{10}{s} \left[\operatorname{H}(1-x) - \frac{1}{2} \mathrm{e}^{-\sqrt{s} (x+1)} + \frac{\operatorname{sgn}(x-1)}{2} \mathrm{e}^{-\sqrt{s} \lvert x-1\rvert}\right]. \end{align} Finally, since $$ \operatorname{\mathcal{L}}^{-1} \left(s \mapsto \frac{\mathrm{e}^{-a \sqrt{s}}}{s}\right) (t) = \operatorname{erfc}\left(\frac{a}{2 \sqrt{t}}\right) = 1 - \operatorname{erf}\left(\frac{a}{2 \sqrt{t}}\right) $$ holds for $a, t \geq 0$, we recover the result you derived using the Fourier transform: \begin{align} u(x,t) &= 5 \left[\operatorname{erf}\left(\frac{x+1}{2 \sqrt{t}}\right) - \operatorname{erf}\left(\frac{x-1}{2 \sqrt{t}}\right)\right] + 10 \left[\operatorname{H}(1-x) - \frac{1 - \operatorname{sgn}(x-1)}{2}\right] \\ &= 5 \left[\operatorname{erf}\left(\frac{x+1}{2 \sqrt{t}}\right) - \operatorname{erf}\left(\frac{x-1}{2 \sqrt{t}}\right)\right] . \end{align}