FEMM 5

Исследование 6 · FEMM 5 · встроенный Python

Нелинейный трансформатор: четыре вещественные неизвестные и полный метод Ньютона

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

В этом примере используются собственная геометрия курса, её однозначная восьмиточечная B–H-кривая, нагрузка 10 Ом и внешнее сопротивление первичной обмотки 0,8 Ом. Источники 100 и 105 В здесь амплитудные: это примерно 70,71 и 74,25 В действующего значения.

Место в курсе.Геометрию, физическую постановку и результаты можно изучать вместе с основным занятием. Серийную программу и обработку разбирайте после введения во встроенный Python.

Объект, физическая модель и вывод соотношений

Соглашения об обозначениях и единицах
  • \(i(t)\), \(u(t)\) — мгновенные значения; \(\widehat I\), \(\widehat U\) — комплексные амплитуды; \(\underline I_{\text{д}}\), \(\underline U_{\text{д}}\) — комплексные действующие значения. Для синусоиды \(\widehat I=\sqrt{2}\,\underline I_{\text{д}}\); действующее значение тока равно \(|\underline I_{\text{д}}|\). Нормировка обязательно учитывается в формулах мощности.
  • \(\Phi=\int_S\mathbf B\cdot\mathrm d\mathbf S\) — поток через ориентированную поверхность; \(\Psi=\sum_k\oint_{C_k}\mathbf A\cdot\mathrm d\mathbf l\) — потокосцепление последовательных витков с согласованными замкнутыми контурами. Только при одинаковом сцепленном потоке каждого витка \(\Psi=N\Phi\). Обе величины измеряются в веберах: число витков безразмерно. Запись «Вб·виток» подчёркивает суммирование.
  • \(\theta_{\text{м}}\), \(\Omega\) (также \(\omega_{\text{мех}}\)) — механические угол и угловая скорость. Для основной гармоники \(\theta_{\text{э}}=p\theta_{\text{м}}\), \(\omega_{\text{э}}=p\Omega\), где \(p\) — число пар полюсов. Геометрический период зубцов задаётся отдельно. В производных по углу используйте радианы.
  • Статическая индуктивность \(L_{\text{ст}}=\Psi/i\) — наклон секущей из начала координат; это определение используется при \(\Psi(0)=0\) и \(i\ne0\). Динамическая индуктивность \(L_{\text{дин}}=\partial\Psi/\partial i\) — наклон касательной (дифференциальная, инкрементальная); положение и остальные токи фиксированы. Для линейной характеристики без смещения обе равны одной индуктивности \(L\). При постоянных магнитах и взаимном возбуждении сначала учитывают ненулевое смещение потокосцепления.
  • Положительные направления, оси и специальные индексы определяются в конкретной задаче. Единицу берите из условия или заголовка таблицы; перед подстановкой в формулу СИ переводите миллиметры в метры, миллигенри в генри и градусы в радианы.
  • Численные результаты обычно показаны с четырьмя значащими цифрами; для сравнения близких значений в таблицах сходимости сохранено больше цифр. Полные значения остаются в исходных расчётных данных. Малые остатки вычисления координат на осях показаны как ноль; малые физические величины сохраняются в показательной записи.

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

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

Линеаризация остатка в окрестности текущего вектора x даёт r(x+\Delta x)\approx r(x)+J\Delta x. Условие нулевого линейного остатка приводит к системе J\Delta x=-r. Каждый столбец J получают двумя полноценными полевыми решениями, изменяя одну из четырёх компонент токов.

Следовательно, построение центрального Якобиана требует восьми возмущённых решений. Плюс базовое поле, пробный шаг и заключительный повтор — именно поэтому число полевых обращений больше числа итераций.

\mathbf x=\begin{bmatrix}\operatorname{Re}\widehat I_1&\operatorname{Im}\widehat I_1&\operatorname{Re}\widehat I_2&\operatorname{Im}\widehat I_2\end{bmatrix}^{\mathsf T}
(Р6.1)
Токи задаются амплитудными фазорами.
\underline r_1=R_1\widehat I_1+\widehat U_{1,\mathrm{поле}}-\widehat U_{\mathrm{ист}},\qquad \underline r_2=\widehat U_{2,\mathrm{поле}}+R_{\mathrm{наг}}\widehat I_2
(Р6.2)
Внешний R_{1} добавляется один раз: исходные обмоточные области этой модели имеют \sigma =0.
J_{jk}\approx\frac{r_j(\mathbf x+h\mathbf e_k)-r_j(\mathbf x-h\mathbf e_k)}{2h},\qquad \mathbf J\Delta\mathbf x=-\mathbf r,\qquad \mathbf x_{\mathrm{нов}}=\mathbf x+\alpha\Delta\mathbf x
(Р6.3)
Множитель 0<\alpha \le 1 выбирают так, чтобы остаток уменьшался.

