"""Полный цикл оптимизации внутри FEMM 5 / MicroPython. Запуск: FEMM → Run Script → этот сохранённый файл после распаковки комплекта. Здесь нет pyFEMM, внешнего контроллера, CPython-мостов или сторонних пакетов. Результаты пишутся в новый каталог results-run-NN; старые не затираются. """ import os import sys import json import math script = globals().get('__file__') or (sys.argv[0] if sys.argv else '') if not script or not os.path.isfile(script): raise RuntimeError('Запустите сохранённый run.py командой Run Script в FEMM 5.') os.chdir(os.path.dirname(os.path.abspath(script))) exec(open('model.py').read()) index = 1 while os.path.exists('results-run-%02d' % index): index += 1 folder = 'results-run-%02d' % index os.mkdir(folder) evaluations, cache = [], {} FORCE_SEARCH = 50.5 # Н: запас 1 % к обязательным 50 Н, а не подмена ограничения. result = {'schema': 'femm5-actuator-optimization-v1', 'parameters': PARAMETERS, 'search_force_N': FORCE_SEARCH, 'BH_T_A_per_m': BH, 'evaluations': evaluations, 'search_1d': {}, 'search_2d': {}, 'verification': [], 'algorithm_log': [], 'scope': 'Локально проверенный минимум медных потерь в принятом семействе при фиксированных I и N; не глобальное доказательство и не тепловой проект аппарата.'} def save(): with open(os.path.join(folder, 'results.json'), 'w') as stream: json.dump(result, stream) headers = ['id', 'stage', 'r_mm', 't_mm', 'force_N', 'power_copper_W', 'J_A_mm2', 'B_sample_max_T', 'mass_copper_kg', 'mass_steel_kg', 'feasible', 'search_feasible', 'mesh_scale', 'outer_radius_mm', 'elements', 'energy_identity_relative'] with open(os.path.join(folder, 'evaluations.csv'), 'w') as stream: stream.write(','.join(headers)+'\n') for row in evaluations: stream.write(','.join(str(row[h]) for h in headers)+'\n') def evaluate(r, t, stage, mesh=1.0, outer=100.0, gap=None): key = (round(r, 6), round(t, 6), mesh, outer, gap) if key in cache: row = cache[key] result['algorithm_log'].append({'stage': stage, 'id': row['id'], 'cached': True}) return row stem = 'field-%03d' % len(evaluations) row = field(r, t, folder, stem, mesh, outer, gap) row['id'], row['stage'] = len(evaluations), stage row['search_feasible'] = row['feasible'] and row['force_N'] >= FORCE_SEARCH evaluations.append(row) cache[key] = row result['algorithm_log'].append({'stage': stage, 'id': row['id'], 'cached': False}) save() print('%03d %s r=%.4f t=%.4f F=%.4f P=%.5f J=%.4f B=%.4f %s' % (row['id'], stage, r, t, row['force_N'], row['power_copper_W'], row['J_A_mm2'], row['B_sample_max_T'], row['search_feasible'])) return row def lower_feasible_radius(t, stage, iterations=8, bracket=None, mesh=1.0): """Восстановление допустимости по r. Применимо после проверки монотонности F. При фиксированном t медные потери строго возрастают по r аналитически. Поэтому ищется левая допустимая граница, а не произвольный минимум штрафа. """ lo, hi = PARAMETERS['r_bounds_mm'] if bracket is None else bracket low = evaluate(lo, t, stage, mesh) high = evaluate(hi, t, stage, mesh) # Сила служит монотонным ограничением; остальные ограничения проверяются явно. if low['force_N'] >= FORCE_SEARCH: return low if low['search_feasible'] else None if high['force_N'] < FORCE_SEARCH: return None for _ in range(iterations): mid = (lo+hi)/2 row = evaluate(mid, t, stage, mesh) if row['force_N'] >= FORCE_SEARCH: hi, high = mid, row else: lo = mid return high if high['search_feasible'] else None # 0. Полностью допустимый исходный проект для честного сравнения. baseline = evaluate(11.0, 4.0, 'baseline') assert baseline['search_feasible'] result['baseline_id'] = baseline['id'] # 1. Независимый одномерный перебор, затем дихотомия границы при t=4 мм. grid1 = [evaluate(6.0+0.5*k, 4.0, 'grid_1d') for k in range(15)] assert all(b['force_N'] > a['force_N'] for a, b in zip(grid1[:-1], grid1[1:])) one = lower_feasible_radius(4.0, 'bisection_1d', 9) assert one is not None grid1_best = min([x for x in grid1 if x['search_feasible']], key=lambda x: x['power_copper_W']) result['search_1d'] = {'fixed_t_mm': 4.0, 'grid_ids': [x['id'] for x in grid1], 'grid_best_id': grid1_best['id'], 'selected_id': one['id'], 'radius_bracket_width_mm': 7.0/512, 'monotonic_force_checked': True} # 2. Независимая сетка двух параметров. Никакая поверхность не подменяет FEMM. grid2 = [evaluate(r, t, 'grid_2d') for t in [3.0, 4.0, 5.0, 6.0, 7.0] for r in [8.0, 9.0, 10.0, 11.0, 12.0, 13.0]] best_grid2 = min([x for x in grid2 if x['search_feasible']], key=lambda x: x['power_copper_W']) # Ограниченный координатный поиск по толщине. На каждом шаге радиус вновь # определяется из НОВЫХ полей дихотомией, восстанавливающей допустимость. # Это двумерная задача с аналитически исключаемой монотонной переменной. best = one path = [{'t_mm': best['t_mm'], 'r_mm': best['r_mm'], 'id': best['id'], 'step_mm': 1.0}] step = 1.0 for iteration in range(8): candidates = [best] for sign in [-1, 1]: t = max(3.0, min(7.0, best['t_mm']+sign*step)) if abs(t-best['t_mm']) > 1e-9: item = lower_feasible_radius(t, 'coordinate_2d', 8) if item is not None: candidates.append(item) selected = min(candidates, key=lambda x: x['power_copper_W']) improved = selected['power_copper_W'] < best['power_copper_W']-0.0002 best = selected path.append({'t_mm': best['t_mm'], 'r_mm': best['r_mm'], 'id': best['id'], 'step_mm': step, 'improved': improved}) if not improved: step *= 0.5 if step < 0.125: break result['search_2d'] = {'grid_ids': [x['id'] for x in grid2], 'grid_best_id': best_grid2['id'], 'path': path, 'selected_id': best['id'], 'thickness_resolution_mm': step, 'radius_resolution_mm': 7.0/256, 'algorithm': 'Ограниченный координатный поиск с восстановлением допустимости по радиусу'} # Локальная доводка на более мелкой сетке. Три близкие толщины и узкая # дихотомия радиуса отделяют разрешение поисковой сетки от полевого решения. polish = [] for dt in [-0.1, 0.0, 0.1]: t = best['t_mm']+dt if 3 <= t <= 7: centre = best['r_mm']-0.55*dt item = lower_feasible_radius(t, 'boundary_polish', 5, [centre-0.15, centre+0.15], 0.6) if item is not None: polish.append(item) assert polish best = min(polish, key=lambda row: row['power_copper_W']) result['search_2d']['polish_ids'] = [row['id'] for row in polish] result['search_2d']['polish_radius_resolution_mm'] = 0.3/32 # Вторая независимая стадия: прямой поиск в локальной двумерной окрестности. # Если проверка находит лучший допустимый сосед, он принимается, а не скрывается. local_path = [best['id']] for local_iteration in range(4): neighbours = [best] for dr, dt in [(-0.1, 0), (0.1, 0), (0, -0.1), (0, 0.1), (-0.05, 0.1), (0.05, -0.1)]: r, t = best['r_mm']+dr, best['t_mm']+dt if 6 <= r <= 13 and 3 <= t <= 7: item = evaluate(r, t, 'local_pattern', 0.6) if item['search_feasible']: neighbours.append(item) selected = min(neighbours, key=lambda row: row['power_copper_W']) local_path.append(selected['id']) if selected['id'] == best['id']: break best = selected result['search_2d']['local_pattern_ids'] = local_path result['search_2d']['selected_id'] = best['id'] result['selected_id'] = best['id'] result['selected'] = best # 3. Проверки не входят в выбор по исходной сетке: иной размер элементов, # иная внешняя область, окрестность размеров, производная коэнергии. for name, mesh, outer in [('fine_mesh', 0.6, 100.0), ('finer_mesh', 0.4, 100.0), ('larger_domain', 0.6, 160.0)]: row = evaluate(best['r_mm'], best['t_mm'], name, mesh, outer) result['verification'].append({'kind': name, 'id': row['id'], 'force_relative_to_search': row['force_N']/best['force_N']-1, 'power_relative_to_search': row['power_copper_W']/best['power_copper_W']-1, 'feasible': row['feasible']}) assert row['feasible'] for dr, dt in [(-0.1, 0), (0.1, 0), (0, -0.1), (0, 0.1), (-0.05, 0.1), (0.05, -0.1)]: if 6 <= best['r_mm']+dr <= 13 and 3 <= best['t_mm']+dt <= 7: row = evaluate(best['r_mm']+dr, best['t_mm']+dt, 'neighbour', 0.6) result['verification'].append({'kind': 'neighbour', 'id': row['id'], 'dr_mm': dr, 'dt_mm': dt, 'feasible': row['feasible'], 'search_feasible': row['search_feasible'], 'power_difference_W': row['power_copper_W']-best['power_copper_W']}) dg = 0.02 minus = evaluate(best['r_mm'], best['t_mm'], 'coenergy_gap_minus', 0.4, 100.0, PARAMETERS['gap_mm']-dg) plus = evaluate(best['r_mm'], best['t_mm'], 'coenergy_gap_plus', 0.4, 100.0, PARAMETERS['gap_mm']+dg) derivative_force = -(plus['coenergy_J']-minus['coenergy_J'])/(2*dg*1e-3) reference = evaluations[result['verification'][1]['id']] result['coenergy_check'] = {'dg_mm': dg, 'minus_id': minus['id'], 'plus_id': plus['id'], 'attraction_from_coenergy_N': derivative_force, 'attraction_WST_N': reference['force_N'], 'relative_difference': derivative_force/reference['force_N']-1} assert abs(result['coenergy_check']['relative_difference']) < 0.02 assert all(x['energy_identity_relative'] < 0.005 for x in evaluations) result['improvement_power_percent'] = 100*(1-best['power_copper_W']/baseline['power_copper_W']) result['field_evaluations'] = len(evaluations) result['passed'] = True save() print('DONE', folder, 'fields', len(evaluations), 'power decrease %', result['improvement_power_percent'])