Численные модели в интроскопии

7. ФОРМИРОВАНИЕ МАТРИЦЫ КОНЕЧНОЭЛЕМЕНТНОГО МЕТОДА

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

F ?? Fi
i

Численные модели в интроскопии
Функционалы Fi суммируются по всем элементам
дискретной сети.
Предполагая, что сеть имеет N узлов, условие минимизации
функционала для Vk -ого узлового потенциала будет
выглядеть следующим образом
N
?F
?F
? ? i ?0
?Vk i?1 ?Vk

Отметим, что хотя сумма и формируется по всем
элементам сети, Fi отличен от нуля лишь в самом i-ом
элементе и равен нулю вне этого элемента

Численные модели в интроскопии
Функционал для элемента равен нулю вне самого элемента.
Производные функционала по отношению к потенциалу
k-го узла отличны от нуля только для элементов a, b, c, f, g и h

N
?F
?Fi ?Fa ?Fb ?Fc ?F f ?Fg ?Fh
??
?
?
?
?
?
?
?Vk i?1 ?Vk ?Vk ?Vk ?Vk ?Vk ?Vk ?Vk

Численные модели в интроскопии
7.1. Вывод матричных выражений в электростатике
Дискретизированный функционал для
двумерной задачи электростатического поля имеет вид
?1 2
?
Fi ? ?? ?E ? ?V ? dS
2
?
Si ?

где для двумерного случая объем равен численному
значению площади i-го элемента Si
(глубина элемента равна единице)

Численные модели в интроскопии
?Fi
Рассчитаем ?Vk , предполагая, что

Vk - потенциал

одного из трех узлов 1, 2 или 3 элемента i
?Fi
? ?1 2
?
?
?
E
?
?
V
?
? dS
?
?Vk ?Vk Si ? 2
?

Рассчитаем первое слагаемое функционала,
используя выражение для напряженности
электрического поля
1 3
E ?? ? ? iql ? jrl ?Vl
D l ?1

(вектор E постоянен в пределах элемента!!!)

Численные модели в интроскопии
? ?1 2?
?
?
E
dS
?
?
?
?
?Vk Si ? 2
?Vk
?

2

1 D ?E 2
?1 2?
? ?E ? ?dS ? ?
2 2 ?Vk
?2
? Si

3

2

?E
? 1 ?
?
?
?
?
i
q
?
j
r
V
l
l l? ?
2 ??
?Vk ?Vk D ? l ?1
?
2
2
3
3
?
? 1 ?
? ?
? ?
?
? q V ? ? ? ? rlVl ? ?
2 ? ? l l
?Vk D ??? l ?1
? ? l ?1
? ??

Численные модели в интроскопии
2
2
3
3
?
? ?
? ?
?E
1 ? ?
? 2
?? ? qlVl ? ? ? ? rlVl ? ?
?Vk D ?Vk ??? l ?1
? ? l ?1
? ??
2

Поскольку k-ый узел является одним из трех
узлов в элементе, то
3
?E 2
1 ? 3
?
? 2 ? 2? qlVl qk ? 2? rlVl rk ?
?Vk D ? l ?1
l ?1
?

Численные модели в интроскопии
В результате для первой составляющей
производной функционала получим
? ?1 2?
1 D ?E 2
?
? ?E ? dS ? ?
?
?Vk Si ? 2
2 2 ?Vk
?

?
? ? q1V1 ? q2V2 ? q3V3 ? qk ? ? r1V1 ? r2V2 ? r3V3 ? rk ? ?
2D
?
?
? ? q1qk ? r1rk ?V1 ? ? q2 qk ? r2 rk ?V2 ? ? q3qk ? r3rk ?V3 ?
2D
?

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

?
? ? q1qk ? r1rk ?
2D

? q2 qk ? r2 rk ?

?V1 ?
? q3qk ? r3rk ?? ??V2 ??
??V3 ??

В численной процедуре при учете элемента
это матричное выражение записывается для всех трех узлов