Сходимость и независимые контроли

Шаг разности h выбирается достаточно большим относительно ошибки поля и достаточно малым для локальной линейности. Повтор с h/2 выполняется из той же начальной точки: запуск из уже найденного решения практически не проверяет Якобиан.

Для продолжения по напряжению решение при 100 В используется как начальная точка для 105 В. Это ускоряет поиск, но не освобождает от заключительной проверки цепных уравнений при новом напряжении.

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

\varepsilon_r=\frac{\|\mathbf r\|_\infty}{\max(U_{\mathrm{ист,амп}},1\,\mathrm{В})},\qquad P=\frac12\operatorname{Re}\!\left(\widehat U\widehat I^{*}\right)
(Р6.4)
Нормировка остатка имеет объявленный размерный масштаб.

Постановка и последовательность расчёта

Таблица Р6.1. Исходные данные
ВеличинаЗначениеЕдиница
Амплитуды напряжения100;105В
Внешнее первичное сопротивление0,8Ом
Нагрузка10Ом
Разностный шаг0,001;0,0005А
Критерий относительного остатка10 × 10⁻⁶
  1. Изучите обозначения портов и выпишите четыре вещественных неизвестных до запуска программы.
  2. Выполните reproduce.py и сопоставьте новое контрольное поле с принятой рабочей точкой.
  3. Полные программы, названные далее, запускайте через run-series.py из корня рабочей папки: выберите SERIES_NUMBER по списку README.html. Диспетчер сам устанавливает подпапку и контекст исходной программы; прежний рабочий каталог FEMM не важен.
  4. Для полного поиска запустите reference/native-circuits.py встроенным Python из каталога reference. Файл native-newton-validation.py повторяет проверку шага из исходной точки.
  5. Просмотрите журнал r(k) и число обращений FEMM. Сравните 100 и 105 В; не заменяйте источник амплитудными 120 В из другого разбора.
  6. Сверьте заключительное поле, половинный шаг и строгую точность поля 10⁻¹⁰.

Численные результаты и их проверка

Показаны сохранённые результаты FEMM 5. Переключатель и курсор исследуют готовую серию; повторный расчёт выполняется программой в установленном FEMM.

Таблица Р6.2. Численные данные
Напряжение, амплитуда, ВRe I₁, АIm I₁, АRe I₂, АIm I₂, АОбращения FEMMНорма остатка ∞, В
1002,657-5,276-4,873-0,19674250,89 × 10⁻⁶
1052,99-7,556-5,102-0,28653256,17 × 10⁻⁶

Контрольные сопоставления

Как воспроизвести в FEMM 5

Скачайте программу и необходимые модели, данные и инструкции из списка ниже. Сохраните указанные имена файлов и структуру подпапок. Откройте программу в редакторе встроенного Python FEMM 5 и запустите её. Внешний CPython для управления FEMM не требуется.

Скачайте все нужные файлы со страницы исследования в новую папку, сохранив имена и структуру подпапок. Выполните reproduce.py встроенным Python FEMM 5: проверяются хеши, заново решаются контрольные поля и пересчитываются численные проверки таблиц. Рабочий каталог определяется по пути сохранённого сценария. Полные программы запускайте через run-series.py; номер SERIES_NUMBER выбирается по списку в README.html.

Программа встроенного Python

Текст программы
# Выполняется только встроенным Python FEMM 5.
# Распакуйте весь ZIP в новую папку и запустите этот файл через меню Scripts.
# Запускайте сохранённый файл целиком: путь передаётся в sys.argv[0].
import os, sys, json, math, hashlib
assert sys.implementation.name == 'micropython', 'Нужен встроенный Python FEMM 5'
_saved_file = globals().get('__file__') or (sys.argv[0] if sys.argv else '')
assert _saved_file and os.path.isfile(_saved_file), 'Запустите сохранённый файл из полного распакованного комплекта'
os.chdir(os.path.dirname(os.path.abspath(_saved_file)))
with open('manifest.json') as f: manifest = json.load(f)

def digest(path):
    h = hashlib.sha256()
    with open(path, 'rb') as f:
        while True:
            b = f.read(65536)
            if not b: break
            h.update(b)
    return ''.join('%02x' % b for b in h.digest())

# До расчёта проверяем переносимость и неизменность всех приложенных входов.
for item in manifest['files']:
    assert digest(item['path']) == item['sha256'], 'Изменён вход ' + item['path']
with open('expected-native.json') as f: cases = json.load(f)
os.mkdir('reproduced')   # Отказ от перезаписи предыдущего результата намеренный.
out = {'implementation': str(sys.implementation), 'runtime': sys.version,
       'case_id': manifest['case_id'], 'fresh_fields': [], 'checks': []}

