Scientific journal
Modern high technologies
ISSN 1812-7320
"Перечень" ВАК
ИФ РИНЦ = 0,940

MATHEMATICAL MODELING OF MAGNETIC FIELD USING BOUNDARY INTEGRAL EQUATION METHOD AND JILES-ATHERTON MODEL

Denisov P.A. 1
1 Platov South-Russian State Polytechnic University (NPI)
Modification of boundary integral equation method was developed to simulate a stationary plane-parallel magnetic field in nonlinear medium with hysteresis, characterized by the use of scalar potentials of simple and double layer in combination with explicit expressions of the Jiles-Atherton model. The proposed method of accounting for hysteresis consists in the extracting from the area bounded by the hysteresis loop, a unambiguous characteristic defined using explicit expressions of the Jiles-Atherton model. As a result of the shifting of the obtained characteristic along the axis of magnetic intensity, the ferromagnetic element is represented as an overlay of a nonlinear ferromagnetic with a characteristic passing through the zero, and a permanent magnet with a magnetization equal to the value of the coercivity. The calculation of the magnetic field distribution is reduced to the finding of scalar secondary sources in the form of fictitious magnetic charges placed at the boundaries of the regions of homogeneity. The potential of a simple-layer was used to simulate a field of known magnetization, and a double-layer potential was used for the ferromagnetic reaction field. The computer program has been developed and the test problem has been calculated in order to validate the proposed method. The results were compared with those obtained using the finite element analysis software Femm 4.2.
mathematical modeling
magnetic field
hysteresis
Jiles-Atherton model
integral equation method
scalar single-layer potential
scalar double-layer potential

Актуальность исследования и развития математических моделей гистерезиса и методов моделирования физических полей связана с широким использованием в технических устройствах элементов, выполненных из магнитных материалов и диэлектриков, имеющих неоднозначные зависимости характеристик. Современный уровень развития техники предполагает все более активное использование компонентов и материалов на пределе их возможностей, что требует широкого применения компьютерного моделирования при проектировании и совершенствовании математических моделей с целью повышения их соответствия исследуемым объектам [1, 2].

Целью исследования является усовершенствование инструментов моделирования магнитного поля с учетом гистерезиса. В рамках поставленной цели выполнено объединение в единый метод ранее разработанных модификаций модели Джилса – Атертона [3] и метода моделирования магнитного поля, построенного на основе граничных интегральных уравнений [4].

Материалы и методы исследования

Согласно модели Джилса – Атертона, в объеме ферромагнитного материала векторы магнитной индукции denis01.wmf, напряженности магнитного поля denis02.wmf и намагниченности denis03.wmf коллинеарны, а их проекции на общую ось удовлетворяют соотношению

denis04.wmf (1)

Связь M(H) задается дифференциальными соотношениями, основанными на безгистерезисной кривой Man(H), для описания которой была выбрана функция Ланжевена

denis05.wmf

где Ms – намагниченность насыщения, a – параметр формы. Также в модели используются коэффициент сцепления k, параметр c, характеризующий выпуклость границ домена, и коэффициент междоменной связи α. Перечисленные параметры можно найти по измеренной петле гистерезиса, используя метод Нелдера – Мида [5] безусловной оптимизации.

В случае, когда междоменным взаимодействием можно пренебречь и положить α = 0, для описания кривой M(H), выходящей из известной начальной точки (H0, M0), можно использовать явные выражения, устойчивые к ошибкам численной реализации [3].

При denis06.wmf рабочая точка будет двигаться вправо, т.е. в сторону возрастания H. Возможны два случая. Если denis07.wmf, то рабочая точка всегда будет находиться под Man(H), а ее движение будет описываться кривой

denis08.wmf (2)

где denis09.wmf. Если denis10.wmf, то траектория M(H) и кривая Man(H) пересекутся в точке (Hcr, Mcr), определяемой из системы уравнений

denis11.wmf, (3)

denis12.wmf. (4)

При denis13.wmf траектория будет описываться кривой

denis14.wmf, (5)

а при H ≥ Hcr – выражением (4), где вместо (H0, M0) необходимо подставить (Hcr, Mcr).

При denis15.wmf рабочая точка будет двигаться влево, т.е. в сторону убывания H. Возможны два случая. Если denis16.wmf, то рабочая точка всегда будет находиться над Man(H), а ее движение будет описываться кривой

