In this lecture, we study how to solve the semi-homogeneous heat problem:
$$
\begin{aligned}
&u_t=\alpha^2 u_{xx}+F(x,t),\ 0<x<L,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
-k_1u_x(0,t)+h_1u(0,t)&=0,\ t>0\\
k_2u_x(L,t)+h_2u(L,t)&=0,
\end{aligned}
\right.\\
&{\rm IC:}\ u(x,0)=\phi(x),\ 0<x<L\ \ \ \ \ (1)
\end{aligned}
$$

We consider a solution of the form
$$u(x,t)=\sum_{n=0}^\infty T_n(t)X_n(x)$$
In addition, we assume that
$$F(x,t)=\sum_{n=0}^\infty F_n(t)X_n(x)$$
Let $A_n:=T_n(0)$. Then
$$\phi(x)=\sum_{n=0}^\infty A_nX_n(x)$$
For each $n=0,1,2,\cdots$, $X_n(x)$ is the solution of the Sturm-Liouville problem:
$$
\begin{aligned}
&X_n^{\prime\prime}=-\lambda_n^2X_n\\
&{\rm BCs:}\ \left\{\begin{aligned}
-k_1X_n'(0)+h_1X_n(0)&=0\\
k_2X_n'(L)+h_2X_n(L)&=0
\end{aligned}
\right.\ \ \ \ \ (2)
\end{aligned}
$$
and $T_n(t)$ satisfies the first order ODE with initial condition:
$$\left\{\begin{aligned}
\dot{T}_n+\lambda_n^2\alpha^2T_n=F_n\\
T_n(0)=A_n.
\end{aligned}
\right.\ \ \ \ \ (3)$$
Then $T_n(t)$ is given by
$$
T_n(t)=\left[A_n+\int_0^t F_n(s)e^{\lambda_n^2\alpha^2 s}ds\right]e^{-\lambda_n^2\alpha^2 t}.
$$
Finally, $F_n(t)$ and $A_n$ are to be determined by
\begin{align*}
F_n(t)&=\frac{1}{L_n}\int_0^L F(x,t)X_n(x)dx\ \ \ \ \ (3)\\
A_n&=\frac{1}{L_n}\int_0^L\phi(x)X_n(x)dx\ \ \ \ \ (4)
\end{align*}
As a summary we have the recipe for solving the semi-homogeneous heat IBVP (1):

1. Find the eigenvalue $\lambda_n$ and the corresponding eigenfunction $X_n$ by solving the Sturm-Liouville problem (2)
2. Find the coefficients $F_n(t)$ and $A_n$ using the formulae (4) and (5) respectively.
3. For each $n=0,1,2,\cdots$, find $T_n(t)$ using the formula (3)
4. Put all things together as $$u(x,t)=\sum_{n=0}^\infty T_n(t)X_n(x)$$

Example. Solve the semi-homogeneous heat IBVP:
\begin{align*}
&u_t=\alpha^2 u_{xx}+b\sin\pi x,\ 0<x<1,\ 0<t<\infty\\
&{\rm BCs:}\ \left\{\begin{aligned}
u(0,t)&=0,\ 0<t<\infty\\
u(1,t)&=0,
\end{aligned}
\right.\\
&{\rm IC:}\ u(x,0)=1,\ 0<x<1
\end{align*}

Solution:

Step 1. The eigenfunctions are
$$X_n(x)=\sin n\pi x,\ n=1,2,3,\cdot.$$

Step 2. \begin{align*}
F_n(t)&=2b\int_0^1\sin\pi x\sin n\pi xdx\\
&=\left\{\begin{array}{ccc}
b & {\rm if} & n=1\\
0 & {\rm otherwise}.\end{array}\right.\\
A_n&=2\int_0^1\sin n\pi xdx\\
&=\frac{2[1-(-1)^n]}{n\pi}
\end{align*}
For $n=1,2,3,\cdots$, $A_{2n}=0$ and $A_{2n-1}=\frac{4}{2(n-1)\pi}$

Step 3. \begin{align*}
T_1(t)&=\left[\frac{4}{\pi}+\int_0^t be^{\alpha^2\pi^2s}ds\right]e^{-\alpha^2\pi^2t}\\
&=\frac{b}{\alpha^2\pi^2}+\left(\frac{4}{\pi}-\frac{b}{a^2\pi^2}\right)e^{-\alpha^2\pi^2t}
\end{align*}
For $n=2,3,\cdots,$
$$T_{2n-1}(t)=\frac{4}{(2n-1)\pi}e^{-\alpha^2(2n-1)^2\pi^2t}$$

Step 4. The solution is therefore
\begin{align*}
u(x,t)=\left[\frac{b}{\alpha^2\pi^2}+\left(\frac{4}{\pi}-\frac{b}{a^2\pi^2}\right)e^{-\alpha^2\pi^2t}\right]\sin\pi x\\
+\sum_{n=2}^\infty \frac{4}{(2n-1)\pi}e^{-\alpha^2(2n-1)^2\pi^2t}\sin(2n-1)\pi x
\end{align*}

The eventual temperature is
$$\lim_{t\to\infty}u(x,t)=\frac{b}{\alpha^2\pi^2}\sin\pi x$$

References:

David Betounes, Partial Differential Equations for Computational Science with Maple and Vector Analysis, TELOS, Springer-Verlag

Consider the heat IBVP:
\begin{align*}
&u_t=\alpha^2 u_{xx},\ 0<x<L,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
&u(0,t)=g_1(t),\ t>0\\
&u_x(L,t)+hu(L,t)=g_2(t),
\end{aligned}
\right.\\
&{\rm IC:}\ u(x,0)=\phi(x),\ 0<x<L.
\end{align*}
Again we seek a solution of the form:
$$u(x,t)=S(x,t)+U(x,t)$$
where $S(x,t)$ satisfies the non-homogeneous BCs:
\begin{align*}
S(0,t)&=g_1(t)\ \ \ \ \ (1)\\
S_x(L,t)+hS(L,t)&=g_2(t)
\end{align*}
Suppose that
$$S(x,t)=A(t)x+B(t)$$
From the BCs (1) we obtain
$$S(x,t)=\frac{g_2(t)-hg_1(t)}{1+hL}x+g_1(t)$$
$U(x,t)$ satisfies the semi-homogeneous heat IBVP:
\begin{align*}
&U_t=\alpha^2 U_{xx}-S_t,\ 0<x<L,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
&U(0,t)=0,\ t>0\\
&U_x(L,t)+hU(L,t)=0,
\end{aligned}
\right.\\
&{\rm IC:}\ U(x,0)=\phi(x)-S(x,0),\ 0<x<L.
\end{align*}
Since $u(x,t)\approx S(x,t)$ for large time, the eventual temperature is $\lim_{t\to\infty}S(x,t)$. The homogeneous BCs are pretty much similar to the ones we studied here but the heat equation is non-homogeneous. Such a problem (non-homogeneous heat equation + homogeneous BCs) is called a semi-homogeneous IBVP. Unfortunately, we cannot solve the semi-homogenous heat IBVP yet. We will discuss how to solve semi-homogeneous IBVPs later.

Example. Solve the heat IBVP:
\begin{align*}
&u_t=\alpha^2 u_{xx}-ae^{-t}x,\ 0<x<1,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
&u(0,t)=1,\ t>0\\
&u(1,t)=ae^{-t},
\end{aligned}
\right.\\
&{\rm IC:}\ u(x,0)=0,\ 0<x<1.
\end{align*}

Solution: Let $u(x,t)=S(x,t)+U(x,t)$
where $S(x,t)=A(t)x+B(t)$ satisfies the BCs. Then
$S(x,t)=(ae^{-t}-1)x+1$. $U(x,t)$ satisfies the homogeneous heat IBVP:
\begin{align*}
&U_t=\alpha^2 U_{xx},\ 0<x<1,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
&U(0,t)=0,\ t>0\\
&U(1,t)=0,
\end{aligned}
\right.\\
&{\rm IC:}\ U(x,0)=-(a-1)x-1,\ 0<x<1.
\end{align*}
We know how to solve this problem and the solution $U(x,t)$ is given by
$$U(x,t)=\sum_{n=1}^\infty A_n e^{-\alpha^2n^2\pi^2t}\sin n\pi x$$
where for each $n=1,2,\cdots$
\begin{align*}
A_n&=2\int_0^1[-(a-1)x-1]\sin n\pi xdx\\
&=2\frac{a(-1)^n-1}{n\pi}
\end{align*}
That is,
$$U(x,t)=2\sum_{n=1}^\infty\frac{a(-1)^n-1}{n\pi}e^{-\alpha^2n^2\pi^2t}\sin n\pi x$$
Therefore, the solution to the original heat problem is
$$u(x,t)=(ae^{-t}-1)x+1+2\sum_{n=1}^\infty\frac{[a(-1)^n-1]}{n\pi}e^{-\alpha^2 n^2\pi^2 t}\sin n\pi x.$$
The eventual temperature is
$$\lim_{t\to\infty}S(x,t)=-x+1$$

References:

David Betounes, Partial Differential Equations for Computational Science with Maple and Vector Analysis, TELOS, Springer-Verlag

Consider the following 1-dimensional heat IBVP:
\begin{align*}
&u_t=\alpha^2 u_{xx},\ 0<x<L,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
u(0,t)&=c_1,\ t>0\\
u(L,t)&=c_2,
\end{aligned}
\right.\\
&{\rm IC:}\ u(x,0)=\phi(x),\ 0<x<L.
\end{align*}
Since the BCs are not homogeneous, the method of separation of variables cannot be used. Let us assume that $u(x,t)$ can be decomposed to
$$u(x,t)=S(x,t)+U(x,t)$$
where $S(x,t)$ satisfies the non-homogeneous BCs i.e.
\begin{align*}
S(0,t)&=c_1\ \ \ \ \ (1)\\
S(L,t)&=c_2
\end{align*}
This implies that $U(x,t)$ satisfies the homogeneous BCs
$$U(0,t)=U(L,t)=0$$
We don’t know the function $S(x,t)$ but we can make a guess. We can first try a possibly simplest one, a linear function in $x$
$$S(x,t)=A(t)x+B(t)$$
If this works, it is our luck day. If not, try something else that is next simplest one, say a quadratic function in $x$. It turns out that our guess works! This means that we can determine the coefficients $A(t)$ and $B(t)$ from the BCs (1) and we find that
$$A(t)=\frac{c_2-c_1}{L},\ B(t)=c_1$$
Thus
$$S(x,t)=\frac{c_2-c_1}{L}x+c_1$$
$U(x,t)$ satisfies the homogeneous IBVP:
\begin{align*}
&U_t=\alpha^2 U_{xx},\ 0<x<L,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
U(0,t)&=0,\ t>0\\
U(L,t)&=0,
\end{aligned}
\right.\\
&{\rm IC:}\ U(x,0)=\phi(x)-S(x,0)=\phi(x)-\frac{c_2-c_1}{L}x-c_1,\ 0<x<L.
\end{align*}
We already know how to solve the homogeneous IBVP (see here) and
$$U(x,t)=\sum_{n=1}^\infty A_ne^{-\alpha^2\frac{n^2\pi^2}{L^2}t}\sin\frac{n\pi}{L}x$$
where
$$A_n=\frac{2}{L}\int_0^L[\phi(x)-\frac{c_2-c_1}{L}x-c_1]\sin\frac{n\pi}{L}xdx,\ n=1,2,\cdots$$

Note that $\lim_{t\to\infty}U(x,t)=0$ and so $\lim_{t\to\infty}u(x,t)=S(x)$. $U(x,t)$ is called the transient and $S(x)$ is called the steady state. The steady state is eventual state for large time.

Example. Solve the IBVP:
\begin{align*}
&u_t=\alpha^2 u_{xx},\ 0<x<1,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
u(0,t)&=c_1,\ t>0\\
u(1,t)&=c_2,
\end{aligned}
\right.\\
&{\rm IC:}\ u(x,0)=0,\ 0<x<1.
\end{align*}

Solution: First we solve
\begin{align*}
&U_t=\alpha^2 U_{xx},\ 0<x<1,\ t>0\\
&{\rm BCs:}\ \left\{\begin{aligned}
U(0,t)&=0,\ t>0\\
U(1,t)&=0,
\end{aligned}
\right.\\
&{\rm IC:}\ U(x,0)=-c_1-(c_2-c_1)x,\ 0<x<1.
\end{align*}
The transient $U(x,t)$ is given by
$$U(x,t)=\sum_{n=1}^\infty A_ne^{-\alpha^2n^2\pi^2 t}\sin n\pi x$$
where
\begin{align*}
A_n&=-2\int_0^1[c_1+(c_2-c_1)]\sin n\pi xdx\\
&=\frac{2}{n\pi}[c_2(-1)^n-c_1].
\end{align*}
Hence the solution is
$$u(x,t)=(c_2-c_1)x+c_1+2\sum_{n=1}^\infty\frac{[c_2(-1)^n-c_1]}{n\pi}e^{-\alpha^2n^2\pi^2 t}\sin n\pi x$$

References:

David Betounes, Partial Differential Equations for Computational Science with Maple and Vector Analysis, TELOS, Springer-Verlag

Let us consider the heat IBVP:
\begin{align*}
u_t=\alpha^2u_{xx},\ 0<x<1,\ t>0\\
u_x(0,t)=0,\ t>0\\
u_x(1,t)+u(1,t)=0,\ t>0\\
u(x,0)=1,\ 0<x<1.
\end{align*}
First we solve the Sturm-Liouville problem:
\begin{align*}
X^{\prime\prime}=kX\ \mbox{with BCs}\\
\left\{
\begin{aligned}
&X'(0)=0\\
&X'(1)+X(1)=0.
\end{aligned}
\right.
\end{align*}
The two cases $k=0$ and $k=\lambda^2>0$ are left to be checked for students. Note that the case $k=0$ leads to a trivial solution and for $k=\lambda^2>0$ case there is no solution. Let $k=-\lambda^2$. Then
$$X(x)=A\cos\lambda x+B\sin \lambda x$$
and
$$X'(x)=-A\lambda\sin\lambda x+B\lambda\cos\lambda x.$$
We obtain $B=0$ from the condition $X'(0)=0$. So $X'(x)=-A\lambda\sin\lambda x$. The BC $X'(1)+X(1)=0$ then implies that the eigenvalues $\lambda$ satisfies
$$\lambda=\cot\lambda\ \ \ \ \ (1)$$

The roots of (1) can be found by the Newton’s method. In Maxima, the Newton’s method can be performed as follows.

First we define $\cot(x)-x$ as a function $f(x)$:

(%i1) f(x):=cot(x)-x;

Newton’s method applied to f(x) is performed by

(%i2) newton(t):=subst(t,x,x-f(x)/diff(f(x),x));

To use Newton’s method we need the initial approximation. We can make a guess from its graph. As you can see the first zero is near 0.5, so we choose $x_0=0.5$.

If you just run

(%i3) newton(0.5);

Maxima will return the output

(%o3)                              0.74865744241706

In order to calculate next approximation $x_1$, simply run

(%i4) newton(%);

(%o4)                              0.85239848545338

For the next approximation $x_2$, again run

(%i5) newton(%);

(%o5)                              0.86029880244504

By running the same command, the next approximation $x_3$ is obtained as $x_3=0.86033358835819$. This value is accurate to eight decimal places. If you think this is good enough, you can stop.

Similarly we obtain $\lambda_n$ for $n=2,3,4,\cdots$:
$$\begin{array}{|c||c|c|}
\hline
n & \lambda_n & \frac{(2n-1)\pi}{2}\\
\hline
1 & .8603358 & 1.570796327\\
\hline
2 & 3.425618 & 4.712388981\\
\hline
3 & 6.437298 & 7.853981635\\
\hline
4 & 9.529334 & 10.99557429\\
\hline
5 & 12.645287 & 14.13716694\\
\hline
\cdots & \cdots & \cdots\\
\hline
\end{array}$$
The table shows that as $n$ increases $\frac{(2n-1)\pi}{2}$ gets closer to $\lambda_n$. So for a practical purpose one may use the approximation
$$\lambda_n\approx \frac{(2n-1)\pi}{2},\ \lambda_n=1,2,3,\cdots$$
The temperature distribution $u(x,t)$ is then given by
$$u(x,t)=\sum_{n=1}^\infty A_ne^{-\alpha^2\lambda_n^2t}\cos\lambda_n x$$
The $L_n$’s and $A_n$’s are computed to be:
\begin{align*}
L_n&=\int_0^1\cos^2\lambda_nxdx\\
&=\frac{1}{2}\int_0^1[\cos(2\lambda_nx)+1]dx\\
&=\frac{1}{2}\left[\frac{\sin 2\lambda_n}{2\lambda_n}+1\right]\\
&=\frac{1}{2}\left[\frac{\sin\lambda_n\cos\lambda_n}{\lambda_n}+1\right]\ (\sin2\lambda_n=2\sin\lambda_n\cos\lambda_n)\\
&=\frac{1}{2}\left[\frac{\sin\lambda_n\cos\lambda_n+\lambda_n}{\lambda_n}\right]
\end{align*}
and
\begin{align*}
A_n&=\frac{1}{L_n}\int_0^1\phi(x)\cos\lambda_nxdx\\
&=\frac{2\lambda_n}{\sin\lambda_n\cos\lambda_n+\lambda_n}\int_0^1\cos\lambda_nxdx\\
&=\frac{2\sin\lambda_n}{\sin\lambda_n\cos\lambda_n+\lambda_n}.
\end{align*}

This calculation is done by assuming that $\cos\lambda_n x$, $n=1,2,\cdots$ are orthogonal functions. In fact they are orthogonal functions if $\lambda_n=\frac{(2n-1)\pi}{2}$, $n=1,2,\cdots$.

Finally we have
$$u(x,t)=\sum_{n=1}^\infty\frac{2\sin\lambda_n}{\sin\lambda_n\cos\lambda_n+\lambda_n}e^{-\alpha^2\lambda_n^2 t}\cos\lambda_n x.$$
For an approximation one may take
$$u(x,t)\approx\sum_{n=1}^\infty\frac{4(-1)^{n-1}}{(2n-1)\pi}e^{-\alpha^2\frac{(2n-1)^2\pi^2}{4}t}\cos\frac{(2n-1)\pi x}{2}\ \ \ \ \ (2)$$

The 3D plot of $u(x,t)$ in (2) is



References:

David Betounes, Partial Differential Equations for Computational Science with Maple and Vector Analysis, TELOS, Springer-Verlag

Sturm-Liouville Problems

The homogeneous boundary conditions of 1D heat conduction problem are given by
\begin{align*}
-\kappa_1u_x(0,t)+h_1u(0,t)&=0,\ t>0\\
\kappa_2u_x(L,t)+h_2u(L,t)&=0,\ t>0
\end{align*}
(See here)

The homogeneous BCs for the second order linear differential equation \begin{equation}\label{eq:ho}X^{\prime\prime}=kX\end{equation} is then
\begin{equation}\label{eq:bc}\begin{aligned}
-k_1X'(0)+h_1X(0)&=0\\
k_2X'(L)+h_2X(L)&=0
\end{aligned}\end{equation}
Finding solutions of the second order linear differential equation \eqref{eq:ho} for $k=0$, $k=\lambda^2$, and $k=-\lambda^2$ that satisfy the BCs \eqref{eq:bc} is called a Sturm-Liouville Problem. Here, we study the Sturm-Liouville Theory with the following example.

Remark. In case of homogeneous heat BVPs, the eventual temperature would be 0 as there is no heat source. So, we see that $k=-\lambda^2<0$ is the only physically relevant case.

Example. [Fixed temperature at both ends]

Consider the heat BVP:
\begin{align*}
u_t&=\alpha^2 u_{xx}\ \mbox{PDE}\\
u(0,t)&=u(1,t)=0\ \mbox{(BCs)}
\end{align*}
From the above BCs, we obtain the BCs for $X(x)$:
$$X(0)=X(1)=0$$
For $k=0$ and $k=\lambda^2>0$ we have a trivial solution $X(x)=0$. For $k=-\lambda^2<0$ $X(x)=A\cos\lambda x+B\sin\lambda x$. With the BCS we find the eigenvalues
$$\lambda_n=n\pi,\ n=1,2,3,\cdots$$
and the corresponding eigenfunctions
$$X_n(x)=\sin n\pi x,\ n=1,2,3,\cdots$$
The $\{X_n: n=1,2,3,\cdots\}$ is a linearly independent set so they form a basis for the solution space which is infinite dimensional. The general solution to the heat BVP is given by
$$u(x,t)=\sum_{n=1}^\infty A_n e^{-n^2\pi^2\alpha^2t}\sin n\pi x$$
There are undetermined coefficients $A_n$ called Fourier coefficients. They can be determined by initial condition (initial temperature).

Orthogonal Functions and Solution of a Homogeneous Heat IBVP

Consider a heat distribution function $u(x,t)$ of the following form
$$u(x,t)=\sum_{n=0}^\infty A_ne^{-\lambda_n^2\alpha^2t}X_n(x)$$
where $X_n$’s are eigenfunctions corresponding to the eigenvalues $\lambda_n$’s respectively. The eigenfunctions $X_n$’s form a basis for the solution space (which is often infinite dimensional) of a given heat IBVP, furthermore they can form an orthonormal basis with respect to the inner product
\begin{equation}\label{eq:innerprod}\langle X_m,X_n\rangle=\int_0^LX_mX_ndx\end{equation}
We say that eigenfunctions $X_m$ and $X_n$ are orthogonal if $\langle X_m,X_n\rangle=0$.

Example. $X_n(x)=\sin n\pi x$, $n=1,2,3,\cdots$ form an orthogonal basis with respect to \eqref{eq:innerprod}, where $0<x<1$:
\begin{align*}
\langle X_m,X_n\rangle&=\int_0^1\sin m\pi x\sin n\pi xdx\\
&=\left\{\begin{aligned}
\frac{1}{2}\ &{\rm if}\ m=n\\
0\ &{\rm if}\ m\ne n.
\end{aligned}\right.
\end{align*}

Remark. [The Gram-Schmidt Orthogonalization Process]
If $\{X_n\}$ is not an orthogonal basis, one can construct an orthogonal basis from $\{X_n\}$ using the inner product (3). The standard process is called the Gram-Schmidt orthogonalization process. Details can be found in many standard linear algebra textbooks.

Now we assume that $\{X_n\}$ is an orthogonal basis for the solution space. Let $L_n:=\langle X_n,X_n\rangle=\int_0^LX_n^2dx$. Let the initial condition be given by
$u(x,0)=\phi(x)$. Then
$$\phi(x)=\sum_{n=0}^\infty A_nX_n$$
Multiply this by $X_m$ and then integrate:
$$\int_0^LX_m\phi(x)dx=\sum_{n=0}^\infty A_n\int_0^LX_nX_mdx$$
By orthogonality we obtain
$$L_mA_m=\int_0^LX_m\phi(x)dx$$
or
$$A_m=\frac{1}{L_m}\int_0^L\phi(x)X_mdx,\ m=0,1,2,\cdots$$

Example. Consider the heat BVP in the previous example with initial condition $\phi(x)=T$, a constant temperature. For $n=1,2,3,\cdots$, $X_n(x)=\sin n\pi x;\ 0<x<1$ so
$$L_n=\int_0^1\sin^2 n\pi x dx=\frac{1}{2}$$
The Fourier coefficients are then computed to be
\begin{align*}
A_n&=2\int_0^1\phi(x)\sin n\pi xdx\\
&=2T\int_0^1\sin n\pi xdx\\
&=\frac{2T}{n\pi}[1-\cos n\pi]\\
&=\frac{2T}{n\pi}[1-(-1)^n].
\end{align*}
$A_n=0$ for $n={\rm even}$ and $A_{2n-1}=\frac{4T}{(2n-1)\pi},\ n=1,2,3,\cdots$. Hence
$$u(x,t)=\sum_{n=1}^\infty\frac{4T}{(2n-1)\pi}e^{-(2n-1)^2\pi^2\alpha^2t}\sin(2n-1)\pi x.$$

References:

David Betounes, Partial Differential Equations for Computational Science with Maple and Vector Analysis, TELOS, Springer-Verlag