Проверка результатов моделирования с помощью метода искусственного решения

Temesgen Kindo 27/07/2015
Share this on Facebook Share this on Twitter Share this on LinkedIn

Как мы можем удостовериться, правильно ли работает средство моделирования? Одним из подходов является метод искусственного решения. Этот процесс включает допущение некого решения, получение исходных условий и других вспомогательных условий согласующихся с первоначальным предположением, решение задачи с этими условиями в качестве входных параметров для пакета моделирования и сравнение полученных результатов с предполагаемым решением. Метод прост в использовании и достаточно универсален. Например, исследователи из Национальной лаборатории Сандия использовали его в некоторых собственных программных разработках.

Проверка и подтверждение

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

Для численного моделирования физической проблемы, мы предпринимаем два шага:

  1. Построение математической модели физической системы. Это то место, где учитываются все факторы, которые оказывают влияние на наблюдаемое поведение и формулируются основные управляющие уравнения. Результатом зачастую является система неявных связей между входными и выходными величинами. Они часто записываются в виде системы дифференциальных уравнений в частных производных с начальными и граничными условиями, которая в совокупности называется как начально-краевая задача (initial boundary value problem — IBVP).
  2. Решение математической модели для получения выражений для выходных величин в виде явных функций от входных значений параметров. Однако, такая замкнутая форма решения недоступна для большинства задач, представляющих практический интерес. В таком случае, мы применяем численные методы для получения приближенных решений, часто с помощью вычислительной техники для решения больших систем, в общем случае, нелинейных алгебраических уравнений и неравенств.

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

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

Сравнение между процессами подтверждения и проверки.

Теперь, давайте углубимся в проблему проверки численного решения начально-краевых задач (IBVP-задач).

Различные подходы к процедуре проверки

Как мы можем удостовериться, обеспечивает ли инструмент моделирования достаточно точное решение IBVP-задачи? Одна из возможностей заключается в том, чтобы выбрать задачу, имеющую точное аналитическое решение и использовать точное решение в качестве критерия правильности. Метод разделения переменных, например, можно использовать для получения решения простой IBVP-задачи. Применимость этого подхода ограничивается тем, что большинство задач, представляющих практический интерес, не имеют точного решения, что и является основанием для использования компьютерного моделирования. Тем не менее, этот подход является полезным при проверке работоспособности алгоритмов и программ.

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

Библиотеки Приложений в среде COMSOL Multiphysics содержат множество проверочных моделей (verification models), которые используют один или оба этих подхода. Они сгруппированы по областям физики.

Скриншот, показывающий каким образом получить доступ к моделям для верификации  среды COMSOL Multiphysics.
Данные модели доступны в Библиотеках Приложений среды COMSOL Multiphysics.

Что делать, если нам потребуется проверить наши результаты при отсутствии точных математических решений и экспериментальных данных? Мы можем обратиться к методу искусственного решения.

Реализация метода искусственного решения

Целью решения начально-краевой задачи (IBVP-задачи) является нахождения явного выражения для решения в терминах независимых переменных, как правило, пространства и времени, с учетом заданных параметров задачи, таких как свойства материала, а также граничные, начальные и исходные условия. Общая форма исходных условий подразумевает описание всех сил, действующих на тело, таких как гравитация в задачах строительной механики и аэро- и гидродинамики, силы сопротивления в задачах массопереноса и источников тепла в задачах термодинамики.

В методе искусственного решения (Method of Manufactured Solutions — MMS), мы решаем задачу от противного , т.е. начинаем с предполагаемого явного выражения для решения. Затем, мы подставляем решение в дифференциальные уравнения и получаем согласованный набор исходных, начальных и граничных условий. Обычно, это включает в себя оценку ряда производных. Вскоре мы увидим, как макрос на языке символической алгебры среды COMSOL Multiphysics может помочь в этом процессе. Аналогичным образом, мы вычисляем предполагаемое решение в момент времени t = 0 и на границах для получения начальных и граничных условий.

Далее идет этап проверки. С учетом только что полученных исходных и вспомогательных условий, мы применяем инструмент моделирования для получения численного решения IBVP-задачи и сравниваем его с изначально предполагавшимся решением, с которого мы “стартовали”.

Продемонстрируем эти шаги на простом примере.

Проверка одномерной задачи теплопроводности

Рассмотрим 1D задачу теплопроводности в стержне длиной L

