Уравнение теплопроводности

Пример численного решения уравнения теплопроводности. Цветом и высотой поверхности передана температура данной точки.

Уравнение теплопроводности — важное уравнение в частных производных, которое описывает распространение тепла в заданной области пространства во времени.

Вид уравнения

Для функции u(x,y,z,t) трёх пространственных переменных (x,y,z) и времени t, уравнение теплопроводности имеет вид

$\frac{\partial u}{\partial t} -\alpha\left(\frac{\partial^2u}{\partial x^2}+\frac{\partial^2u}{\partial y^2}+\frac{\partial^2u}{\partial z^2}\right)=0$

Для произвольной системы координат:

$\frac{\partial u}{\partial t} - \alpha \nabla^2 u=0$

где α — положительная константа, а Δ или ∇² — оператор Лапласа.

Способы решения уравнений теплопроводности

Пример численного решения уравнения теплопроводности. Цветом и высотой поверхности передана температура данной точки.

Метод разделения переменных (Метод Фурье)

Однородное уравнение теплопроводности с однородными граничными условиями

Рассмотрим следующую задачу

$\begin{array}{l} u_t=a^2 u_{xx},\quad 0<x<l,\;0<t\leqslant T \\ u(x,\;0)=\varphi(x);\quad 0\leqslant x\leqslant l \\ \left.\begin{array}{l} u(0,\;t)=0, \\ u(l,\;t)=0. \\ \end{array}\right\}\quad 0\leqslant t\leqslant T \\ \end{array}$

Требуется найти функцию $u(x,\;t)$ для $\forall(x,\;t):0\leqslant x\leqslant l,\;0\leqslant t\leqslant T$.

Представим искомую функцию в виде произведения

$u(x,\;t)=X(x)T(t).$

Затем предполагаемую форму решения подставим в исходное уравнение, получим

$X(x)T'(t)=a^2 X''(x)T(t).$

Разделим выражение на $a^2 X(x)T(t)$:

$\frac{1}{a^2}\frac{T'(t)}{T(t)}=\frac{X''(x)}{X(x)}=-\lambda,\;\lambda=\mathrm{const}.$

Так как в левой части уравнения у нас находится функция зависящая только от $t$, а в правой — только от $x$, то, фиксируя любое значение $x$ в правой части, получаем, что для любого $t$ значение левой части уравнения постоянно. Таким же образом можно убедиться, что и правая часть постоянна, то есть равна некой константе $-\lambda$ (минус взят для удобства). Таким образом, мы получаем два обыкновенных линейных дифференциальных уравнения:

$\begin{array}{l} X''(x)+\lambda X(x) = 0, \\ T'(t)+a^2\lambda T(t)=0. \end{array}$

Обратим внимание на граничные условия исходной задачи и подставим в них предполагаемый вид уравнения, получим:

$\begin{array}{l} u(0,\;t)=X(0)T(t)=0, \\ u(l,\;t)=X(l)T(t)=0, \end{array}$

откуда $X(0)=X(l)=0$ ($T(t)\ne 0$, так как в противном случае мы имели бы решение $u(x,\;t)=0$, а мы ищем только нетривиальные решения).

С учетом полученных граничных условий мы получаем задачу Штурма — Лиувилля:

$\begin{array}{l} X''(x)+\lambda X(x)=0; \\ X(0)=0, \\ X(l)=0. \\ \end{array}$

Её решение сводится к решению линейного дифференциального уравнения и рассмотрению трёх случаев:

1. $\lambda<0.$

   В этом случае общий вид решения будет следующим:

   $X(x)=C_1 e^{\sqrt{-\lambda}x}+C_2 e^{-\sqrt{-\lambda}x}.$

   Подставив граничные условия, мы убедимся, что решение будет $X(x)\equiv 0$, а мы ищем только нетривиальные решения, следовательно, этот случай не подходит.

2. $\lambda=0.$

   Общий вид решения

   $X(x)=C_1 x+C_2.$

   Несложно убедиться, что этот вариант нам также не подходит.

3. $\lambda>0.$

   Общий вид решения

   $X(x)=C_1\cos(\sqrt\lambda x)+C_2\sin(\sqrt\lambda x).$

   Подставим граничные условия:

   $\begin{array}{l} X(0)=C_1=0, \\ X(l)=C_2\sin(\sqrt\lambda l)=0. \end{array}$

   Так как мы ищем только нетривиальные решения, $C_2=0$ нам не подходит, следовательно

   $\begin{array}{l} \sin(\sqrt\lambda l)=0, \\ \sqrt\lambda l=\pi n,\quad n=1,\;2,\;\ldots \\ \end{array}$

   $\lambda_n=\left(\frac{\pi n}{l}\right)^2,\quad n=1,\;2,\;\ldots$

   Отсюда

   $X_n(x)=C_n\sin\left(\frac{\pi n}{l}x\right),\;\quad n=1,\;2,\;\ldots$

C учетом найденных $\lambda$, выведем общее решение линейного дифференциального уравнения

$T'(t)+a^2\left(\frac{\pi n}{l}\right)^2 T(t)=0.$

Должен получиться ответ

$T_n(t)=D_n\exp\left(-a^2\left(\frac{\pi n}{l}\right)^2 t\right),\quad D_n=\mathrm{const}.$

Теперь всё готово для того, чтобы записать решение исходной задачи:

$u_n(x,\;t)=X_n(x)T_n(t)=C_n\sin\left(\frac{\pi n}{l}x\right)\exp\left(-a^2\left(\frac{\pi n}{l}\right)^2 t\right),\quad n=1,\;2,\;\ldots$

В результате у нас получилось бесконечное количество частных решений уравнения. Все эти частные решения линейно независимы, то есть линейная комбинация любого количества решений равна нулю, только если все коэффициенты при них равны нулю. Поэтому логично предположить, что суммируя все частные решения по $n$ от единицы до бесконечности, мы получим общее решение исходной задачи.

$u(x,\;t)=\sum\limits_{n=1}^\infty u_n(x)=\sum\limits_{n=1}^\infty C_n\sin\left(\frac{\pi n}{l}x\right)\exp\left(-a^2\left(\frac{\pi n}{l}\right)^2 t\right).$

Осталось определить значение константы $C$ (зависящей от $n$) из начального условия

$u(x,\;0)=\varphi(x).$

Для того, чтобы определить значение $C_n$, необходимо разложить функцию $\varphi(x)$ в ряд Фурье:

$\begin{array}{l} \varphi(x)=\sum\limits_{n=1}^\infty A_n\sin\left(\dfrac{\pi n}{l}x\right), \\ A_n=\dfrac{2}{l}\displaystyle\int\limits_0^l \varphi(\xi)\sin\left(\dfrac{\pi n}{l}\xi\right)\,d\xi. \end{array}$

Получаем:

$\begin{array}{l} u(x,\;0)=\sum\limits_{n=1}^\infty C_n\sin\left(\dfrac{\pi n}{l}x\right)=\sum\limits_{n=1}^\infty A_n\sin\left(\dfrac{\pi n}{l}x\right), \\ C_n=A_n=\dfrac{2}{l}\displaystyle\int\limits_0^l \varphi(\xi)\sin\left(\dfrac{\pi n}{l}\xi\right)\,d\xi. \end{array}$

Откуда общее решение:

$u(x,\;t)=\sum\limits_{n=1}^\infty \left(\dfrac{2}{l}\int\limits_0^l \varphi(\xi)\sin\left(\dfrac{\pi n}{l}\xi\right)\,d\xi\right)\sin\left(\dfrac{\pi n}{l}x\right)\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right).$

В курсе математической физики доказывается, что полученный ряд удовлетворяет всем условиям данной задачи, то есть функция $u(x,\;t)$ дифференцируема (и ряд сходится равномерно), удовлетворяет уравнению в области определения и непрерывна в точках границы этой области.

Неоднородное уравнение теплопроводности с однородными граничными условиями

Рассмотрим способ решения неоднородного уравнения:

$\begin{array}{l} u_t=a^2 u_{xx}+f(x,\;t),\quad 0<x<l,\;0<t\leqslant T \\ u(x,\;0)=0;\quad 0\leqslant x\leqslant l \\ \left.\begin{array}{l} u(0,\;t)=0, \\ u(l,\;t)=0. \\ \end{array} \right\}\quad 0\leqslant t\leqslant T \end{array}$

Пусть

$\begin{array}{l} u_n(x,\;t)=X_n(x)T_n(t), \\ f_n(x,\;t)=X_n(x)F_n(t), \\ X_n(x)=\sin\left(\dfrac{\pi n}{l}x\right). \end{array}$

Тогда, пользуясь очевидным соотношением $X''_n(x)=-\left(\frac{\pi n}{l}\right)^2 X_n(x)$, перепишем исходное уравнение как:

$\begin{array}{l} X_n(x)T'_n(t)=-a^2\left(\dfrac{\pi n}{l}\right)^2 X_n(x)T_n(t)+X_n(x)F_n(t), \\ T'_n(t)=-a^2\left(\dfrac{\pi n}{l}\right)^2 T_n(t)+F_n(t). \end{array}$

Решим последнее линейное неоднородное уравнение методом вариации постоянной. Сначала найдём общее решение однородного линейного уравнения

$\begin{array}{l} T'_n(t)=-a^2\left(\dfrac{\pi n}{l}\right)^2 T_n(t), \\ T_n(t)=D\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right). \end{array}$

В общем решении заменим постоянную $D$ на переменную $D(t)$ и подставим в исходное уравнение.

$\begin{array}{l} T_n(t)=D(t)\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right), \\ D'_n(t)\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)-a^2\left(\dfrac{\pi n}{l}\right)^2\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)D_n(t)=-a^2\left(\dfrac{\pi n}{l}\right)^2\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)D_n(t)+F_n(t), \\ D'_n(t)\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)=F_n(t), \\ D_n(t)=\displaystyle\int F_n(t)\exp\left(a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)\,dt+A_n, \\ T_n(t)=A_n\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)+\exp\left(-a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)\displaystyle\int F_n(t)\exp\left(a^2\left(\dfrac{\pi n}{l}\right)^2 t\right)\,dt. \end{array}$

Из начального условия получаем:

$\begin{array}{l} u_n(x,\;0)=X_n(x)T_n(0)=0, \\ T_n(0)=0. \end{array}$

С учетом условия для $T$, получаем

$T_n(t)=\int\limits_0^t \exp\left(-a^2\left(\frac{\pi n}{l}\right)^2 (t-\tau)\right)F_n(\tau)\,d\tau.$

Так как

$f_n(x,\;t)=X_n(x)F_n(t)=\sin\left(\frac{\pi n}{l}x\right)F_n(t),$

то $F_n(t)$, очевидно, является коэффициентом ряда Фурье, и равен

$F_n(t)=\frac{2}{l}\int\limits_0^l f(\xi,\;t)\sin\left(\frac{\pi n}{l}\xi\right)\,d\xi.$

В результате, общая формула такова:

$u(x,\;t)=\sum\limits_{n=1}^\infty X_n(x)T_n(t)=\sum\limits_{n=1}^\infty \left[\int\limits_0^t \exp\left(-a^2\left(\frac{\pi n}{l}\right)^2 (t-\tau)\right)\left\{\frac{2}{l}\int\limits_0^l f(\xi,\;\tau)\sin\left(\frac{\pi n}{l}\xi\right)\,d\xi\right\}\,d\tau\right]\sin\left(\frac{\pi n}{l}x\right).$

Общая первая краевая задача

Во многих случаях удаётся решить уравнение теплопроводности с неоднородными краевыми и начальным условиями

$\begin{array}{l} u_t=a^2 u_{xx}+f(x,\;t), \\ u(x,\;0)=\varphi(x), \\ u(0,\;t)=\mu_1(t), \\ u(l,\;t)=\mu_2(t) \end{array}$

с помощью методов, описанных выше и следующего несложного приёма. Представим искомую функцию в виде суммы:

$\begin{array}{l} u(x,\;t)=\tilde u(x,\;t)+U(x,\;t), \\ \tilde u(x,\;0)=u(x,\;0)-U(x,\;0)=\varphi(x)-U(x,\;0), \\ \tilde u(0,\;t)=0, \\ \tilde u(l,\;t)=0. \end{array}$

Найдём функцию $U(x,\;t)$:

$\begin{array}{l} U(x,\;t)=Ax+b, \\ U(0,\;t)=b=\mu_1(t), \\ U(l,\;t)=Al+\mu_1=\mu_2\Rightarrow A=\dfrac{\mu _2(t)-\mu_1(t)}{l}, \\ U(x,\;t)=\dfrac{\mu_2(t)-\mu_1(t)}{l}x+\mu_1(t). \end{array}$

Таким образом, исходная задача свелась к следующей:

$\begin{array}{l} \tilde u_t=a^2\tilde u_{xx}+f(x,\;t)-\dfrac{\mu'_2(t)-\mu'_1(t)}{l}x-\mu'_1(t), \\ \tilde u(x,\;0)=\varphi(x)-\dfrac{\mu_2(0)-\mu_1(0)}{l}x-\mu_1(0), \\ \tilde u(0,\;t)=0, \\ \tilde u(l,\;t)=0. \end{array}$

После того, как мы найдём функцию $\tilde u(x,\;t)$, искомую функцию найдём по формуле

$u(x,\;t)=\tilde u(x,\;t)+\frac{\mu_2-\mu_1}{l}x+\mu_1.$

Литература

• Cannon, John Rozier (1984), «The One–Dimensional Heat Equation», vol. 23 (1st ed.), Encyclopedia of Mathematics and Its Applications, Reading-Menlo Park–London–Don Mills–Sidney–Tokyo/ Cambridge–New York–New Rochelle–Melbourne–Sidney: Addison-Wesley Publishing Company/Cambridge University Press, сс. XXV+483, ISBN 978-0-521-30243-2, <http://books.google.com/?id=XWSnBZxbz2oC&printsec=frontcover#v=onepage&q=> .

• Crank, J.; Nicolson, P. & Hartree, D. R. (1947), "«A Practical Method for Numerical Evaluation of Solutions of Partial Differential Equations of the Heat-Conduction Type»", Proceedings of the Cambridge Philosophical Society Т. 43: 50–67, DOI 10.1017/S0305004100023197 

• Einstein, Albert (1905), "«Über die von der molekularkinetischen Theorie der Wärme geforderte Bewegung von in ruhenden Flüssigkeiten suspendierten Teilchen»", Ann. Phys. Leipzig 17 Т. 322 (8): 549–560, DOI 10.1002/andp.19053220806 

• Evans, L.C. (1998), «Partial Differential Equations», American Mathematical Society, ISBN 0-8218-0772-2 

• John, Fritz (1991), «Partial Differential Equations» (4th ed.), Springer, ISBN 978-0-387-90609-6 

• Wilmott, P.; Howison, S. & Dewynne, J. (1995), «The Mathematics of Financial Derivatives:A Student Introduction», Cambridge University Press 

• Carslaw, H. S. & Jaeger, J. C. (1959), «Conduction of Heat in Solids» (2nd ed.), Oxford University Press, ISBN 978-0-19-853368-9 

• Thambynayagam, R. K. M. (2011), «The Diffusion Handbook: Applied Solutions for Engineers», McGraw-Hill Professional, ISBN 978-0-07-175184-1 

• Perona, P & Malik, J. (1990), "«Scale-Space and Edge Detection Using Anisotropic Diffusion»", IEEE Transactions on Pattern Analysis and Machine Intelligence Т. 12 (7): 629–639 

• Unsworth, J. & Duarte, F. J. (1979), "«Heat diffusion in a solid sphere and Fourier Theory»", Am. J. Phys. Т. 47 (11): 891–893, DOI 10.1119/1.11601 

Ссылки

	Уравнение теплопроводности на Викискладе?

• Вывод уравнения теплопроводности
• Linear heat equations: Particular solutions and boundary value problems — from EqWorld