Уравнение диффузии

Уравнение диффузии или уравнение теплопроводности представляет собой частный вид дифференциального уравнения в частных производных. Бывает нестационарным и стационарным.

Нестационарное уравнение

Нестационарное уравнение диффузии классифицируется как параболическое дифференциальное уравнение. Оно описывает распространение растворяемого вещества вследствие диффузии или перераспределение температуры тела в результате теплопроводности.

Одномерный случай

В случае одномерного диффузионного процесса с коэффициентом диффузии (теплопроводности) D уравнение имеет вид: $$\frac{\partial}{\partial t} c(x,t)=\frac{\partial}{\partial x} D \frac{\partial}{\partial x}{c(x,t)} + f(x,t).$$

При постоянном D приобретает вид: $$\frac{\partial}{\partial t} c(x,t)=D\frac{\partial^2}{\partial x^2} {c(x,t)} + f(x,t)$$

где c(x,t) – концентрация диффундирующего вещества, a f(x,t) - функция, описывающая источники вещества (тепла).

Трёхмерный случай

В трёхмерном случае уравнение приобретает вид: $$\frac{\partial}{\partial t} c(\vec{r},t) = (\nabla, D \nabla c(\vec{r},t)) + f (\vec{r},t) ,$$

а при постоянном D приобретает вид: $$\frac{\partial}{\partial t} c(\vec{r},t) = D \Delta c(\vec{r},t)+ f (\vec{r},t),$$

где $\Delta=\frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2}+\frac{\partial^2}{\partial z^2}$ — оператор Лапласа.

В одномерном случае фундаментальное решение однородного уравнения (при начальном условии, выражаемом дельта-функцией c(x,0)=δ(x) и граничном условии c(∞,t)=0) есть $$c(x,t)=\sqrt{\frac{1}{4\pi Dt}} e^{-\frac{x^2}{4Dt}}$$

В этом случае c(x,t) можно интерпретировать как плотность вероятности того, что одна частица, находившаяся в начальный момент времени в исходном пункте, через время t перейдёт в пункт с координатой x. Тогда средний квадрат пройденного пути есть $$= \int_{-\infty}^\infty x^2 c(x,t)dx = 2Dt$$

Такому параболическому закону подчиняется, например, рост окисных плёнок на поверхностях металлов при высоких температурах. Необходимым условием при этом является достаточное количество кислорода в окружающей среде.

В случае произвольного начального распределения $c(x,0)$ общее решение уравнения диффузии представляется в интегральном виде как свёртка: $$c(x,t) = \int_{-\infty}^\infty c(x',0) \frac{1}{ \sqrt{4\pi D t}} \exp\left(-\frac{(x-x')^2}{4 D t}\right)\,dx'$$

Частные замечания

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

Типичные ошибки

Наряду с собственными фантазиями авторов, существует априорный метод вывода уравнений для всех фундаментальных взаимодействий .

В частности, показано, что уравнения Максвелла изначально квантовые. Однако, они пропорциональны первой степени постоянной тонкой структуры $\alpha$, которая в системе единиц СИ может быть определена так: $$\begin{equation} \label{alpha} \boxed{\alpha=\frac{e^2}{\ 4 \pi \varepsilon_0 \hbar c}=\frac{e^2}{2 \varepsilon_0 h c}}\end{equation},$$ где $\ e$ — элементарный электрический заряд,

$\hbar=h/2\pi$ — постоянная Дирака (или приведённая постоянная Планка),

$\ c$ — скорость света в вакууме,

$\varepsilon_0$ — электрическая постоянная (раньше — диэлектрическая проницаемость вакуума).

Поэтому явно постоянная тонкой структуры $\alpha$ в уравнения Максвелла и не входит. А в пространственно-временное уравнение диффузии $\alpha$ входит.

Далее надо отметить, что априорная теория всего рассматривает два типа взаимодействий: пространственные и временные. Поэтому временное уравнение требует применения метода типа предиктор-корректор.

Стационарное уравнение

В случае, когда ставится задача по нахождению установившегося распределения плотности или температуры (например, в случае, когда распределение источников не зависит от времени), из нестационарного уравнения выбрасывают члены уравнения, связанные с временем. Тогда получается стационарное уравнение теплопроводности, относящееся к классу эллиптических уравнений. Его общий вид: $$-(\nabla, D \nabla c(\vec{r})) = f (\vec{r}) ,$$

