FEMM 5

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

Нелинейная IPM и MTPA: выбор тока по среднему моменту

У машины с внутренними магнитами поток возбуждения сочетается с магнитной анизотропией ротора. Момент возникает как от взаимодействия тока с полем магнитов, так и от различия магнитных сопротивлений вдоль осей d и q.

MTPA означает максимальный момент при заданной амплитуде вектора тока. В учебной IPM сильное шунтирование магнитов перемычками даёт небольшой момент; это особенность конкретной геометрии, а не показатель качества промышленных IPM. Оптимизировать нужно средний момент по периоду, а не случайный максимум в одном положении.

Место в курсе.Геометрию, физическую постановку и результаты можно изучать вместе с основным занятием. Серийную программу и обработку разбирайте после введения во встроенный 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\). При постоянных магнитах и взаимном возбуждении сначала учитывают ненулевое смещение потокосцепления.
  • Положительные направления, оси и специальные индексы определяются в конкретной задаче. Единицу берите из условия или заголовка таблицы; перед подстановкой в формулу СИ переводите миллиметры в метры, миллигенри в генри и градусы в радианы.
  • Численные результаты обычно показаны с четырьмя значащими цифрами; для сравнения близких значений в таблицах сходимости сохранено больше цифр. Полные значения остаются в исходных расчётных данных. Малые остатки вычисления координат на осях показаны как ноль; малые физические величины сохраняются в показательной записи.

От фазных величин к осям d–q

Амплитудно-инвариантное преобразование Парка задаёт токи d и q проекциями трёх фаз на вращающиеся оси. При симметрии iA+iB+iC=0 оно сохраняет величину iд²+iк²=2/3Σiф². В этом комплекте физическая ось d имеет смещение 90^\circ относительно программного отсчёта; оно явно сохранено.

Электромагнитный момент следует из производной коэнергии по положению ротора при фиксированных фазных токах. После преобразования в установившуюся синхронную систему получается выражение через Ψд, Ψк, iд, iк. В линейном приближении подстановка Ψд=Ψм+Lд iд и Ψк=Lк iк разделяет магнитный и реактивный вклады.

При насыщении постоянные Lд и Lк теряют универсальность, поэтому в этом исследовании момент берётся непосредственно из FEMM и усредняется по одинаковым механическим положениям. Формула d–q остаётся независимым согласующим контролем.

i_d=\frac23\sum_{j=0}^{2}i_j\cos(\theta_{\mathrm{эл}}-2\pi j/3),\qquad i_q=-\frac23\sum_{j=0}^{2}i_j\sin(\theta_{\mathrm{эл}}-2\pi j/3)
(Р10.1)
Латинские d и q — общепринятые обозначения осей.
T=\frac32p(\Psi_d i_q-\Psi_q i_d)=\frac32p\left[\Psi_{\mathrm{м}}i_q+(L_d-L_q)i_d i_q\right]
(Р10.2)
Второе равенство справедливо для линейной разделённой модели потокосцеплений.

Как проводится поиск максимума

На окружности iд=Icosβ, iк=Isinβ перебирается \beta. Для каждого \beta ротор и фазные токи согласованно проходят механический период 30^\circ, обусловленный симметрией четырёхполюсной 12-пазовой машины. Затем вычисляют средний момент и размах пульсаций.

После грубого поиска оптимум уточняют по \beta, проверяют 32 угловыми положениями, отрицательной половиной окружности и отдельным сгущением сетки. Полный переносимый план содержит 1888 полей IPM: 1056 в исходной и уточнённой сетках токовых углов, 736 в отрицательной половине окружности и 96 сеточных контролей. Это отдельные полевые задачи, не один статический расчёт.

При 5 А лучший найденный \beta равен 100^\circ, при 10 А — 115^\circ. Возрастание отклонения от 90^\circ показывает изменение оптимального соотношения d–q, но не доказывает глобальную оптимальность при любых токах, напряжениях и температурах.

\overline T(I,\beta)=\frac1N\sum_{j=0}^{N-1}T(I,\beta,\theta_j),\qquad \beta_{\mathrm{опт}}=\arg\max_\beta\overline T(I,\beta)
(Р10.3)
Критерий выбран для фиксированного I; ограничение напряжения и ослабление поля рассматриваются отдельно.

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

Таблица Р10.1. Исходные данные
ВеличинаЗначениеЕдиница
Амплитуды тока5;10А
Пары полюсов2
Пазы12
Механический период нагруженного момента30{}^\circ
Положение физической d-оси90{}^\circ
  1. Сначала разберите соглашение об осях: ошибочный сдвиг на 90^\circ поменяет смысл MTPA.
  2. Запустите reproduce.py; он заново решит точку оптимального \beta и проверит токовую окружность и потокосцепления.
  3. Полные программы, названные далее, запускайте через run-series.py из корня рабочей папки: выберите SERIES_NUMBER по списку README.html. Диспетчер сам устанавливает подпапку и контекст исходной программы; прежний рабочий каталог FEMM не важен.
  4. В reference/pm-data.json сравните средние моменты для всей сетки \beta и двух токов, затем исходные строки fine_angular_checks.
  5. Для полной серии запускайте reference/rebuild-complete-series.py последовательно до завершения контрольного плана; один запуск делает не более 40 полевых задач. Подходящая оболочка run-complete-batch.py устанавливает рабочую папку.
  6. Сравните средний момент WST и d–q; отдельно сравните размах пульсаций и уточнение сетки.

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

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

Таблица Р10.2. Численные данные
Амплитуда тока, Аβ от физической d-оси, °Средний момент WST, Н·мРазмах, Н·м
51000,023350,01283
5900,022950,01273
101150,050840,03025
10900,046810,0282

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

Как воспроизвести в 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. [33]
  2. [59]

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