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

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

Предварительные понятия

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

От устройства к числу

Инженер редко решает уравнения Лапласа ради них самих. Ему нужно знать ёмкость кабеля, силу притяжения якоря, температуру обмотки. Между реальным устройством и этими числами лежит цепочка решений, каждое из которых вносит свою погрешность и свои допущения. Численное моделирование — это не «нажать кнопку и получить картинку», а осознанное прохождение этой цепочки с пониманием, что именно теряется на каждом шаге.

Цепочка выглядит так:

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

Ошибки, сделанные на шагах 2–3, никакое сгущение сетки не исправит: если выбрана неверная модель, её точное решение так же бесполезно, как и грубое. Поэтому начнём с модельных решений, а уже затем перейдём к математике дискретизации.

Концептуальная и математическая модель

Расчётная область и её усечение

Поле, строго говоря, занимает всё пространство. Расчётная область конечна, и её внешнюю границу приходится придумывать. Для электростатической задачи о двух электродах окружающий «воздушный ящик» обычно берут в несколько характерных размеров системы и задают на нём либо $\varphi=0$ (условие Дирихле), либо $\partial\varphi/\partial n=0$ (условие Неймана). Оба условия искусственны. Практический приём: решить задачу с обоими типами условий на внешней границе и увеличивать область до тех пор, пока результаты не перестанут заметно различаться. Если расхождение велико, область мала, а не сетка груба.

Упрощение геометрии

Чертёж изделия содержит детали, несущественные для поля, но дорогие для сетки: мелкие скругления, тонкие покрытия, резьбу. Решение об их удалении зависит от искомой величины. Скругление радиусом 0,1 мм не изменит ёмкость конденсатора с зазором 10 мм, но определит максимальную напряжённость на кромке электрода — для расчёта пробивного напряжения его убирать нельзя. Тонкий слой лака толщиной в десятки микрометров сетка разрешит только ценой огромного числа элементов; его заменяют поверхностным условием (контактным сопротивлением в тепловой задаче, поверхностным импедансом — в задаче о вихревых токах). Симметрия сокращает область в два, четыре и более раз: на плоскости симметрии ставится условие Неймана или Дирихле в зависимости от чётности поля относительно этой плоскости.

Материальные данные

Проницаемость, проводимость, теплопроводность берут из справочников, протоколов материалов или измерений. Все они известны с ограниченной точностью: проводимость меди меняется на 0,4 % на каждый кельвин, кривая $B$ – $H$ электротехнической стали различается от партии к партии на десятки процентов в области насыщения, а коэффициент теплоотдачи $\alpha$ при естественной конвекции известен в лучшем случае с точностью 20–30 %. Погрешность исходных данных часто превышает погрешность численного метода, и бессмысленно добиваться четвёртого знака от сетки, когда третий не гарантирован материалом.

Граничные условия как модельные решения

Граничное условие — не свойство уравнения, а утверждение о том, что происходит за пределами области. Условие $\varphi=U$ на электроде означает, что электрод идеально проводящий; условие $\partial T/\partial n=0$ на поверхности — что она идеально теплоизолирована; конвективное условие $-\lambda\,\partial T/\partial n=\alpha(T-T_\infty)$ заменяет всю гидродинамику окружающего воздуха одним коэффициентом. Каждое из них — упрощение внешнего мира, и его следует проверять так же критически, как геометрию и материалы.

Метод конечных элементов: идея

Пусть требуется решить уравнение $-\nabla\cdot (\varepsilon\,\nabla \varphi)=\rho$ в области $\Omega$ с границей $\Gamma$ . Точное решение — функция бесконечного числа «степеней свободы». Метод конечных элементов (МКЭ) заменяет её функцией, определяемой конечным набором чисел — значениями в узлах сетки.

Аппроксимация и функции формы

