Трёхмерная краевая задача теплопроводности

Рассмотрим трёхмерную несимметричную задачу теплопроводности в шаре радиуса $R$. Уравнение теплопроводности (1.8$\displaystyle \frac{\partial T(M, t)}{\partial t}$ $\displaystyle = a^2 \cdot \Delta T(M, t)$ $\displaystyle + f(M, t),$) в сферических координатах принимает вид

где $r$ $\in (0,R)$, $\theta$ $\in (0,\pi)$, $\phi$ $\in (0,2\pi)$, $R$ — радиус, $T(r, \theta, \phi, t)$ — температура, $a^2$ — коэффициент температуропроводности, $f(r, \theta, \phi, t)$ — функция плотности тепловых источников.

Краевые условия Дирихле (1.13$\displaystyle T(M, t)$ $\displaystyle = \Phi(M, t), \quad M$ $\displaystyle \in S.$) в трёхмерном случае имеют вид

где $T_0(r, \theta, \phi)$ — начальное условие, $\Phi(\theta, \phi, t)$ — граничное условие, которое, разумеется, периодично по $\phi$; на полюсах ($\theta$ $= 0$ и $\theta$ $= \pi$) вместо граничных условий ставится условие ограниченности решения.

Для избавления от неоднородности введём $T(r, \theta, \phi, t)$ $= \widehat{T}(r, \theta, \phi, t)$ $+ U(\theta, \phi, t)$, выберем функцию $U(\theta, \phi, t)$ $= \Phi(\theta, \phi, t)$, в соответствии с (2.2$\displaystyle \widehat{T}(M, t)$ $\displaystyle + U(M, t)$ $\displaystyle = \Phi(M, t), \quad \widehat{T}(M, t)$ $\displaystyle = 0 \;$ $\displaystyle \Rightarrow\; U(M, t)$ $\displaystyle = \Phi(M, t), \quad M$ $\displaystyle \in S.$). Тогда уравнение (2.71$\displaystyle \frac{\partial T(r, \theta, \phi, t)}{\partial t}$ $\displaystyle = a^2 \cdot \frac{\partial^2 T(r, \theta, \phi, t)}{\partial r^2}$ $\displaystyle + \frac{2 \cdot a^2}{r} \cdot \frac{\partial T(r, \theta, \phi, t)}{\partial r}$ $\displaystyle + \frac{a^2}{r^2 \cdot \sin(\theta)} \cdot \frac{\partial}{\partial \theta} \left( \sin(\theta) \cdot \frac{\partial T(r, \theta, \phi, t)}{\partial \theta} \right)$ $\displaystyle + \frac{a^2}{r^2 \cdot \sin^2 \theta} \cdot \frac{\partial^2 T(r, \theta, \phi, t)}{\partial \phi^2}$ $\displaystyle + f(r, \theta, \phi, t),$) и краевые условия (2.72$\displaystyle T(r, \theta, \phi, 0)$ $\displaystyle = T_0(r, \theta, \phi),$$\displaystyle T(R, \theta, \phi, t)$ $\displaystyle = \Phi(\theta, \phi, t),$$\displaystyle |T(r, 0, \phi, t)|$ $\displaystyle < \infty, \quad |T(r, \pi, \phi, t)|$ $\displaystyle < \infty,$$\displaystyle T(r, \theta, \phi, t)$ $\displaystyle = T(r, \theta, \phi + 2\pi, t),$) примут вид

где

Краевые условия Неймана (1.14$\displaystyle \lambda \cdot \frac{\partial T(M, t)}{\partial \overrightarrow{n}}$ $\displaystyle = \Phi(M, t), \quad M$ $\displaystyle \in S,$) в трёхмерном случае имеют вид

где $\frac{\displaystyle \partial T(r, \theta, \phi, t)}{\displaystyle \partial \mathbf{n}} \bigg|_{r=R}$ — нормальная производная температуры на сфере радиуса $R$; направление нормали совпадает с направлением координаты $r$.

