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

29/01/2014

Интегрирование — один из важнейших математических инструментов, особенно в численном моделировании. Например, дифференциальные уравнения в частных производных обычно выводятся из интегральных уравнений сохранения. Когда возникает необходимость численного решения уравнения в частных производных, интегрирование также играет важную роль. В этой статье приведен обзор методов и подходов интегрирования, доступных в COMSOL Multiphysics, а также конкретные примеры их использования.

Важность интегралов

В COMSOL используется метод конечных элементов, который преобразует описывающее некоторый процесс уравнение в частных производных в интегральное уравнение — другими словами, в слабую форму (weak form). При детальном и глубоком изучении формулировок, используемых в интерфейсах COMSOL, вы обнаружите, что множество граничных условий реализованы через интегралы. В качестве наиболее характерных примеров можно привести условия Total heat flux (Общий тепловой поток) или Floating potential (Плавающий потенциал). Вычисление интегралов также играет ключевую роль в процессе постобработки результатов, поскольку COMSOL рассчитывает большое количество вспомогательных величин через интегралы, например энергию электрического поля, скорость потока или общий тепловой поток. Разумеется, пользователи вольны использовать интегрирование в COMSOL в своих целях, и в этой статье мы объясним вам, как это делать максимально эффективно.

Вычисление интегралов в узле Derived Values

Интеграл общего вида имеет форму

\int_{t_0}^{t_1}\int_{\Omega}F(u)\ \mathrm{d A} \mathrm{d} t

где [t_0,t_1] — это временной интервал, \Omega — это пространственная область, а F(u) — это произвольное выражение, включающее зависимую переменную u и произвольные функции от нее, в том числе производные по пространству, времени, а также любой другой величине.

Наиболее удобный способ вычисления интегралов — использование узла Derived Values (Расчет выражений) в разделе Results (Результаты) ленты Ribbon или дерева модели (Лента Ribbon отсутствует в том случае, если ваш компьютер работает не под управлением ОС Windows®).
Как добавить операцию вычисления выражения
Добавление операций расчета пространственных интегралов по объему, поверхности или линии в узле Derived Values (Расчет выражений)

Вы можете обратиться к любому доступному решению, выбрав соответствующий набор данных (data set). В поле Expression (Выражение) вводится подынтегральная функция, включающая зависимые или производные переменные. Для данных расчета во временной области пространственный интеграл вычисляется на каждом временном шаге. В качестве альтернативы, в окне Settings (Настройки) узла Data Series Operations (Операции с массивами данных) можно выбрать опцию Integral (Интегрирование), что позволит вычислить общий пространственно-временной интеграл.
Окно настроек для вычисления интеграла по поверхности
Пример настроек вычисления интегралов по поверхности (Surface Integration) с дополнительным вычислением интеграла по времени в разделе Data Series Operations.

Оператор Average (Усреднение) — еще одна операция в разделе Derived Values, связанная с вычислением интегралов. Оператор вычисляет интеграл и делит его на объем, площадь или длину выбранной области. Операция Averageв узле Data Series Operations аналогично вводит деление на продолжительность временного диапазона. Операторы узла Derived Values — важный инструмент, однако их можно использовать только во время постобработки, а значит с их помощью можно рассчитать далеко не любой интеграл. Именно поэтому в COMSOL представлены другие более мощные и гибкие инструменты для вычисления интегралов. Мы продемонстрируем их работу на представленном ниже модельном примере.

Вычисление пространственного и временного интегралов для демо-модели из области теплопередачи

Рассмотрим простую модель теплопередачи: двумерный единичный квадрат из алюминия в (x,y)-плоскости. Температура верхней и правой сторон постоянна и равна комнатной (293.15 K), в то время как для левой и нижней границы задан общий входящий тепловой поток (General inward heat flux), составляющий 5 000 W/m^2. Стационарное и нестационарное решение (в момент времени 100 секунд) представлены на иллюстрациях ниже.

Стационарное решение
Стационарное решение, нажмите на изображение для увеличения.
Нестационарное решение
Нестационарное решение (для момента времени 100 секунд), нажмите на изображение для увеличения.

Вычисление пространственного интеграла с использованием операторов узла Component Coupling

Операторы узла Component Coupling (Сопряжение компонентов) используются в тех случаях, когда, например, в одном выражении объединяются несколько интегралов, или интегралы требуются в процессе вычислений, или требуется множество контурных интегралов. Операторы данного узла определяются в разделе Definitions (Определения). На этом этапе режультат использования оператора не просчитывается, а указываются только их название и выборки областей.

Добавление операторов из узла Component Coupling
Добавление операторов через узел Component Couplings

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

\int_{\Omega}T(x,y)\ \mathrm{d}x\mathrm{d}y = 301.65

В пакете COMSOL оператор вычисления интеграла по умолчанию получает имя intop1.