Область $\Omega$ разбивается на элементы (на плоскости — треугольники или четырёхугольники) с узлами $\ur_i$ , $i=1,\dots,n$ . Приближённое решение ищется в виде $\varphi(\ur)\approx\varphi_h(\ur)=\sum_{i=1}^{n}u_i\,N_i(\ur),$ где $u_i$ — неизвестные узловые значения, а $N_i(\ur)$ — функции формы: заданные кусочно-полиномиальные функции, отличные от нуля только в элементах, примыкающих к узлу $i$ . Индекс $h$ обозначает характерный размер элемента.

Для линейного треугольника с вершинами $(x_1,y_1)$ , $(x_2,y_2)$ , $(x_3,y_3)$ и площадью $S$ функции формы линейны по координатам: $N_i=\frac{a_i+b_i x+c_i y}{2S},\qquad a_1=x_2y_3-x_3y_2,\quad b_1=y_2-y_3,\quad c_1=x_3-x_2,$ остальные коэффициенты получаются циклической перестановкой индексов $1\to2\to3\to1$ . Удвоенная площадь $2S=(x_2-x_1)(y_3-y_1)-(x_3-x_1)(y_2-y_1)$ положительна при обходе вершин против часовой стрелки.

Свойства функций формы линейного треугольника

Свойство Кронекера: $N_i(\ur_j)=\delta_{ij}$ — функция равна единице в своём узле и нулю в остальных. Поэтому коэффициент $u_i$ есть именно значение $\varphi_h$ в узле $i$ .
Разбиение единицы: $N_1+N_2+N_3=1$ в каждой точке элемента. Следовательно, постоянное поле $\varphi=\mathrm{const}$ воспроизводится точно.
Градиент постоянен в элементе: $\nabla N_i=(b_i,\ c_i)/(2S)$ , а потому $\nabla \varphi_h=\sum_i u_i\nabla N_i$ — один и тот же вектор во всех точках треугольника.

Из третьего свойства следует важное: напряжённость $\uE_h=-\nabla \varphi_h$ в линейных элементах кусочно-постоянна и терпит разрыв при переходе через рёбра. Сам потенциал $\varphi_h$ непрерывен, поскольку на общем ребре двух треугольников он определяется одними и теми же двумя узловыми значениями.

Слабая формулировка

Приближение $\varphi_h$ — ломаная поверхность; её вторые производные не существуют, и подставить её в уравнение Пуассона напрямую невозможно. Выход — заменить уравнение равносильным ему интегральным соотношением, в котором производных на порядок меньше.

Слабая (интегральная) форма уравнения Пуассона

Умножим уравнение $-\nabla\cdot (\varepsilon\nabla \varphi)=\rho$ на произвольную достаточно гладкую пробную функцию $v(\ur)$ и проинтегрируем по области: $-\int_\Omega v\,\nabla\cdot (\varepsilon\nabla \varphi)\,\dV=\int_\Omega\rho\,v\,\dV .$ Воспользуемся тождеством $\nabla\cdot (v\,\varepsilon\nabla \varphi)=v\,\nabla\cdot (\varepsilon\nabla \varphi)+\varepsilon\,\nabla v\cdot\nabla \varphi$ и теоремой Остроградского–Гаусса из урока об интегральных теоремах: $\int_\Omega\varepsilon\,\nabla \varphi\cdot\nabla v\,\dV=\int_\Omega\rho\,v\,\dV+\oint_\Gamma v\,\varepsilon\frac{\partial\varphi}{\partial n}\,\dS .$ Это и есть слабая форма. В ней от $\varphi$ требуется лишь существование первых производных — ломаная поверхность $\varphi_h$ подходит.

Поверхностный интеграл раскрывает природу граничных условий. Пусть граница состоит из части $\Gamma_D$ , где потенциал задан ( $\varphi=\varphi_0$ ), и части $\Gamma_N$ , где задана нормальная составляющая смещения ( $-\varepsilon\,\partial\varphi/\partial n=D_n$ ).

На $\Gamma_D$ значение решения известно заранее, поэтому пробные функции выбирают равными нулю на $\Gamma_D$ , и соответствующая часть интеграла исчезает. Такое условие называют главным (существенным): его накладывают явно на узловые значения.
На $\Gamma_N$ в интеграл подставляется известная величина $D_n$ . Такое условие называют естественным: оно входит в правую часть само, без каких-либо ограничений на неизвестные. Однородное условие Неймана $\partial\varphi/\partial n=0$ вообще не требует никаких действий — оно выполняется «по умолчанию» на любой границе, где ничего не задано.