Для избавления от неоднородности введём $T(r, \theta, \phi, t)$ $= \widehat{T}(r, \theta, \phi, t)$ $+ U(r, \theta, \phi, t)$, выберем функцию $U(r, \theta, \phi, t)$ $= \Phi(\theta, \phi, t) \cdot r$ — она удовлетворяет условию (2.3$\displaystyle \frac{\partial \widehat{T}(M, t)}{\partial \overrightarrow{n}}$ $\displaystyle + \frac{\partial U(M, t)}{\partial \overrightarrow{n}}$ $\displaystyle = \Phi(M, t), \quad \frac{\partial \widehat{T}(M, t)}{\partial \overrightarrow{n}}$ $\displaystyle = 0 \;$ $\displaystyle \Rightarrow\; \frac{\partial U(M, t)}{\partial \overrightarrow{n}}$ $\displaystyle = \Phi(M, t), \quad M$ $\displaystyle \in S.$): $\frac{\displaystyle \partial U(r, \theta, \phi, t)}{\displaystyle \partial \mathbf{n}} \bigg|_{r=R}$ $= \Phi(\theta, \phi, t)$. Тогда уравнение (2.71$\displaystyle \frac{\partial T(r, \theta, \phi, t)}{\partial t}$ $\displaystyle = a^2 \cdot \frac{\partial^2 T(r, \theta, \phi, t)}{\partial r^2}$ $\displaystyle + \frac{2 \cdot a^2}{r} \cdot \frac{\partial T(r, \theta, \phi, t)}{\partial r}$ $\displaystyle + \frac{a^2}{r^2 \cdot \sin(\theta)} \cdot \frac{\partial}{\partial \theta} \left( \sin(\theta) \cdot \frac{\partial T(r, \theta, \phi, t)}{\partial \theta} \right)$ $\displaystyle + \frac{a^2}{r^2 \cdot \sin^2 \theta} \cdot \frac{\partial^2 T(r, \theta, \phi, t)}{\partial \phi^2}$ $\displaystyle + f(r, \theta, \phi, t),$) и краевые условия (2.74$\displaystyle T(r, \theta, \phi, 0)$ $\displaystyle = T_0(r, \theta, \phi),$$\displaystyle \frac{\displaystyle \partial T(r, \theta, \phi, t)}{\displaystyle \partial \mathbf{n}} \bigg|_{r=R}$ $\displaystyle = \Phi(\theta, \phi, t),$$\displaystyle |T(r, 0, \phi, t)|$ $\displaystyle < \infty, \quad |T(r, \pi, \phi, t)|$ $\displaystyle < \infty,$$\displaystyle T(r, \theta, \phi, t)$ $\displaystyle = T(r, \theta, \phi + 2\pi, t),$) примут вид

где

Уравнение (2.8$\displaystyle \Delta \Psi(M)$ $\displaystyle + \gamma^2 \cdot \Psi(M)$ $\displaystyle = 0.$) примет вид

Представим функцию $\Psi(r, \theta, \phi)$ в виде произведения трёх функций $\Psi(r, \theta, \phi)$ $= R(r) \cdot \Theta(\theta) \cdot \Upsilon(\phi)$, подставим в уравнение (2.76$\displaystyle \frac{\displaystyle \partial^2 \Psi(r, \theta, \phi)}{\displaystyle \partial r^2}$ $\displaystyle + \frac{\displaystyle 2}{\displaystyle r} \cdot \frac{\displaystyle \partial \Psi(r, \theta, \phi)}{\displaystyle \partial r}$ $\displaystyle + \frac{\displaystyle 1}{\displaystyle r^2 \cdot \sin(\theta)} \cdot \frac{\displaystyle \partial}{\displaystyle \partial \theta} \left( \sin(\theta) \cdot \frac{\displaystyle \partial \Psi(r, \theta, \phi)}{\displaystyle \partial \theta} \right)$ $\displaystyle + \frac{\displaystyle 1}{\displaystyle r^2 \cdot \sin^2 \theta} \cdot \frac{\displaystyle \partial^2 \Psi(r, \theta, \phi)}{\displaystyle \partial \phi^2}$ $\displaystyle + \gamma^2 \cdot \Psi(r, \theta, \phi)$ $\displaystyle = 0.$) и проведём преобразования

Выражение в первой скобке зависит только от $r$, во второй — только от углов, поэтому, поделив на $R(r) \cdot \Theta(\theta) \cdot \Upsilon(\phi)$, обе части можно приравнять константе $\gamma_1^2$ — по аналогии с разделением переменных для времени и геометрии. Уравнение для радиуса примет вид

и продолжим преобразования

Здесь выражение в скобках зависит только от $\theta$, а последнее слагаемое — только от $\phi$, поэтому, домножив на $\frac{\displaystyle \sin^2(\theta)}{\displaystyle \Theta(\theta) \cdot \Upsilon(\phi)}$, снова разделим переменные константой $m^2$ и получим уравнения для полярного и азимутального углов