Окно настроек оператора интегрирования
Окно настроек оператора интегрирования.
Расчет результата интегрирования через оператор
Расчет результата интегрирования через оператор.

Теперь давайте рассмотрим, как оператор интегрирования может использоваться непосредственно в процессе расчета модели. С его помощью мы могли бы, например, выяснить, какая нагревательная мощность потребуется для получения средней температуры 303.15 К, то есть температуры, на 10 К превышающей температуру окружающей среды. Прежде всего нам необходимо вычислить разницу между требуемым и действительным средними значениями. Среднее значение вычисляется путем деления интеграла от T на интеграл от постоянной функции 1, который равен площади области. Нетрудно догадаться, что вычисление подобного вида легко выполнить с помощью представленного в COMSOL оператора Average (Усреднение), см. комментарии выше. По умолчанию данный оператор получает название aveop1. Обратите внимание, что среднее значение для области в нашем примере совпадает с интегралом, т.к. область имеет единичную площадь. Соответствующая разность равна

303.15-\int_{\Omega}T(x,y)\mathrm{d} x\mathrm{d} y = 1.50

Далее нам необходимо найти значение General heat flux (Общий тепловой поток) на левой и нижней границах, при котором была бы достигнута требуемая средняя температура. Для этого мы введем дополнительную степень свободы под названием q_hot и дополнительное ограничение в качестве глобального уравнения (Global Equation). Общий входящий тепловой поток (General inward heat flux) тогда перезадается через переменную q_hot.

Добавление дополнительной степени свободы и глобального уравнения
Добавление дополнительной степени свободы (переменной) и глобального уравнения, для неяного подбора средней температуры, равной 303.15 K.

Решение данной сопряженной системы с помощью стационарного исследования дает значение q_{hot}=5 881.30 W/m^2. Т.е. полученное значение можно задать в качестве граничного условия для общего входящего теплового потока, чтобы средняя температура в рассматриваемой области стала равна 303.15 К.

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

В своих обращениях в службу поддержки пользователи часто задают один и тот же вопрос: как рассчитать неопределенный пространственный интеграл? Для этой цели нам также пригодится оператор интегрирования, задаваемый через Component Couplings. Нахождение неопределенного интеграла — операция, обратная дифференцированию. Неопределенный интеграл позволяет вычислять площади произвольных областей, ограниченных графиками функций. Одна из самых важных прикладных задач — вычисление вероятностей в статистическом анализе. Для того чтобы это продемонстрировать, мы зафиксируем y=0 и обозначим неопределенный интеграл от T(x,0) как u(x). Это значит, что \frac{\partial u}{\partial x}=T(x,0). Тогда неопределенный интеграл имеет вид

u(\bar x) = \int_0^{\bar x}T(x,0)\mathrm{d} x

Здесь мы используем \bar x, чтобы отличать переменную интегрирования от внешней переменной. В отличие от представленных выше интегралов, результатом интегрирования является функция, а не скалярная величина. Нам необходимо указать для ПО, что для каждого значения \bar x\in[0,1] соответствующее значение u(\bar x) вычисляется при помощи интеграла. В среде COMSOL это можно легко сделать всего за три шага. Во-первых, потребуется логическое выражение для переписывания интеграла в виде

u(\bar x) = \int_0^1T(x,0)\cdot(x\leq\bar x)\ \mathrm{d} x

Во-вторых, нам понадобится оператор вычисления интеграла, который будет действовать на нижней границе области из примера. Давайте обозначим его как intop2. В-третьих, мы должны отличать переменную интегрирования от внешней переменной. Принятые обозначения для такого случая: x называется источником (source), а \bar xточкой назначения (destination). При использовании операторов интегрирования доступен встроенный оператор dest, который позволяет явно оглашать, что соответствующее выражение не относится к переменным интегрирования. Точнее, это значит, что в COMSOL \bar x=dest(x). Объединив логическое выражение с оператором dest, мы получим выражение вида T*(x<=dest(x)), которое является именно тем входным выражением, которое требуется для intop2. Объединив все вместе, мы можем вычислить неопределенный интеграл, воспользовавшись выражением intop2(T*(x<=dest(x))). Результат данной операции можно проиллюстрировать следующим графиком:

График неопределенного интеграла
Как построить график неопределенного интеграла с помощью оператора интегрирования, оператора dest и логического выражения.

В пакете COMSOL дополнительно доступны еще два оператора вычисления интеграла, а именно общая проекция (general projection) и линейная проекция (linear projection). Эти операторы можно использовать для получения множества контурных интегралов в любом направлении в области. Другими словами, вычисление интеграла производится только вдоль одного измерения. В результате мы получаем функцию размерности на единицу меньше, чем размерность области. Для двухмерного примера результатом будет одномерная функция, которая может быть рассчитана на любой границе. Более подробная информация об использовании данных операторов будет представлена в одной из следующих публикаций в нашем компоративном блоге.