Это объясняет известный практический факт: свободная граница сетки ведёт себя как линия симметрии поля (нормальная составляющая $\uD$ равна нулю).

Метод Галёркина

В слабую форму подставляется приближение $\varphi_h=\sum_j u_jN_j$ , а в качестве пробных функций поочерёдно берутся сами функции формы $v=N_i$ . Получается система $n$ линейных уравнений относительно $n$ неизвестных $u_j$ : $\sum_{j=1}^{n}K_{ij}u_j=f_i,\qquad K_{ij}=\int_\Omega\varepsilon\,\nabla N_i\cdot\nabla N_j\,\dV,\qquad f_i=\int_\Omega\rho\,N_i\,\dV+\int_{\Gamma_N}(-D_n)N_i\,\dS .$ Матрицу $\mathbf K$ по традиции, идущей из механики, называют матрицей жёсткости, вектор $\mathbf f$ — вектором нагрузки, систему записывают как $\mathbf K\mathbf u=\mathbf f$ [Сегерлинд; Зенкевич].

Метод Галёркина для уравнения Пуассона равносилен минимизации энергетического функционала $\tfrac12\int\varepsilon|\nabla \varphi_h|^2\dV-\int\rho\varphi_h\dV$ на множестве функций вида $\sum u_jN_j$ [Сильвестер]. Иными словами, МКЭ находит среди всех ломаных поверхностей ту, которая ближе всех к точному решению по энергии. Этот факт будет ключом к пониманию скоростей сходимости в уроке о сходимости.

Элементная матрица и сборка

Интеграл $K_{ij}$ распадается на сумму интегралов по элементам, причём в элемент $e$ вклад дают только три его узла. Поскольку градиенты функций формы постоянны, интегрирование сводится к умножению на площадь: $K^{e}_{ij}=\int_{e}\varepsilon\,\nabla N_i\cdot\nabla N_j\,\dS =\varepsilon\,S\,\frac{b_ib_j+c_ic_j}{(2S)^2} =\frac{\varepsilon\,(b_ib_j+c_ic_j)}{4S},\qquad i,j=1,2,3 .$ Вектор нагрузки от постоянной в элементе плотности заряда распределяется между вершинами поровну: $f^e_i=\int_e\rho N_i\,\dS=\rho S/3$ , так как интеграл линейной функции по треугольнику равен произведению площади на её среднее значение, а среднее значение $N_i$ по треугольнику равно $1/3$ .

Матрица жёсткости одного треугольника

Треугольник с вершинами $1\,(0,0)$ , $2\,(1,0)$ , $3\,(0,1)$ в среде с $\varepsilon=1$ (единицы условные). Площадь $S=\tfrac12$ , $2S=1$ , $4S=2$ .

Коэффициенты: $b_1=y_2-y_3=-1$ , $c_1=x_3-x_2=-1$ ; $b_2=y_3-y_1=1$ , $c_2=x_1-x_3=0$ ; $b_3=y_1-y_2=0$ , $c_3=x_2-x_1=1$ . Так как $2S=1$ , функции формы $N_1=1-x-y$ , $N_2=x$ , $N_3=y$ (при этом $a_1=1$ , $a_2=a_3=0$ ). Проверка свойства Кронекера: $N_1(0,0)=1$ , $N_1(1,0)=0$ , $N_1(0,1)=0$ ; сумма $N_1+N_2+N_3=1$ .