Решение уравнения (2.79$\displaystyle \frac{\displaystyle d^2 \Upsilon(\phi)}{\displaystyle d\phi^2}$ $\displaystyle + m^2 \cdot \Upsilon(\phi)$ $\displaystyle = 0.$) мы уже находили для двумерного случая, поэтому сошлёмся на (2.54$\displaystyle \Upsilon_{1m}(\theta)$ $\displaystyle = C_{1m} \cdot \sin(m \cdot \theta), m$ $\displaystyle \in (1..\infty),$$\displaystyle \Upsilon_{2m}(\theta)$ $\displaystyle = C_{2m} \cdot \cos(m \cdot \theta), m$ $\displaystyle \in (0..\infty).$). А для решения уравнения (2.78$\displaystyle \frac{\displaystyle 1}{\displaystyle \sin(\theta)} \cdot \frac{\displaystyle d}{\displaystyle d\theta} \left( \sin(\theta) \cdot \frac{\displaystyle d \Theta(\theta)}{\displaystyle d\theta} \right)$ $\displaystyle + \left[ \gamma_1^2 - \frac{\displaystyle m^2}{\displaystyle \sin^2(\theta)} \right] \cdot \Theta(\theta)$ $\displaystyle = 0,$) введём замену $z$ $= \cos(\theta)$: при $\theta$ $\in (0, \pi)$ имеем $z$ $\in (-1, 1)$ и $dz$ $= - \sin(\theta) \cdot d\theta$. Тогда уравнение примет вид

Это уравнение разобрано в приложении «Уравнение Лежандра», его решение имеет вид (E.6$\displaystyle \Theta_{km}(z)$ $\displaystyle = P_k^{(m)}(z)$ $\displaystyle = (1 - z^2)^{\frac{\displaystyle m}{\displaystyle 2}} \cdot \frac{\displaystyle d^m}{\displaystyle dz^m} P_k(z),$)

где $\gamma_{1k}^2$ $= k \cdot (k + 1)$ — собственные значения, $P_k^{(m)}(z)$ — присоединённые полиномы Лежандра, причём нетривиальные решения существуют только при $m$ $\leq k$. Норма присоединённого полинома Лежандра имеет вид (F.6$\displaystyle \int_{-1}^1 \left[ P_k^{(m)}(z) \right]^2 \,dz$ $\displaystyle = \frac{\displaystyle 2}{\displaystyle 2 \cdot k + 1} \cdot \frac{\displaystyle (k+m)!}{\displaystyle (k-m)!}, \quad m$ $\displaystyle \leq k.$). Теперь вернёмся к углу $\theta$, сделав обратную подстановку

Вообще говоря, это решение определено с точностью до постоянного множителя — положим его равным единице: он всё равно вбирается в коэффициенты разложения. Теперь вернёмся к уравнению (2.77$\displaystyle r^2 \cdot \frac{\displaystyle d^2 R(r)}{\displaystyle dr^2}$ $\displaystyle + 2 \cdot r \cdot \frac{\displaystyle d R(r)}{\displaystyle dr}$ $\displaystyle + \left( r^2 \cdot \gamma^2 - \gamma_1^2 \right) \cdot R(r)$ $\displaystyle = 0,$) и перепишем его с учётом $\gamma_1^2$ $= k \cdot (k + 1)$

Это уравнение разобрано в приложении «Уравнение Бесселя в сферических координатах»: подстановка $R(r)$ $= \frac{\displaystyle \widehat{R}(r)}{\displaystyle \sqrt{r}}$ сводит его к уравнению Бесселя с полуцелым индексом, и решение имеет вид (G.2$\displaystyle R(r)$ $\displaystyle = \frac{\displaystyle 1}{\displaystyle \sqrt r} \cdot J_{k + 1/2}(\gamma \cdot r),$). Следует заметить, что решение ограничено в точке $r$ $= 0$: функция Бесселя $J_{k + 1/2}(\gamma \cdot r)$ стремится к нулю не медленнее, чем $\sqrt{r}$.

Таким образом, решения уравнения (2.76$\displaystyle \frac{\displaystyle \partial^2 \Psi(r, \theta, \phi)}{\displaystyle \partial r^2}$ $\displaystyle + \frac{\displaystyle 2}{\displaystyle r} \cdot \frac{\displaystyle \partial \Psi(r, \theta, \phi)}{\displaystyle \partial r}$ $\displaystyle + \frac{\displaystyle 1}{\displaystyle r^2 \cdot \sin(\theta)} \cdot \frac{\displaystyle \partial}{\displaystyle \partial \theta} \left( \sin(\theta) \cdot \frac{\displaystyle \partial \Psi(r, \theta, \phi)}{\displaystyle \partial \theta} \right)$ $\displaystyle + \frac{\displaystyle 1}{\displaystyle r^2 \cdot \sin^2 \theta} \cdot \frac{\displaystyle \partial^2 \Psi(r, \theta, \phi)}{\displaystyle \partial \phi^2}$ $\displaystyle + \gamma^2 \cdot \Psi(r, \theta, \phi)$ $\displaystyle = 0.$) имеют вид