? q1q1 ? r1r1 q1q2 ? r1r2
? ?
q
q
?
r
r
q
q
?
r
r
2
1
2
1
2
2
2
2
?
2D
?? q3q1 ? r3r1 q3q2 ? r3r2

q1q3 ? r1r3 ? ?V1 ?
q2 q3 ? r2 r3 ?? ??V2 ??
q3q3 ? r3r3 ?? ??V3 ??

Численные модели в интроскопии
Матрица коэффициентов, умножаемая на вектор
потенциалов, является симметричной.
Это так называемая элементная матрица, или
матрица элемента, обычно обозначается через

?
? qn qk ? rn rk ?
S ? n, k ? ?
2D
где n и k - номера строк и столбцов в матрице

S ? n, k ?

Численные модели в интроскопии
Рассчитаем второе слагаемое производной функционала,
предполагая, что объемная плотность заряда ?
постоянна в элементе

?
?Vk

?
? ? ?V ? dS ?? ?
?
?Vk
Si

??V ? dS

Si

Заметим, что потенциал есть сумма трех составляющих
3

V ? x, y ? ?? ? pl ? ql x ? rl y ?Vl
l ?1

Численные модели в интроскопии
Обозначая через коэффициент ?1 , умножаемый
на потенциал V1 , то есть

1
?1 ? x, y ? ? ? p1 ? q1 x ? r1 y ?
D
получим

3

V ? x, y ? ?? ? lVl ??1V1 ? ? 2V2 ? ? 3V3
l ?1

x ? x1 ,

y ? y1

?

V(x1,y1 ) ?V1

V(x1,y1 ) ?1?V1 ? 0 ?V2 ? 0 ?V3

?1 ? x1 , y1 ? ?1, ? 2 ? x1 , y1 ? ?0, ? 3 ? x1 , y1 ? ?0

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

x ? x2 ,

y ? y2 ? ?1 ? x2 , y2 ? ?0, ? 2 ? x2 , y2 ? ?1, ? 3 ? x2 , y2 ? ?0

x ? x3 ,

y ? y3 ?

?1 ? x3 , y3 ? ?0, ? 2 ? x3 , y3 ? ?0, ? 3 ? x3 , y3 ? ?1

Таким образом,
?1 есть двумерная (линейная) функция,
равная 1 для x ? x1 , y ? y1 и нулю в других узлах

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

Численные модели в интроскопии
Возвращаясь к рассчитываемому выражению, получим

?
?Vk

?
? ? ?V ? dS ?? ?
?
?Vk
Si

??V ?dS ?

Si

?
?? ? ? ??1V1 ? ? 2V2 ? ? 3V3 ? dS ?? ? ?? k dS
?Vk
Si
Si
где k - узел треугольника. С учетом формы коэффициента

?k
как показано, интеграл равен объему тетраэдра
D1
1
? ? ?? k dS ?? ?
?? ?D
23
6
Si

Численные модели в интроскопии
Повторяя подобные расчеты для всех трех узлов элемента
и добавляя их к матрице, окончательно получим
? ?Fi ?
?
?
?
V
? 1 ? ? S11
? ?Fi ? ?? S
? ?V2 ? ? 21
? ?F ? ?? S31
? i?
?? ?V3 ??

S12
S 22
S32

S13 ? ?V1 ?
?? ?
D? ?
?
?
?
S 23 ? ??V2 ? ? ? ? ?
6
?? ? ??
S33 ?? ??V3 ??

Это выражение является базовым строительным блоком
для формирования системы уравнений
в методе конечных элементов

Численные модели в интроскопии
7.2. Матрица элемента в задаче со стационарными токами
Это явление характеризуется энергетическим
функционалом вида
?1
?
Fi ? ?? ?E 2 ? dS
2
?
Si ?

По аналогии, расчет

?Fi ?Vk

дает матрицу эквивалентного вида

?
? qn qk ? rn rk ?
S ? n, k ? ?
2D

