3.2 Конечно-элементная дискретизация уравнений

Сначала построим дискретную модель для задачи

 = 0 в объеме , (3.24)

на границе сверхпроводника (3.25)

₊ – _– = I на поверхности разреза (3.26)

Минимизируя функционал задачи

(3.27)

на множестве функций вида

(3.28)

(_i – коэффициенты, N_i – базисные функции – полиномы), придем к системе уравнений

, (3.30)

Рис. 3.5.

где . Матрица {S_ij} является симметричной и положительно определенной.

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

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

.

Следуя стандартному методу конечных элементов [93], проведем объединение этих систем. Если S не является линией разреза, то в узлах 2 и 3 должно выполняться условие непрерывности, т.е.

(3.30)

В этом случае получается система (с учетом обозначений _i  _i^(l))

(3.31)

Если на S происходит скачок, то каждому узлу, лежащему на S, необходимо поставить в соответствие два значения потенциала _i+ и _i– . Одно из этих значений будет относиться к элементам, лежащим по одну сторону от S, другое – к элементам, лежащим по другую сторону относительно S. В нашем случае условие скачка можно записать в виде

(3.32)

При объединении элементов для каждого узла на S следует оставить только одну независимую переменную, например, , а и заменить соответственно выражениями ₂ – I, ₃ – I.

В итоге получим систему

(3.33)

Таким образом, по сравнению с (3.31) данная система отличается только правой частью. Решений она не имеет (объясните, почему?). Для получения нужного решения следует заметить, что дополнительное условие (3.26) накладывается только на разность значений потенциала, а не на саму функцию. Тем самым, решение задачи определяется с точностью до константы. Поэтому необходимо в одном каком-либо узле зафиксировать потенциал, т.е. придать ему вполне определенное значение, которое играет роль начала отсчета. Например, поставим в узле 4 условие

₄ = 0.

Тогда система примет вид

(3.34)

Теперь система имеет отличный от нуля определитель. В случае, если I = 0 (скачка нет), ее решение – тождественный нуль, если I  0 (скачок есть), получим искомое ненулевое решение задачи.

По итогам решения системы конечно-элементных уравнений найдем значение  в произвольной точке элементов ₁ и ₂ согласно (3.28):

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

Добавим к элементам 1 и 2 (рис. 3.5) еще два элемента, как показано на рис. 3.6. Выпишем второе уравнение глобальной системы; это уравнение соответствует узлу 2, который лежит на S и является общим для всех четырех элементов. С этим узлом связаны потенциалы ₂₊ и _2– , причем

Вводя обозначение ₂  ₂₊ , получим следующее уравнение

.

Сформировав остальные уравнения и положив, например, ₄ = 0, получим новую систему, решение которой даст набор чисел {_i}. При вычислении потенциала внутри элементов 1 и 3 в качестве узлового значения узла 2 следует брать ₂ , а внутри элементов 2 и 4 – величину (₂ – I).

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

Осесимметричная формулировка отличается от плоской только способом вычисления матричных элементов :

.

Рис. 3.6.

В рассмотренной постановке ток I должен быть известным. Однако такая модель непригодна для описания случая «замороженного» потока, когда поток Ф₀ через отверстие замкнутого сверхпроводника фиксирован, а в самом сверхпроводнике устанавливается такой ток, который препятствует изменению потока Ф₀. В такой постановке ток I является изменяемым параметром, который по-прежнему учитывается как скачок потенциала.

Получим вариационную фор­мулировку такой задачи. Пусть функционал имеет вид

(3.35)

где I  (₊ – _–)|_Г2 – скачок потенциала на поверхности разреза Г₂; в различных точках этой поверхности значения потенциала могут заметно отличаться, но разность ₊ и _–, вычисленных в одной точке, остается постоянной. Поэтому слагаемое Ф₀I в приведенном выражении есть число, как того и требует понятие функционала. Для первой вариации получим

I| _Г2 есть вариация величины скачка потенциала на Г₂. Так как эта величина неизвестна, то I  0. Перепишем F, используя формулу Грина

Здесь Г₂₊ и Г_2– – разные стороны линии (поверхности) разреза, а Г – остальная часть границы (не включающая разрез). Очевидно, что направления вектора нормали к поверхностям Г₂₊ и Г_2– противоположны, тем самым

.

Вариация I = const на Г₂, так как скачок потенциала I постоянен в каждой точке поверхности Г₂. Таким образом,

Итак, для того, чтобы функционал F достигал экстремума, необходимо:

в области  вне сверхпроводника, (3.36а)

на поверхности Г' сверхпроводника (3.36б)

(3.36в)

(знак «+» или «–» определяется выбором направления нормали к поверхности Г₂).

Если наложить на потенциал главное граничное условие в виде

на внешней границе области ,

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

Рассмотрим особенности учета условия сохранения потока в реализации МКЭ. Как известно, стандартная конечно-элементная процедура состоит в вычислении матричных элементов и вектора правых частей для каждого конечного элемента (т. н. локальной системы линейных уравнений) и добавлении их в глобальную систему, относящуюся ко всей конечно-элементной сетке. Отметим, что в рассматриваемой формулировке поверхность разреза Г является внутренней по отношению к  и тем самым должна находиться внутри конечно-элемент­ной сетки, но при этом она не может проходить по внутренней части элементов, а может проходить только по границам элементов. Возьмем два двумерных конечных элемента, примыкающих к поверхности Г с разных сторон, как показано на рис.3.7.

Если бы условия (3.36) не было (в функционале нет второго слагаемого), то локальные системы этих элементов имели вид

г

Рис. 3.7.

де S^(e) – матрица жесткости элемента с номером e; индексы i, j, k совпадают с глобальными номерами узлов элемента e. Глобальная система уравнений, соответствующая этим элементам, выглядит так:

Если Г является поверхностью разреза с условием (3.36), то в узлах 1 и 3 значения  справа и слева не равны друг другу, а связаны соотношением

₁₊ – _1– = I, ₃₊ – _3– = I.

Так что по одной узловой переменной можно исключить, например ₁₊ и ₃₊ . Чтобы правильно построить глобальную систему, запишем выражения для элементных функционалов

Заметим, что каждому узлу на Г соответствует только одна независимая переменная – та, что вычисляется по левую сторону. Но при этом появилась дополнительная переменная, которая ставится в соответствие не узлам, а всему разрезу в целом – сила тока I. Суммарный функционал представится в виде

F = F⁽¹⁾ + F⁽²⁾ + Ф₀I.

Минимизируя F по переменным _1– , ₂ , _3– , ₄ и I, получим систему

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

,

,

и

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

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

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

Как видно, такой подход требует решения дополнительно k(k+1)/2 задач (k – число токонесущих сверхпроводников; здесь k = 2) для определения коэффициентов само- и взаимной индукции из (3.5). (Укажите способы вычисления этих коэффициентов).