Элементы матрицы по формуле $K_{ij}=(b_ib_j+c_ic_j)/2$ : $K_{11}=\frac{1+1}{2}=1,\quad K_{12}=\frac{-1+0}{2}=-\tfrac12,\quad K_{13}=\frac{0-1}{2}=-\tfrac12,\quad K_{22}=\frac{1+0}{2}=\tfrac12,\quad K_{23}=\frac{0+0}{2}=0,\quad K_{33}=\frac{0+1}{2}=\tfrac12 .$ $\mathbf K^e=\begin{pmatrix}1&-\tfrac12&-\tfrac12\\-\tfrac12&\tfrac12&0\\-\tfrac12&0&\tfrac12\end{pmatrix}.$ Матрица симметрична, сумма элементов каждой строки равна нулю: при одинаковых узловых потенциалах ( $\mathbf u=(1,1,1)$ ) получается $\mathbf K^e\mathbf u=0$ — постоянный потенциал не создаёт поля. Диагональные элементы положительны. Энергия $\tfrac12\mathbf u^{\mathsf T}\mathbf K^e\mathbf u$ при $\mathbf u=(0,1,0)$ равна $\tfrac14$ : действительно, $\nabla \varphi_h=(1,0)$ , $|\nabla \varphi_h|^2=1$ , и $\tfrac12\cdot1\cdot S=\tfrac14$ .

Сборка. Глобальная матрица получается суммированием элементных: вклад $K^e_{ij}$ прибавляется к элементу глобальной матрицы на пересечении строки и столбца, соответствующих глобальным номерам узлов $i$ и $j$ элемента. Поскольку каждый узел связан лишь с несколькими соседями, в строке глобальной матрицы отлично от нуля лишь 5–10 элементов из тысяч — матрица разреженная, и хранятся только ненулевые элементы.

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

Решение системы

Полученная матрица симметрична и положительно определена. Для её решения применяют два класса методов.

Прямые методы (разложение Холецкого, $LU$ -разложение) дают решение за конечное число операций с точностью до ошибок округления. Их недостаток — заполнение: в ходе разложения нули разреженной матрицы превращаются в ненулевые элементы, и для больших трёхмерных задач памяти не хватает. Для плоских задач с сотнями тысяч узлов прямые методы остаются предпочтительными.

Итерационные методы строят последовательность приближений, сходящуюся к решению. Для симметричных положительно определённых матриц основным является метод сопряжённых градиентов: на каждом шаге требуется лишь умножение матрицы на вектор, память не растёт, а точность задаётся порогом невязки $\|\mathbf K\mathbf u-\mathbf f\|$ . Скорость сходимости зависит от обусловленности матрицы, которая ухудшается при сгущении сетки и при сильном контрасте материалов (сталь — воздух), поэтому метод применяют с предобусловливанием.

Решатель находит решение дискретной системы с заданной точностью — скажем, $10^{-8}$ по невязке. Это не означает, что решение краевой задачи получено с такой точностью: погрешность дискретизации, то есть разница между $\varphi_h$ и $\varphi$ , обычно на много порядков больше и не зависит от того, насколько точно решена система.

Нелинейные задачи

В магнитостатике с ферромагнетиком проницаемость $\mu$ зависит от индукции $B$ , которая сама подлежит определению. Матрица становится функцией решения: $\mathbf K(\mathbf u)\mathbf u=\mathbf f$ . Два основных подхода:

Метод простых итераций (Пикара): задаются начальным $\mu$ , решают линейную задачу, по найденному $B$ в каждом элементе определяют новое $\mu(B)$ из кривой намагничивания, собирают матрицу заново и повторяют до стабилизации. Сходится надёжно, но медленно, особенно при глубоком насыщении; часто требует релаксации — усреднения новых и предыдущих значений проницаемости в элементах.
Метод Ньютона: уравнение линеаризуется относительно поправки с учётом производной $\dd\mu/\dd B$ (дифференциальной проницаемости). Сходимость квадратичная вблизи решения — число верных знаков удваивается на каждом шаге, — но при неудачном начальном приближении метод может расходиться; на практике первые шаги делают по Пикару, затем переключаются на Ньютона.

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

Численное интегрирование

Для линейных треугольников интегралы в $K^e_{ij}$ вычисляются точно. При переменных в элементе свойствах ( $\varepsilon(\ur)$ , $\lambda(T)$ ), криволинейных элементах или элементах высшего порядка подынтегральные функции перестают быть полиномами низкой степени, и интегралы берут численно — квадратурами Гаусса: значение функции вычисляют в нескольких специально выбранных точках элемента и суммируют с весами. Квадратура из $m$ точек на отрезке точна для полиномов степени $2m-1$ ; недостаточное число точек — ещё один, обычно малый, источник погрешности.