Численные модели в интроскопии
7.3. Задача магнитного поля: скалярный потенциал
?H
?
Fi ? ?? ?BdH ?dS
Si ? 0
?

Оценим производную функционала
H
?
?Fi
? ?
??
? ?BdH ?dS
?Vk S i ?Vk ? 0
?

с учетом зависимости H от Vk , получим
H
? ?H
?Fi
? ?
?H
? ? ? ?BdH ?
dS ? ?B
dS
?Vk Si ?H ? 0
?Vk
Si
? ?Vk

Численные модели в интроскопии
Для конечных элементов первого порядка

B

и
H
постоянны внутри каждого элемента
?Fi D ?H D
?H ?D ?H 2
? B
? ?H
?
?Vk 2 ?Vk 2
?Vk
4 ?Vk
?H 2
Расчет ?Vk приводит к известному выражению,

которое дает аналогичное матричное уравнение
?
? qn qk ? rn rk ?
S ? n, k ? ?
2D

Численные модели в интроскопии
7.4. Задача магнитного поля: векторный потенциал
Локальный энергетический функционал:
?H
?
Fi ? ?? ?HdB ? J ?A ?dS
Si ? 0
?

Выражение для ?Fi ?Ak получается суммированием частных
производных двух составляющих
Первая составляющая аналогична полученным выражениям
1
Отмечая, что H=?B, с учетом ? ? получим
?

?
? qn qk ? rn rk ?
S ? n, k ? ?
2D

Численные модели в интроскопии
?
?Ak

Производная

подобна производной

?? J ?AdS

Si

?
?Vk

?? ? ?V dS

Si

Используя ранее полученные результаты,
получим для ?Fi ?Ak следующее матричное выражение
? ?Fi ?
?
?
?
A
? 1 ? ? S11
? ?Fi ? ?? S
? ?A ? ? 21
? ?F2 ? ?? S 31
? i?
?? ?A3 ??

S12
S 22
S 32

S13 ? ? A1 ?
?J ?
D? ?
?
?
?
S 23 ? ?? A2 ? ? ? J ?
6
?
?
?
?? J ??
S 33 ? ? A3 ?

Численные модели в интроскопии
7.5. Поле постоянного магнита
Предположим постоянство магнитной проницаемости.
Тогда соответствующий функционал будет иметь вид
1
1
2
Fi ? ? ? ? B ? Br ? dS ? ? ?Bi2 dS
Si 2
Si 2

Предположим также, что в i-ом элементе магнита
отсутствуют токи
?Fi ?Fi ?Bi D ? ?Bi2 ?Bi D?
?Bi
?
?
?
2 Bi
?Ak ?Bi ?Ak 2 2 ?Bi ?Ak
4
?Ak

Численные модели в интроскопии
В двумерном случае используется только одна
компонента А, перпендикулярная плоскости 0xy.
Рассматривая А как скалярную величину, можем записать
?Bi
?B i
Bi
?B i
?Ak
?Ak

потому что вектора

Bi

и

?B i
?Ak

(производная вектора

по отношению к скаляру) имеют то же направление
1 3
B i ?rotA ? B r ? ? ? irl ? jql ? Al ? iBrx ? jBry
D l ?1

Численные модели в интроскопии
B r ?iBrx ? jBry

?B i 1
? ? irk ? jqk ?
?Ak D

?Fi D?
?B i
?
2B i
?
?Ak
4
?Ak
D?
?
2

?1 3
? ?1
?
?
?
?
?
i
r
?
j
q
A
?
i
B
?
j
B
?
i
r
?
j
q
l
l
rx
ry ? ?
k
k ?
?D ? l
?
? l ?1
? ?D

После алгебраических преобразований получим
?Fi D? ? 1 3
? ?1
?
?
?
?
?
?
i
r
?
j
q
A
?
i
B
?
j
B
?
i
r
?
j
q
?
l
l
l
rx
ry ? ?
k
k ? ?
?
?Ak
2 ? D l ?1
?
? ?D
? 3
?
? Bry qk ? Brx rk ?
?
?
?
r
r
?
q
q
A
?
?
l k
l k
l
2 D l ?1
2