Постановка краевых задач

Задача с начальными условиями (задача Коши) о распределении температуры на бесконечной прямой

Если рассматривать процесс теплопроводности в очень длинном стержне, то в течении небольшого промежутка времени влияние температур на границах практически отсутствует, и температура на рассматриваемом участке зависит лишь от начального распределения температур.

Найти решение уравнения теплопроводности в области $- \infty \le x \le \infty$ и $t \ge t_0$, удовлетворяющее условию $u(x,t_0 ) = \varphi (x){\rm{ }} ~( - \infty < x < + \infty )$, где $\varphi (x)$ — заданная функция.

Первая краевая задача для полубесконечного стержня

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

Найти решение уравнения теплопроводности в области $- \infty \le x \le \infty$ и $t \ge t_0$, удовлетворяющее условиям $\begin{cases} u(x,t_0 ) = \varphi (x){\rm{ }} \> (0 < x < \infty ) \\ u(0,t) = \mu (t){\rm{ }} \> (t \ge t_0 ) \\ \end{cases}$, где $\varphi (x)$ и $\mu (t)$ — заданнные функции.

Краевая задача без начальных условий

Если момент времени который нас интересует достаточно удален от начального, то имеет смысл принебречь начальными условиями, поскольку их влияние на процесс с течением времени ослабевает. Таким образом, мы приходим к задаче, в которой заданы краевые условия и отсутствуют начальные.

Найти решение уравнения теплопроводности в области $0 \le x \le l$ и $- \infty < t$, удовлетворяющее условиям $\begin{cases} u(0,t) = \mu _1 (t) \\ u(l,t) = \mu _2 (t) \\ \end{cases}$, где $\mu _1 (t)$ и $\mu _2 (t)$ — заданнные функции.

Краевые задачи для ограниченного стержня

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

$u_t = a^2 u_{xx} + f(x,t),~~0 < x < l,~~0 < t \le T$ — уравнение теплопроводности. Если $f(x,t) = 0$, то такое уравнение называют однородным, в противном случае — неоднородным.

$u(x,0) = \varphi (x),0 \le x \le l$ — начальное условие в момент времени $t = 0$, температура в точке $x$ задается функцией $\varphi (x)$

$\begin{cases} u(0,t) = \mu _1 (t),0 \le t \le T \\ u(l,t) = \mu _2 (t),0 \le t \le T \\ \end{cases}$ — краевые условия. Функции $\mu _1 (t)$ и $\mu _2 (t)$ задают значение температуры в граничных точка 0 и $l$ в любой момент времени $t$.

В зависимости от рода краевых условий, задачи для уравнения теплопроводности можно разбить на три типа. Рассмотрим общий случай при $\alpha _i^2 + \beta _i^2 \ne 0,~~(i = 1,2)$.

$\begin{cases} \alpha _1 u_x (0,t) + \beta _1 u(0,t) = \mu _1 (t) \\ \alpha _2 u_x (l,t) + \beta _2 u(l,t) = \mu _2 (t) \\ \end{cases}$

Если $\alpha _i=0,~~ (i=1,2)$, то такое условие называют условием первого рода, если $\beta _i = 0,~~(i = 1,2)$ — второго рода, а если $\alpha _i$ и $\beta _i$ отличны от нуля, то условием третьего рода. Отсюда получаем задачи для уравнения теплопроводности — первую, вторую и третью краевую.

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

Метод разделения переменных

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

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

$\begin{array}{l} ~~u_t = a^2 u_{xx} ~~0 < x < l,0 < t \le T \\ ~~u(x,0) = \varphi (x)~~0 \le x \le l \\ \left\{ u(0,t) = 0; u(l,t) = 0 \right\} ~~0 \le t \le T \\ \end{array}$
Требуется найти функцию $~u(x,t)$ для $\forall (x,t):~0 \le x \le l,0 \le t \le 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 = 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,~~n = 1,2...\\ \end{array}$$ $$\lambda _n = \left( {\frac{{\pi n}}{l}} \right)^2,~~n = 1,2...$$ Отсюда $$X_n (x) = C_n\sin \left( {\frac{{\pi n}}{l}x} \right),~~n = 1,2...$$