Вычисление пространственного интеграла посредством дополнительного физического интерфейса

Наиболее гибким способом вычисления пространственных интегралов является техника с добавлением дополнительного PDE-интерфейса. Давайте вспомним пример выше с неопределенным интегралом и предположим, что мы хотим вычислить неопределенный интеграл не только для y=0. Данная задача может быть сформулирована в виде дифференциального уравнения в частных производных

\frac{\partial u}{\partial x}=T(x,y)

с граничным условием типа Дирихле u=0 на левой границе. Расчет такого уравнения проще всего реализовать в физическом (математическом) интерфейсе Coefficient Form PDE (Дифференциальное уравнение в частных производных, коэффициентная форма записи), который потребует следующих настроек:

PDE-интерфейс для вычисления пространственного интеграла
Вычисление пространственного интеграла посредством дополнительного PDE-интерфейса.

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

Вычисление временного интеграла посредством встроенных операторов

Мы уже упоминали узел Data Series Operations, который можно использовать для вычисления интеграла по времени. Другой крайне полезный способ вычисления интеграла по времени обеспечивается встроенными операторами timeint и timeavg для вычисления интеграла или среднего значения по времени, соответственно. Данные инструменты доступны при постобработке результатов, используются для вычисления интеграла от любого выражения (зависящего от времени) на заданном временном интервале. В нашем примере мы рассчитаем среднее значение температуры в диапазоне от 90 до 100 секунд, то есть:

\frac{1}{10}\int_{90}^{100}T(x,y,t)\ \mathrm{d} t

На поверхностном графике ниже представлен результирующий интеграл, являющийся функцией пространственных переменных (x,y):

Использование timeavg - оператора вычисления интеграла по времени
Использование оператора timeavg – оператора вычисления интеграла по времени.

Схожие операторы существуют для вычисления интегралов на сферических зонах, а именно ballint, circint, diskint и sphint.

Вычисление временного интеграла посредством дополнительного физического интерфейса

В случае если временные интегралы нужно использовать непосредственно в модели в процессе расчета, вам будет необходимо задать их как дополнительные зависимые переменные. Аналогично представленному выше примеру с интерфейсом Coefficient Form PDE, это можно сделать, добавив ODE-интерфейс из раздела Mathematics. Предположим, например, что на каждом временном шаге требуется вычислять интеграл от величины общего теплового потока на промежутке от старта до текущего момента, который показывает накопленную энергию. Переменная для общего теплового потока рассчитывается в COMSOL автоматически и называется ht.tfluxMag. Интеграл может быть вычислен как дополнительная зависимая переменная с помощью узла Distributed ODE (Распределенное обыкновенное дифференциальное уравнение) интерфейса Domain ODEs and DAEs. Правой частью (источниковым членом) для доменного ОДЕ должна выступать подынтегральная функция, что и показано на иллюстрации ниже.

Вычисление временного интеграла через ODE
Использование дополнительного ODE-интерфейса для вычисления интеграла по времени.

В чем польза подобной техники? Полученный интеграл можно повторно использовать в других физических интерфейсах, поля в которых могут зависеть от накопленной в системе энергии. Более того, полученный резултат будет мгноменно доступен для всех видов постобработки, что удобнее и быстрее, чем использование встроенных операторов. Рекомендуем ознакомится с моделью Carbon Deposition in Hetereogeneous Catalysis (Образование сажевых отложений при гетерогенном катализе), в которой ОДЕ в области используется для вычисления пористости катализатора при наличии химических реакций в виде нестационарной полевой переменной.

Вычисление интеграла от аналитических функций и выражений

До сих пор мы демонстрировали, каким образом вычислять интеграл от искомых переменных в процессе расчета или при постобработке. Но не касались случая взятия интегралов от аналитических функций или выражений. Для этой операции в среде COMSOL доступен встроенный оператор integrate(expression, integration variable, lower bound, upper bound).

Выражение может представлять собой любую одномерную функцию, например sin(x). При этом допускается включение дополнительных переменных, например sin(x*y). Второй параметр определяет, по какой переменной вычисляется интеграл. Например, integrate(sin(x*y),y,0,1) выдает функцию переменной x, потому что интегрирование выполняется только по переменной y. Обратите внимание, что данный оператор также может использоваться для работы с аналитическими функциями, заданными в узле Definitions (Определения) текущего компонента.

Добавление аналитической функции
Добавление аналитической функции.
Вычисление интеграла от аналитической функции
Вычисление интеграла от аналитической функции.

Материалы для дальнейшего изучения

Рубрики блога


Комментарии (0)

Оставить комментарий
Войти | Регистрация
Загрузка...
РУБРИКАТОР БЛОГА COMSOL
РУБРИКИ
ТЕГИ