I. Функционал метода Бубнова–Галёркина

Запишем ещё раз функционал для уравнения Эйлера (5.4$\displaystyle I(\upsilon)$ $\displaystyle {} = (L[\upsilon], \upsilon)$​$\displaystyle {} - 2 \cdot (f, \upsilon).$)

Линейный оператор $L[\upsilon]$ для уравнения параболического типа, интересующего нас, имеет вид

Скалярное произведение $(L[\upsilon], \upsilon)$ имеет вид

Скалярное произведение $(f, \upsilon)$ имеет вид

Известна формула Остроградского для оператора Лапласа

где $S$ — граница подмногообразия $M$, а $dS$ — её векторный элемент, направленный по внешней нормали. Проинтегрируем по частям второй интеграл

Таким образом, исходный функционал можно записать в виде

Дадим приращение $\widehat{\upsilon}$​${} + \epsilon \cdot \upsilon$, где $\widehat{\upsilon}$ — точное решение, а $\epsilon$ — некоторое малое число. Полное приращение функционала можно записать в виде

Пренебрегая членами второго порядка малости, пропорциональными $\epsilon^2$, получим

Если вынести $\epsilon$ за скобки и принять во внимание, что $\epsilon$ может быть отрицательным, а вот приращение функционала всегда положительно, так как $\widehat{\upsilon}$ — точное решение, то получим в качестве обязательного условия

Учитывая, что $\widehat{\upsilon}$​${} \approx \upsilon$, и сокращая общий множитель, мы возвращаемся к исходному уравнению. Это означает, что численное решение краевой задачи сводится к минимизации функционала

Множитель $2$ при линейном члене $\int_M f \cdot \upsilon dM$ сохраняем сознательно: именно он обеспечивает, что условие минимума $\delta I$​${} = 0$ возвращает исходное уравнение, а не уравнение с лишним множителем $1/2$ (при дифференцировании квадратичные члены дают множитель $2$, и линейный член должен иметь такой же).

Перейдём к разделению на стационарную и нестационарную задачи. Принципиально важно проводить это разделение и последующую дискретизацию по времени на уровне уравнения, а не функционала: производная по времени не является самосопряжённым оператором, поэтому прямая подстановка $\upsilon$​${} = \psi(M) \cdot \phi(t)$ в функционал привела бы к неверным коэффициентам. Запишем слабую форму уравнения $L[u]$​${} = f$, полученную выше интегрированием по частям

которая должна выполняться для произвольной пробной функции $\upsilon$.

Стационарное уравнение с граничными условиями Дирихле и/или Неймана

Если $\frac{\displaystyle \partial u}{\displaystyle \partial t}$​${} = 0$, временной член исчезает, и в силу симметрии пространственного оператора слабая форма (I.4$\displaystyle \int_M \frac{\displaystyle \partial u}{\displaystyle \partial t} \cdot \upsilon \,dM$​$\displaystyle {} + a^2 \cdot \int_M \nabla u \cdot \nabla \upsilon \,dM$​$\displaystyle {} - a^2 \cdot \int_S \upsilon \cdot \nabla u \,dS$​$\displaystyle {} = \int_M f \cdot \upsilon \,dM,$) эквивалентна минимизации функционала

Нестационарное уравнение с граничными условиями Дирихле и/или Неймана

Производную по времени дискретизируем неявной схемой Эйлера с шагом $\Delta t$, применяя её к слабой форме уравнения: $\frac{\displaystyle \partial u}{\displaystyle \partial t}$​${} \approx \frac{\displaystyle u_n - u_{n-1}}{\displaystyle \Delta t}$, где $u_n$ — решение на текущем временном слое, а $u_{n-1}$ — на предыдущем. Подставив это в (I.4$\displaystyle \int_M \frac{\displaystyle \partial u}{\displaystyle \partial t} \cdot \upsilon \,dM$​$\displaystyle {} + a^2 \cdot \int_M \nabla u \cdot \nabla \upsilon \,dM$​$\displaystyle {} - a^2 \cdot \int_S \upsilon \cdot \nabla u \,dS$​$\displaystyle {} = \int_M f \cdot \upsilon \,dM,$) и домножив на $\Delta t$, получаем уравнение для шага по времени

Это уравнение, в свою очередь, является условием минимума функционала

Таким образом, нестационарное уравнение сведено к последовательному решению стационарных задач: на каждом шаге по известному $\upsilon_{n-1}$ находится $\upsilon_n$.

Матричная форма

Дискретизируем пробную функцию по узлам сетки: $\upsilon$​${} = \sum_i q_i \cdot \phi_i(M)$, где $\phi_i(M)$ — базисные функции, а $\overrightarrow{q}$​${} = (q_0, \ldots, q_n)^T$ — вектор узловых значений. Подстановка превращает каждый интеграл функционала в квадратичную или линейную форму по $\overrightarrow{q}$, а их коэффициенты собираются в матрицы:

где $K$ — матрица жёсткости (интеграл от произведения градиентов), $D$ — матрица демпфирования (интеграл от произведения пробных функций), $\overrightarrow{F}$ — вектор нагрузки (интеграл от источника); матрицы $K$ и $D$ симметричны. Граничный интеграл ${} -a^2 \int_S \upsilon \cdot \nabla \upsilon \,dS$ в матрицу не сворачивается: это поверхностный член, который в дискретной форме входит только в уравнения граничных узлов (у внутреннего узла пробная функция на границе равна нулю) и определяется самим граничным условием.

Дискретную систему получаем по методу Бубнова–Галёркина: подставляем $\upsilon$​${} = \sum_i q_i \phi_i$ в слабую форму (I.4$\displaystyle \int_M \frac{\displaystyle \partial u}{\displaystyle \partial t} \cdot \upsilon \,dM$​$\displaystyle {} + a^2 \cdot \int_M \nabla u \cdot \nabla \upsilon \,dM$​$\displaystyle {} - a^2 \cdot \int_S \upsilon \cdot \nabla u \,dS$​$\displaystyle {} = \int_M f \cdot \upsilon \,dM,$) и поочерёдно берём пробные функции $\phi_j$. Для стационарного уравнения это даёт систему

где $\overrightarrow{s}$ — вектор узловых граничных потоков. Сам поток здесь ничем не задаётся — что с ним делать, определяет тип граничного условия. Учёт граничных условий — задача непростая, и подробно она разобрана в отдельных приложениях: «Учёт граничных условий Дирихле» и «Учёт граничных условий Неймана».

Для нестационарного уравнения (неявная схема Эйлера, см. (I.6$\displaystyle \int_M \upsilon_n^2 \,dM$​$\displaystyle {} + \Delta t \cdot a^2 \cdot \int_M \nabla \upsilon_n \cdot \nabla \upsilon_n \,dM$​$\displaystyle {} - 2 \cdot \int_M \left[ \Delta t \cdot f + \upsilon_{n-1} \right] \cdot \upsilon_n \,dM$​$\displaystyle {} - \Delta t \cdot a^2 \cdot \int_S \upsilon_n \cdot \nabla \upsilon_n \,dS$​$\displaystyle {} \rightarrow \min.$)) система на шаге по времени принимает вид