Численные модели в интроскопии
В матричной форме это выражение можно переписать
в виде
? ?Fi ?
?
?
? ?A1 ? ? S11
? ?Fi ? ?? S
? ?A ? ? 21
? ?F2 ? ?? S 31
? i?
?? ?A3 ??

S12
S 22
S 32

? Bry q1 ? Brx r1 ?
S13 ? ? A1 ?
??
?
?
?
?
S 23 ? ?? A2 ? ? ? Bry q2 ? Brx r2 ?
2
? Bry q3 ? Brx r3 ?
S 33 ?? ?? A3 ??
?
?

Tk

Численные модели в интроскопии
7.6. Матрица элемента в задаче с электрическим
векторным потенциалом
Этот случай описывается локальным энергетическим
функционалом вида
? 1 2?
Fi ? ??
J ? dS
2?
?
Si ?

Расчет ?Fi ?Tk , где Tk - электрический векторный потенциал
в точке k, может быть проведен подобно тому,
как он был получен в задаче с магнитным
векторным потенциалом,
в предположении постоянной величины ?

Численные модели в интроскопии
При этих допущениях локальный функционал
рассчитывается по формуле
? 1 2?
Fi ? ???
B ?? dS
2? ?
Si ?

Следовательно, условная запись S ? n, k ? для выражения
также будет аналогична полученной ранее, то есть
S ? n, k ? ?

1
? qn qk ? rn rk ?
2 D?

Численные модели в интроскопии
7.7 Формирование глобальной матрицы системы
Минимизация энергетических функционалов позволяет
получить элементные матрицы для различных физических
задач. Выражения ?Fi ?Vk
являются матричными уравнениями и
могут быть записаны в общей форме следующим образом
? ?Fi ?
?
?
?
V
? 1 ? ? S11
? ?Fi ? ?? S
? ?V ? ? 21
? ?F2 ? ?? S 31
? i?
?? ?V3 ??

S12
S 22
S 32

S13 ? ?V1 ? ? Q1 ?
S 23 ?? ???V2 ?? ? ??Q2 ??
S 33 ?? ??V3 ?? ?? Q3 ??

Численные модели в интроскопии

? Q?

- это так называемый "вектор-источник" конкретной
задачи, который, в зависимости от постановки, может быть
электрическим зарядом, током или постоянным магнитом.
Уравнение, рассчитанное для i-го элемента, является
частью суммы
N
?F
?F
? ? i ?0
?Vk i?1 ?Vk

где N - число элементов в конечно-элементной сети.
Чтобы рассчитать эту сумму, надо объединить вклады
всех элементов в общую глобальную матрицу системы,
в которой учтены все узлы сети

Численные модели в интроскопии
Объединение локальной матрицы в глобальную.
Вклад элемента с номерами узлов 7, 10 и 14

Численные модели в интроскопии
После того, как все N элементов просчитаны и объединены,
образуется глобальная система

? SS ? ??V ? ?? Q?
Перед тем, как приступить к решению, необходимо учесть
граничные условия Дирихле - это может быть сделано
с помощью метода, основанного на следующем правиле:
• если значение Vm известно, диагональный элемент
строки с номером m устанавливается в единицу;
• все другие коэффициенты в строке обнуляются;
• в строку m вектора ? Q ? вводится значение Vm

0 ?V1 ? 0 ?V2 ? ... ? 1?Vm ? 0 ?Vm ?1 ? ... ? 0 ?VK ?Vm






Чтобы не видеть здесь видео-рекламу достаточно стать зарегистрированным пользователем.
Чтобы не видеть никакую рекламу на сайте, нужно стать VIP-пользователем.
Это можно сделать совершенно бесплатно. Читайте подробности тут.