denis17.wmf, (6)

где denis18.wmf. Если denis19.wmf, то при denis20.wmf траектория определяется согласно (5), а при H ≤ Hcr – выражением (6), где вместо (H0, M0) необходимо подставить (Hcr, Mcr), определяемую из уравнений (3),(4).

Для общего случая при α ≠ 0, для определения M(H) можно использовать итерационный процесс:

denis21.wmf, denis22.wmf, (7)

denis23.wmf, (8)

где функции M(Hj) определяются согласно (2)–(6). Результатом является функция denis24.wmf, которая может быть подставлена в (1) для определения индукции.

Вычислительные правила (1)–(8) позволяют для любой начальной точки (H0, B0) выделить из петли гистерезиса однозначную зависимость B(H), представляющую собой непрерывную, возрастающую, гладкую всюду, кроме точки (H0, B0), функцию (рис. 1, а). Такая зависимость соответствует характеристике ферромагнитного материала при любом монотонном изменении значения напряженности H магнитного поля относительно начального заданного значения H0 и может быть использована для моделирования магнитного поля постоянных токов и постоянных магнитов. Каждая выделенная линия B(H) имеет одну точку пересечения с осью H – остаточную индукцию Br, и с осью B – коэрцитивную силу Hc.

den1a.tif den1b.tif

а) б)

Рис. 1. Модель Джилса – Атертона магнитного гистерезиса: 1 – начальная кривая намагничивания, 2 – петля гистерезиса, 3 – выделенная однозначная зависимость, 4 – модель Джилса – Атертона: а) однозначная характеристика B(H) для начальной точки (H0, B0), б) петля гистерезиса и начальная кривая для модели Джилса – Атертона

В доступных программах анализа физических полей Femm и Elcut для моделирования кривой размагничивания B(H), не проходящей через начало координат, используется подход, предполагающий смещение этой кривой вправо на величину коэрцитивной силы [6, 7]. Магнитная проницаемость μ в этом случае определяется по формуле

denis25.wmf. (9)

Смещение кривой размагничивания в [6, 7] достигается размещением по границам магнита эффективных поверхностных токов с поверхностной плотностью ±Hc. Установлено, что аналогичный эффект достигается заданием остаточной намагниченности, значение которой denis26.wmf. Такой подход использован в работе для расчета поля для ферромагнетика с выделенной однозначной характеристикой B(H).

Для анализа распределения плоскопараллельного магнитного поля была использована модификация метода граничных интегральных уравнений, предполагающая одновременное использование скалярных потенциала простого слоя (ППС) и потенциала двойного слоя (ПДС) [4]. Согласно выбранному методу магнитная индукция ищется в виде суперпозиции:

denis27.wmf

где denis28.wmf – безвихревое поле реакции ферромагнетиков, а denis29.wmf – поле постоянных магнитов, создаваемое в однородной среде с магнитной проницаемостью μ0.

Постоянный магнит представляется в виде совокупности элементарных однородно намагниченных постоянных магнитов, границы каждого из которых покрываются слоем фиктивных магнитных зарядов простого слоя с известной плотностью σ = μ0(Mr)n, где (Mr)n – нормальная составляющая остаточной намагниченности. Индукция, наводимая магнитом с индексом PM в точке Q, определяется по формуле

denis30.wmf,

где denis31.wmf вычисляется интегрированием по границам denis32.wmf магнита:

denis33.wmf,

denis34.wmf,

denis35.wmf – ППС, denis36.wmf – расстояние от точки интегрирования P до точки наблюдения Q, denis37.wmf вне магнита PM, а внутри denis38.wmf,

denis39.wmf

где denis40.wmf, denis41.wmf – проекции вектора denis42.wmf на оси координат x и y, соответственно.

Каждый ферромагнетик разбивается на элементарные области однородности магнитной проницаемости μ, а поле реакции denis43.wmf моделируется введением скалярного ПДС, двойной слой зарядов которого размещается на границах denis44.wmf областей однородности:

denis45.wmf,

denis46.wmf,

denis47.wmf,

где denis48.wmf – нормаль в точке P, τ – неизвестная плотность ПДС, которая определяется из интегрального уравнения второго рода

denis49a.wmf

denis49b.wmf, (10)