def z(v):
    if isinstance(v, dict): return complex(v['re'], v['im'])
    if isinstance(v, (list, tuple)): return complex(v[0], v[1])
    return complex(v)

def serial(v):
    if isinstance(v, complex): return [v.real, v.imag]
    if isinstance(v, (list, tuple)): return [serial(t) for t in v]
    if isinstance(v, dict): return {k: serial(t) for k,t in v.items()}
    return v

def check(label, actual, expected, tolerance):
    error = abs(z(actual) - z(expected))
    out['checks'].append({'name':label,'actual':serial(actual),
                         'expected':serial(expected),'error':error,
                         'tolerance':tolerance,'passed':error <= tolerance})

for j, spec in enumerate(cases):
    for cmd in ('mo_close', 'mi_close'):
        try: globals()[cmd]()
        except Exception: pass
    opendocument(spec['model'])
    mi_saveas('reproduced/control-' + str(j) + '.fem')
    if 'frequency_Hz' in spec:
        mi_probdef(spec['frequency_Hz'], 'millimeters', 'planar', 1e-8,
                   spec.get('depth_mm',100), 20)
    for name,value in spec.get('set_circuits',{}).items():
        mi_modifycircprop(name, 1, z(value))
    if spec.get('rotate_deg'):
        mi_selectgroup(2); mi_moverotate(0,0,spec['rotate_deg'],4); mi_clearselected()
    if 'boundary' in spec:
        name,index,value = spec['boundary']; mi_modifyboundprop(name,index,value)
    mi_saveas('reproduced/control-' + str(j) + '.fem')
    mi_analyze(1); mi_loadsolution(); mo_smooth('off')
    row = {'model':spec['model'],'elements':mo_numelements(),'nodes':mo_numnodes(),
           'problem':serial(mo_getprobleminfo()),'circuits':{}}
    for name, values in spec.get('circuits',{}).items():
        actual = list(mo_getcircuitproperties(name)); row['circuits'][name] = serial(actual)
        for k,expected in enumerate(values):
            check('Поле '+str(j)+' '+name+' '+('I','U','Psi')[k], actual[k], expected,
                  max(1e-9, 1e-5*abs(z(expected))))
    mo_clearblock(); mo_groupselectblock()
    for key,expected in spec.get('integrals',{}).items():
        actual=mo_blockintegral(int(key)); row['integral_'+key]=serial(actual)
        check('Объёмный интеграл '+key,actual,expected,max(1e-7,1e-5*abs(z(expected))))
    mo_clearblock()
    if 'gap_torque' in spec:
        q=spec['gap_torque']; actual=mo_gapintegral(q['name'],0)
        if isinstance(actual,(list,tuple)): actual=actual[0]
        row['gap_torque_Nm']=serial(actual)
        check('Момент границы AGE',actual,q['expected'],q['tolerance'])
    if 'torque' in spec:
        q=spec['torque'];mo_groupselectblock(q['group']);actual=mo_blockintegral(22);mo_clearblock()
        row['torque_Nm']=serial(actual);check('Момент WST',actual,q['expected'],q['tolerance'])
    if 'magnet_probe' in spec:
        # Исходные графические минимумы получены со сглаживанием поля.
        # Настройку сохраняем явно, чтобы сравнивать одинаковые величины.
        mo_smooth('on')
        q=spec['magnet_probe']
        for p in q['extrema']:
            x,y=p['minimum_point_mm'];angle=math.atan2(y,x)
            sign=1 if p['pole']%2==0 else -1;v=mo_getpointvalues(x,y)
            h=(z(v[1]).real*sign*math.cos(angle)+z(v[2]).real*sign*math.sin(angle)-q['Br_T'])/(4*math.pi*1e-7*q['mu'])
            check('Физический H вдоль намагничивания полюса '+str(p['pole']),h,p['Hparallel_min_A_m'],.02)
    out['fresh_fields'].append(row)
    with open('reproduced/native-result.json','w') as f: json.dump(out,f)
    mo_close();mi_close()

# Повторяем предметные арифметические проверки сохранённых серий.
with open('postprocess.py') as f: exec(f.read(),globals())
out['passed']=all(q['passed'] for q in out['checks'])
with open('reproduced/native-result.json','w') as f: json.dump(out,f)
print('AUDIT_CASE',out['case_id'],'FIELDS',len(out['fresh_fields']),
      'CHECKS',len(out['checks']),'PASSED',out['passed'])
assert out['passed'], 'Есть несогласованные численные проверки; см. native-result.json'

Как трактовать результат

Расчётная среда: FEMM 5.0 06Oct2026; встроенный MicroPython; без внешнего Python/COM/TCP

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

Источники и дальнейшее чтение

  1. [15]
  2. [33]

← Все расчётные исследования