A_c\rho C_p \frac{\partial T}{\partial t} + \frac{\partial}{\partial x}(-A_ck\frac{\partial T}{\partial x}) = A_cQ, t \in (0,t_f), x \in (0,L)

с начальным условием

T(x,0) = f(x)

и температурным режимом, поддерживаемым на обоих концах по закону, задаваемым функциями

T(0,t) = g_1(t), \quad T(L,t) = g_2(t).

Константами A_c, \rho, C_p и k являются площадь поперечного сечения, плотность материала, теплоемкость и коэффициент теплопроводности, соответственно. Источник тепла задается как Q.

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

u(x,t) = 500 K + \frac{x}{L}(\frac{x}{L}-1)\frac{t}{\tau} K,

где \tau есть характерное время, которое для данного примера выберем равным одному часу. Введем новую переменную u для предполагаемой температуры, чтобы отличить ее от расчетной температуры T.

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

В случае однородных материала и свойств поперечного сечения, мы можем объявить A_c, \rho, C_p и k как параметры. В общем неоднородном случае для этого потребуются переменные величины, как это проделано для зависящих от времени граничных условий. Обратите внимание на использование оператора d(), одного из встроенных операторов дифференцирования в среду COMSOL Multiphysics, как показано на скриншоте ниже.

Автоматизация оценки частных производных в моделировании.
Макрос на языке символической алгебры среды COMSOL Multiphysics может автоматизировать оценку частных производных.

Мы выполняем эти символические вычисления с той оговоркой, что мы доверяем символической алгебре. Иначе, любые ошибки наблюдаемые позднее проистекали бы из символических вычислений, а не численного решения. Разумеется, мы можем построить график вычисляемого вручную выражения для Q рядом с результатом символических вычислений, показанных выше, для проверки макроса.

Далее, мы вычисляем начальные и граничные условия. Начальное условие – это предполагаемое решение, вычисленное в момент времени t = 0.

f(x) = u(x,0) = 500 K.

Значения температуры на концах стержня есть g_1(t) = g_2(t) = 500 K.

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

Добавление начальных, граничных условий и внешних источников в дерева моделей среды COMSOL Multiphysics.
Добавление начальных значений, граничных условий и внешних источников, получаемых исходя из предполагаемого решения.

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

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

Решение, вычисленное при использовании искусственного решения с линейными элементами (слева) и квадратичными элементами (справа).

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

Метод MMS предоставляет нам необходимую гибкость для проверки различных частей кода. Ради простоты, в приведенном выше примере, мы умышленно оставили многие части IBVP-задачи непроверенными. На практике, каждый элемент в уравнении должен быть проверен в наиболее общем виде. Например, для проверки того, насколько точно программа обрабатывает области с неоднородным поперечным сечением, мы должны определить пространственно-переменную область до получения исходных условий. То же можно сказать и в отношении других коэффициентов уравнения, например для свойств материала.

Аналогичная проверка должна производиться для всех граничных и начальных условий. Если, к примеру, мы хотим на левом конце стержня задать поток тепла взамен температуры, мы сначала оцениваем поток, соответствующий искусственному решению, т. е. -n\cdot(-A_ck \nabla u) где n есть единичная внешняя нормаль. Для предполагаемого решения в данном примере, направленный внутрь поток на левом конце становится равным \frac{A_ck}{L}\frac{t}{\tau}*1K.

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

u(x,t) = 500 K + (\frac{x}{L})^2\frac{t}{\tau} K.

Обратите внимание, что во время проверки, мы проверяем только, корректно ли решаются уравнения. Мы не затрагиваем вопроса о том, соответствует ли решение физической ситуации.

Важно помнить, что как только мы выдумываем новое искусственное решение, мы должны пересчитать исходные, начальные и граничные условия в соответствии предполагаемым решением. Разумеется, при использовании средств символических вычислений в среде COMSOL Multiphysics, мы освобождены от утомительных вычислений!

Скорость сходимости

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

Например, для стационарной модификаци изадачи, стандартная оценка конечного- элементной погрешности для погрешности измеряемой в норме Соболева

m-порядка, есть

||u-u_h||_m \leq C h^{p+1-m}||u||_{p+1},

где u и u_h являются точным и конечно-элементным решениями, h является максимальным размером элемента, p -порядок полиномиальной аппроксимации (функция формы). Для m = 0, это дает оценку погрешности