Постобработка

Решение системы — набор узловых потенциалов. Всё остальное получается из него дифференцированием и интегрированием, и здесь кусочная природа $\varphi_h$ проявляется в полной мере.

Локальное поле. В каждом линейном элементе $\uE_h=-\sum u_i\nabla N_i$ — постоянный вектор. На ребре между элементами касательная составляющая $\uE_h$ непрерывна (потенциал на ребре одинаков с обеих сторон), а нормальная составляющая $\uD_h$ терпит скачок, хотя в точной задаче при отсутствии поверхностного заряда она непрерывна. Величина скачка — непосредственная мера локальной погрешности. Поле в узле обычно получают усреднением по примыкающим элементам; пространственная «ступенчатость» карты $|\uE|$ на грубой сетке — не физика, а артефакт аппроксимации.

Интегральные величины оказываются существенно точнее локального поля. Энергия $W=\tfrac12\int\varepsilon|\uE_h|^2\dV=\tfrac12\mathbf u^{\mathsf T}\mathbf K\mathbf u$ вычисляется из той же квадратичной формы, которую минимизирует метод Галёркина, и локальные ошибки поля в ней частично компенсируются. Ёмкость $C=2W/U^2$ , индуктивность $L=2W/I^2$ , потери $P=\int\sigma|\uE_h|^2\dV$ наследуют эту точность. Напротив, заряд, вычисленный как поток $\uD_h$ через контур вокруг электрода, зависит от выбора контура и использует разрывное поле напрямую, поэтому сходится медленнее. Количественно эти различия разобраны в уроке о сходимости.

Распространённые заблуждения

«Мельче сетка — всегда лучше». Сгущение уменьшает только погрешность дискретизации. Оно не исправит неверно выбранную область, упрощённую не там геометрию, ошибочные свойства материала или неподходящее граничное условие. Кроме того, слишком мелкая сетка ухудшает обусловленность матрицы и увеличивает влияние ошибок округления.
«Цветная картинка — значит, правильный ответ». Любая сетка даёт гладкую на вид карту потенциала; метод Галёркина гарантирует существование решения при любом качестве дискретизации. Правильность проверяется сравнением с аналитическими эталонами, сгущением сетки и балансами, а не внешним видом.
«Решатель выдаёт точное решение». Решатель с точностью до округления решает дискретную систему. Разница между дискретной и исходной задачами — погрешность дискретизации — остаётся и не видна по невязке системы.
«Поле в узле известно так же точно, как потенциал». Поле — производная; при линейных элементах оно теряет один порядок точности по сравнению с потенциалом и разрывно между элементами.
«Граница без условий — ошибка в постановке». Напротив: на ней автоматически выполняется однородное условие Неймана. Ошибка в том, чтобы не осознавать этого и считать такую границу «открытой».

Вопросы и упражнения

Для треугольника с вершинами $(0,0)$ , $(2,0)$ , $(0,2)$ и $\varepsilon=1$ найдите матрицу жёсткости. Сравните с матрицей разобранного примера и объясните, почему матрица не зависит от масштаба треугольника в плоской задаче.
Для того же треугольника при $\rho=3$ (условные единицы) вычислите вектор нагрузки. Что произойдёт с суммой компонент $f_i$ , если площадь элемента удвоить при той же $\rho$ ?
В слабой форме задача о проводнике с током $-\nabla\cdot (\sigma\nabla \varphi)=0$ имеет граничные условия: на электродах $\varphi=\pm U/2$ , на боковой поверхности — ничего не задано. Какое физическое условие при этом выполняется на боковой поверхности и почему оно соответствует изолятору?
Узловые потенциалы треугольника из примера равны $u=(0,\ 1,\ 1)$ В, сторона равна 1 м. Найдите $\uE_h$ в элементе и энергию $\tfrac12\mathbf u^{\mathsf T}\mathbf K^e\mathbf u$ при $\varepsilon=\varepsilon_0$ .
Решатель сообщает невязку $10^{-10}$ . Один инженер говорит, что ёмкость найдена с точностью $10^{-10}$ ; другой — что сетка дала 5 % погрешности. Могут ли оба быть правы одновременно? Что следует сделать, чтобы это выяснить?
Кривая $B$ – $H$ стали задана таблицей. Опишите словами один шаг метода Пикара для плоской задачи о магнитной цепи: какие величины вычисляются в элементах, что меняется в матрице и по какому признаку итерации останавливают.

