PlanetPhysics/Example of Solving the Heat Equation

From testwiki
Revision as of 03:28, 13 September 2020 by imported>MaintenanceBot (Formatting)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search

Let a thin square-formed plate of heat conducting homogeneous material be in the xy-plane with sides on the x-axis (isolated), on the line\, y=π (held at the constant temperature\, u=C), and on the vertical lines\, x=0\, and\, x=π (both held at the constant temperature\, u=0).\, Determine the [[../BoltzmannConstant/|temperature]] [[../Bijective/|function\,]] (x,y)u(x,y)\, on the plate, when the faces of the plate are isolated.

The [[../HeatEquation/|equation of the heat flow]] in this stationary case is

2uu'xx+u'yy=0

under the [[../PiecewiseLinear/|boundary]] conditions u(0,y)=0,u(π,y)=0,u(x,π)=C,u'y(x,0)=0. We first try to separate the variables, i.e. seek the solution of (1) of the form u(x,y):=X(x)Y(y). Then we get u'x=XY,u'xx=XY,u'y=XY,u'yy=XY, and thus (1) gets the form

XY+XY=0

and the boundary conditions X(0)=X(π)=0,X(x)=CY(π),Y(0)=0. We separate the variables in (2): XX=YY This equation is not possible unless both sides are equal to a same negative constant k2, which implies for\, X=k2X\, the solution X:=C1coskx+C2sinkx and for\, Y=k2Y\, the solution Y:=D1coshky+D2sinhky. The two first boundary conditions give\, 0=X(0)=C1,\, 0=X(π)=0+C2sinkπ,\, and since\, C20,\, we must have\, sinkπ=0,\, i.e. 0<k:=n=1,2,3, Therefore X(x):=C2sinnx,Y(y)nD1sinhny+nD2coshny. The fourth boundary condition now gives that\, 0=Y(0)=nD2;\, thus\, D2=0\, and\, Y(y):=D1coshny.\, So (1) has infinitely many solutions

un:=C2D1sinnxcoshny=Ansinnxcoshny

with\, n+\, and they all satisfy the boundary conditions except the third.\, Because of the linearity of (1), also the sum u:=n=1Ansinnxcoshny of the functions (3) satisfy (1) and those boundary conditions, provided that this series converges.\, The third boundary condition requires that C=u(x,π)=n=1Ansinnxcoshnπ=n=1(Ancoshnπ)sinhnx on the interval\, 0xπ.\, But this is the Fourier sine series of the constant function\, xC\, on the half-interval\, [0,π],\, whence Ancoshnπ=2π0πCsinnxdx=2Cnπ(1(1)n)n+. The even n's here give 0 and the odd give A2m+1:=4C(2m+1)πcosh(2m+1)π(m=0,1,2,)

Thus we obtain the solution u(x,y):=4Cπm=0sin(2m+1)xcosh(2m+1)y(2m+1)cosh(2m+1)π. It can be shown that this series converges in the whole square of the plate.

\subsection*{Visualization of the solution}

\begin{figure}[!htb]

\includegraphics{heat-surface.eps} \caption{Surface plot of the solution u(x,y), for C=1}

\end{figure}

\begin{figure}[!htb]

\includegraphics{heat-color.eps} \caption{Color-coded plot of the temperature u(x,y)}

\end{figure}

Remark. The function u has been approximated in the plot by computing a partial sum of the true infinite-series solution.\, However, there is substantial numerical error in the approximate solution near\, y=π,\, evident in the small oscillations observed in the surface plot, that should not be there in theory.\, This phenomenon is actually inevitable given that the boundary conditions are actually discontinuous at the corners\, (0,π)\, and\, (π,π).

More precisely, observe that when\, y=π,\, the formula for\, u(x,y)\, reduces to the Fourier series

4Cπ(sinx+sin3x3+sin5x5+)

for the discontinuous function on\, [π,π]:

x{C,0<x<πC,π<x<0

That means the Fourier expansion will necessarily be subject to the Gibbs phenomenon .\, Of course, the series also cannot converge absolutely; in other words, the terms of the series decay too slowly in [[../AbsoluteMagnitude/|magnitude]], adversely affecting the numerical solution.

 \item Python program to compute\, u(x,y)\, and produce the two figures

Template:CourseCat