||e|| = ||u-u_h|| = (\int_{\Omega}(u-u_h)^2d\Omega)^{\frac{1}{2}} \leq C h^{p+1}||u||_{p+1},

где C – константа независимая от сетки разбиения.

Возвращаясь к методу искусственного решения, это означает, что решение с линейным элементом (p = 1) должно показывать сходимость второго порядка при измельчении сетки. Если мы построим график нормы ошибки по отношению к размеру сетки в логарифмическом масштабе, тангенс угла наклона должен асимптотически приближаться к 2. Если этого не происходит, то мы должны будем проверить код или точность и регулярность входных параметров, таких как материальные и геометрические свойства. Как показано на рисунках ниже, численное решение сходится с теоретически ожидаемой скоростью.

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

Слева: Использование операторов интегрирования для определения норм. Оператор intop1 определяется для интегрирования по области. Справа: График ошибки в зависимости от размера сетки, построенный в логарифмических осях, показывает второй порядок сходимости в норме L_2 (m = 0) для линейных элементов, что согласуется с теоретической оценкой.

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

Нелинейный задачи и взаимосвязанные междисциплинарные задачи

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

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

Единственность

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

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

Попробуйте Самостоятельно

Встроенная функциональность в виде символических вычислений в среде COMSOL Multiphysics позволяет легко реализовать метод MMS для проверки кода программ, а также для образовательных целей. Несмотря на то, что мы всесторонне тестируем коды наших программ, мы приветствуем любые замечания со стороны наших пользователей. В данном топике был представлен универсальный инструмент, который вы можете использовать для проверки различных физических интерфейсов. Вы также можете проверить свои собственные реализации при использовании моделирования на основе пользовательских уравнений или Построителя интерфейсов для новой физики (Physics Builder) в среде COMSOL Multiphysics. Если у вас есть какие-либо вопросы по данной методике , пожалуйста, не стесняйтесь написать нам!

Информационные ресурсы

  • Для всестороннего обсуждения метода искусственного решения, включая его сильные и слабые стороны, посмотрите данный Отчет из Национальной лаборатории Сэндиа (Sandia). В Отчете подробно излагается система “слепых” тестов, в которых один автор заложил ряд ошибок в код программы незаметно от второго автора, который выискивает заложенные “мины” с использованием метода, описанного в этом топике.
  • Для расширенной дискуссии по проверке и подтверждению в контексте научных вычислений, ознакомьтесь с
    • W. J. Oberkampf and C. J. Roy, Проверка и подтверждение в научных расчетах (Verification and Validation in Scientific Computing), Cambridge University Press, 2010
  • Оценка стандартной погрешности в методе конечных элементов имеется в таких книгах, как
    • Thomas J. R. Hughes, Метод конечных элементов: Линейный статический и динамический конечно-элементный анализ (The Finite Element Method: Linear Static and Dynamic Finite Element Analysis), Dover Publications, 2000
    • B. Daya Reddy, Вводный функциональный анализ: с приложениями к краевым задачам и конечным элементам (Introductory Functional Analysis: With Applications to Boundary Value Problems and Finite Elements), Springer-Verlag, 1997

Загрузка комментариев...

Темы публикаций


Теги

3D печать Cерия "Гибридное моделирование" Введение в среду разработки приложений Видео Волновые электромагнитные процессы Глазами пользователя Графен Интернет вещей Кластеры Моделирование высокочастотных электромагнитных явлений на различных пространственных масштабах Модуль AC/DC Модуль MEMS Модуль Акустика Модуль Волновая оптика Модуль Вычислительная гидродинамика Модуль Геометрическая оптика Модуль Динамика многих тел Модуль Композитные материалы Модуль Коррозия Модуль Механика конструкций Модуль Миксер Модуль Нелинейные конструкционные материалы Модуль Оптимизация Модуль Плазма Модуль Полупроводники Модуль Радиочастоты Модуль Роторная динамика Модуль Теплопередача Модуль Течение в трубопроводах Модуль Химические реакции Модуль Электрохимия Модуль аккумуляторов и топливных элементов Охлаждение испарением Пищевые технологии Рубрика Решатели Серия "Геотермальная энергия" Серия "Конструкционные материалы" Серия "Электрические машины" Серия “Моделирование зубчатых передач” Сертифицированные консультанты Технический контент Указания по применению физика спорта