Ответы к пунктам 1, 2 и 4

1. Здесь $2S=4$ , $S=2$ , $b_1=-2$ , $c_1=-2$ , $b_2=2$ , $c_2=0$ , $b_3=0$ , $c_3=2$ . По формуле $K_{ij}=(b_ib_j+c_ic_j)/(4S)=(b_ib_j+c_ic_j)/8$ : $K_{11}=8/8=1$ , $K_{12}=-4/8=-\tfrac12$ , $K_{13}=-\tfrac12$ , $K_{22}=\tfrac12$ , $K_{23}=0$ , $K_{33}=\tfrac12$ — та же матрица. Коэффициенты $b_i$ , $c_i$ пропорциональны линейному размеру $\ell$ , их произведения — $\ell^2$ , площадь — тоже $\ell^2$ ; отношение не зависит от масштаба. Это отражает тот факт, что в плоской задаче «проводимость» области между двумя точками (на единицу глубины) безразмерна относительно масштаба.

2. $f_i=\rho S/3=3\cdot2/3=2$ для каждого узла; сумма $\sum f_i=\rho S=6$ — полный заряд элемента. При удвоении площади сумма удвоится: $\sum f_i=\rho S$ всегда равна заряду элемента.

4. $\nabla \varphi_h=\sum u_i\nabla N_i=0\cdot(-1,-1)+1\cdot(1,0)+1\cdot(0,1)=(1,\ 1)$ В/м, $\uE_h=(-1,\ -1)$ В/м, $|\uE_h|=\sqrt2$ В/м. Энергия: $\mathbf K^e\mathbf u=(-1,\ \tfrac12,\ \tfrac12)$ , $\mathbf u^{\mathsf T}\mathbf K^e\mathbf u=0+\tfrac12+\tfrac12=1$ , откуда $W=\tfrac12\varepsilon_0\cdot1=4{,}43\cdot10^{-12}$ Дж/м (на единицу глубины). Проверка: $\tfrac12\varepsilon_0|\uE_h|^2S=\tfrac12\varepsilon_0\cdot2\cdot\tfrac12=\tfrac12\varepsilon_0$ — совпадает.

Итоги

Численное моделирование — цепочка: физическая система → концептуальная модель → математическая модель → дискретная модель → решение → постобработка → интерпретация. Ошибки модели не исправляются сеткой.
Усечение области, упрощение геометрии, материальные данные и граничные условия — модельные решения, каждое со своей погрешностью; граничное условие описывает то, что находится вне области.
МКЭ ищет решение в виде $\sum u_iN_i$ ; для линейного треугольника $N_i=(a_i+b_ix+c_iy)/(2S)$ , функции формы обладают свойством Кронекера и образуют разбиение единицы, их градиент постоянен в элементе.
Слабая форма получается умножением на пробную функцию и интегрированием по частям; условия Дирихле — главные, Неймана — естественные. Метод Галёркина даёт систему $\mathbf K\mathbf u=\mathbf f$ с $K^e_{ij}=\varepsilon(b_ib_j+c_ic_j)/(4S)$ , равносильную минимизации энергии.
Решатель точно решает дискретную систему; погрешность дискретизации остаётся. Поле кусочно-постоянно и разрывно; интегральные величины точнее локальных.

Далее — размерность моделей: как трансляционная и вращательная симметрия сводят трёхмерную задачу к плоской или осесимметричной и как правильно пересчитывать результаты на реальное устройство.