C учетом найденных $~\lambda$, выведем общее решение линейного дифференциального уравнения $$T '(t) + a^2 \left( {\frac{{\pi n}}{l}} \right)^2 T (t) = 0 .$$ Должен получиться ответ $$T_n (t) = D_n e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} ,~~D_n = const$$

Теперь все готово для того, чтобы записать решение исходной задачи: $$u_n (x,t) = X_n (x)T_n (x) = C_n\sin \left( {\frac{{\pi n}}{l}x} \right)e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} ,~~n = 1,2...$$ В результате у нас получилось бесконечное количество частных решений уравнения. Все эти частные решения линейно-независимы, т.е. линейная комбинация любого количества решений равна нулю, только если все коэффициенты при них равны нулю. Поэтому логично предположить, что суммируя все частные решения по 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)e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} }$$ Осталось определить значение константы 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( {\frac{{\pi n}}{l}x} \right)} \\ A_n = \frac{2}{l}\int\limits_0^l {\varphi (\xi )\sin \left( {\frac{{\pi n}}{l}\xi } \right)} d\xi \\ \end{array}$$ Получаем: $$\begin{array}{l} u(x,0) = \sum\limits_{n = 1}^\infty {C_n \sin \left( {\frac{{\pi n}}{l}x} \right)} = \sum\limits_{n = 1}^\infty {A_n \sin \left( {\frac{{\pi n}}{l}x} \right)} \\ C_n = A_n = \frac{2}{l}\int\limits_0^l {\varphi (\xi )\sin \left( {\frac{{\pi n}}{l}\xi } \right)} d\xi \\ \end{array}$$ Откуда общее решение: $$u(x,t) = \sum\limits_{n = 1}^\infty {\left( {\frac{2}{l}\int\limits_0^l {\varphi (\xi )\sin \left( {\frac{{\pi n}}{l}\xi } \right)} d\xi } \right)\sin \left( {\frac{{\pi n}}{l}x} \right)e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} }$$ В курсе математической физики доказывается, что полученный ряд удовлетворяет всем условиям данной задачи, т.е. функция $~u(x,t)$ дифференцируема (и ряд сходится равномерно), удовлетворяет уравнению в области определения и непрерывна в точках границы этой области.

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

Рассмотрим способ решения неоднородного уравнения: $$\begin{array}{l} ~ u_t = a^2 u_{xx} 0 + f(x,t),~~0< x < l,0 < t \le T \\ ~ u(x,0) = 0,~~0 \le x \le l \\ \left. \begin{array}{l} u(0,t) = 0 \\ u(l,t) = 0 \\ \end{array} \right\}~0 \le t \le 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( {\frac{{\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( {\frac{{\pi n}}{l}} \right)^2 X_n (x)T_n (t) + X_n (x)F_n (t) \\ T'_n (t) = - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 T_n (t) + F_n (t) \\ \end{array} .$$ Решим последнее линейное неоднородное уравнение методом вариации постоянной. Сначала найдем общее решение однородного линейного уравнения $$\begin{array}{l} T'_n (t) = - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 T_n (t) \\ T_n (t) = De^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} \\ \end{array}$$ В общем решении заменим постоянную D на переменную D(t) и подставим в исходное уравнение. $$\begin{array}{l} T_n (t) = D(t)e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} \\ D'_n (t)e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} D_n (t) = - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} D_n (t) + F_n (t) \\ D'_n (t)e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} = F_n (t) \\ D_n (t) = \int {F_n (t)} e^{a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} dt \\ T_n (t) = A_n e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} + e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} \int {F_n (t)} e^{a^2 \left( {\frac{{\pi n}}{l}} \right)^2 t} 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 {e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 \left( {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 {e^{ - a^2 \left( {\frac{{\pi n}}{l}} \right)^2 \left( {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 = \frac{{\mu _2 (t) - \mu _1 (t)}}{l} \\ U(x,t) = \frac{{\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) - \frac{{\mu _2 ^\prime (t) - \mu _1 ^\prime (t)}}{l}x - \mu _1 ^\prime (t) \\ \tilde u(x,0) = \varphi (x) - \frac{{\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$$