где точка Q пробегает границы всех элементарных областей разбиения ферромагнетиков denis50.wmf, denis51.wmf, denis52.wmf и denis53.wmf – магнитные проницаемости со стороны положительного и отрицательного направлений нормали denis54.wmf в точке Q соответственно.

На основе приведенных формул в работе получен алгоритм моделирования магнитного поля, предполагающий выполнение следующих этапов:

1) замена рассматриваемой области расчета кусочно-однородной средой, в том числе разбиение ферромагнетиков на области однородности SE и задание в каждой из них начальных значений магнитной проницаемости;

2) размещение первичных источников магнитного поля (плотностей ППС на границах элементарных постоянных магнитов);

3) определение плотностей скалярного ПДС на границах элементарных областей SE путем решения интегрального уравнения (10);

4) расчет усредненного значения магнитной индукции в каждой области SE:

denis55.wmf,

где denis56.wmf – внешняя нормаль, а denis57.wmf и denis58.wmf – магнитные проницаемости с внутренней и внешней сторон границы LE в точке denis59.wmf;

5) корректировка магнитной проницаемости в каждой элементарной области SE с применением параметра релаксации ν:

denis60.wmf,

где μk и μk+1 – старое и новое значения магнитной проницаемости соответственно, μ(Hcp) рассчитывается согласно (9), а значение Hcp определяется из уравнения denis61.wmf с помощью метода половинного деления отрезка (дихотомии). Выполнение алгоритма завершается, когда для каждой элементарной области однородности SE выполняется denis62.wmf, где ε – заранее заданное максимально допустимое значение погрешности определения магнитной проницаемости. Иначе переход на шаг в).

Результаты исследования и их обсуждение

Рассмотрена задача моделирования плоскопараллельного магнитного поля уединенного ферромагнитного образца C45 прямоугольной формы с остаточной намагниченностью, направленной вдоль оси y (рис. 2, а). Результаты сравнивались с решением, полученным в Femm 4.2. Для материала C45 на тороидальном образце с использованием прибора MagHyst была получена петля гистерезиса с начальной кривой намагничивания (рис. 1, б), определены значения остаточной индукции и коэрцитивной силы: Br = 0,949 Тл, Hc = –682,9 А/м. В Femm была экспортирована кривая размагничивания со смещением вправо на величину denis63.wmf.

По имеющейся петле гистерезиса были получены параметры модели Джилса – Атертона: Ms = 1460000 А/м, a = 1230 А/м, k = 870 А/м, c = 0,17, α = 0,0023, количество повторений n в (8) принято равным 6 (рис. 1, б). Построенная модель характеризовалась следующими значениями остаточной индукции и коэрцитивной силы: Br = 0,939 Тл, Hc = –687,2 А/м. В качестве однозначной характеристики B(H) была выбрана кривая размагничивания модели Джилса – Атертона. В модели была принята кусочно-постоянная аппроксимация плотности ПДС, количество граничных элементов оказалось равным 336.

den2a.tif den2b.tif

а) б)

Рис. 2. Расчет магнитного поля уединенного образца C45: 1 – результат, полученный предложенным методом, 2 – график, импортированный из FEMM 4.2: а) область расчета: ферромагнетик, разбитый на 21 подобласть, б) графики y – составляющей магнитной индукции вдоль оси y

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

Заключение

В работе предложена модификация метода граничных интегральных уравнений, предназначенная для моделирования стационарного плоскопараллельного магнитного поля в нелинейных средах с гистерезисом и основанная на выделении из области, ограниченной петлей гистерезиса, однозначной характеристики, задаваемой с использованием явных выражений модели Джилса – Атертона. Комбинированное использование скалярных ППС и ПДС совмещает преимущества обоих способов размещения зарядов, а именно:

- устойчивость решения в области, занимаемой ферромагнетиком, обеспечиваемую использованием ПДС;

- непрерывность и конечность значений плотности фиктивных магнитных зарядов ПДС, что упрощает ее аппроксимацию;

- простоту моделирования поля постоянных магнитов, обеспечиваемую применением ППС.

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

Результаты работы получены при поддержке ФГБОУ «Фонд содействия развитию малых форм предприятий в научно-технической сфере» в рамках НИР «Разработка программных модулей моделирования размагничивания постоянных магнитов для проектирования высокоэффективных электроприводов мотор-колес» (11690ГУ/2017).