Подставим решения (2.82$\displaystyle \begin{aligned} &\Psi_{1km}(r, \theta, \phi) = A_{km} \cdot \frac{\displaystyle 1}{\displaystyle \sqrt{r}} \cdot J_{k + 1/2}(\gamma \cdot r) \cdot P_k^{(m)}(\cos(\theta)) \cdot \sin(m \cdot \phi),\\ &1 \leq m \leq k, \quad k \geq 1, \end{aligned}$$\displaystyle \begin{aligned} &\Psi_{2km}(r, \theta, \phi) = A_{km} \cdot \frac{\displaystyle 1}{\displaystyle \sqrt{r}} \cdot J_{k + 1/2}(\gamma \cdot r) \cdot P_k^{(m)}(\cos(\theta)) \cdot \cos(m \cdot \phi),\\ &0 \leq m \leq k, \quad k \geq 0, \end{aligned}$) в граничные условия (2.73$\displaystyle \frac{\displaystyle \partial \widehat{T}(r, \theta, \phi, t)}{\displaystyle \partial t}$ $\displaystyle = a^2 \cdot \frac{\displaystyle \partial^2 \widehat{T}(r, \theta, \phi, t)}{\displaystyle \partial r^2}$ $\displaystyle + \frac{2 \cdot a^2}{r} \cdot \frac{\displaystyle \partial \widehat{T}(r, \theta, \phi, t)}{\displaystyle \partial r}$ $\displaystyle + \frac{a^2}{r^2 \cdot \sin(\theta)} \cdot \frac{\displaystyle \partial}{\displaystyle \partial \theta} \left( \sin(\theta) \cdot \frac{\displaystyle \partial \widehat{T}(r, \theta, \phi, t)}{\displaystyle \partial \theta} \right)$ $\displaystyle + \frac{a^2}{r^2 \cdot \sin^2 \theta} \cdot \frac{\displaystyle \partial^2 \widehat{T}(r, \theta, \phi, t)}{\displaystyle \partial \phi^2}$ $\displaystyle + \widehat{f}(r, \theta, \phi, t),$$\displaystyle \widehat{T}(r, \theta, \phi, 0)$ $\displaystyle = T_0(r, \theta, \phi)$ $\displaystyle - U(\theta, \phi, 0),$$\displaystyle \widehat{T}(R, \theta, \phi, t)$ $\displaystyle = 0,$$\displaystyle |\widehat{T}(r, 0, \phi, t)|$ $\displaystyle < \infty, \quad |\widehat{T}(r, \pi, \phi, t)|$ $\displaystyle < \infty,$$\displaystyle \widehat{T}(r, \theta, \phi, t)$ $\displaystyle = \widehat{T}(r, \theta, \phi + 2\pi, t),$) и получим

Понятно, что решение имеет смысл только тогда, когда $A_{km}$ $\neq 0$, что возможно только тогда, когда $J_{k + 1/2}(\gamma \cdot R)$ $= 0$. Обозначим корни уравнения $J_{k + 1/2}(\mu)$ $= 0$ через $\mu_{ks}$ — тогда собственные значения примут вид

Нулевой и отрицательные корни не годятся: при $\mu$ $= 0$ собственная функция $J_{k + 1/2}(0 \cdot r)$ $\equiv 0$ — тождественный нуль, не являющийся собственной функцией (индекс $k$ $+ 1/2$ $> 0$), а отрицательные корни новых решений не дают. Поэтому берём только положительные корни: $\mu_{ks}, s$ $\in (1..\infty)$.

Таким образом, собственные значения и собственные функции имеют вид

где $A_{kms}$ — произвольный постоянный множитель: собственная функция определена с точностью до него. Положим $A_{kms}$ $= 1$ — этот множитель всё равно вбирается в коэффициенты разложения, как и в общем решении.

Для разложения функций в ряд Фурье по $\Psi_{1kms}(r, \theta, \phi)$ и $\Psi_{2kms}(r, \theta, \phi)$ необходимо вычислить нормы $\|\Psi_{1kms}(r, \theta, \phi)\|^2$ и $\|\Psi_{2kms}(r, \theta, \phi)\|^2$, вес для сферических координат $\rho(r, \theta, \phi)$ $= r^2 \cdot \sin(\theta)$.