Note
Go to the end to download the full example code.
Search focusing on achieving temperature rise for mobility constraint.
import importlib.util
from pathlib import Path
import numpy as np
# Import setup_vedef_problem from the numbered file
spec = importlib.util.spec_from_file_location(
"problem_module", Path(__file__).parent / "1_setup_vedef_problem.py"
)
problem_module = importlib.util.module_from_spec(spec)
spec.loader.exec_module(problem_module)
setup_vedef_problem = problem_module.setup_vedef_problem
def evaluate_case(V, e, d_t, f, E_p):
"""Evaluate a parameter case."""
try:
prob = setup_vedef_problem(target_power=25000.0, max_energy_density=1e9, as_subsystem=False)
# Fixed parameters
prob.set_val("pulse_duration", 1e-6, units="s")
prob.set_val("TP_QM", 160.0 / 86400.0, units="kg/s")
prob.set_val("ThermalDiameterModel.gas_density", 0.717, units="kg/m**3")
prob.set_val("ThermalDiameterModel.gas_heat_capacity", 2200.0, units="J/(kg*K)")
prob.set_val("ThermalDiameterModel.thermal_conductivity", 0.08, units="W/(m*K)")
prob.set_val("InitialTemperatureModel.gas_density", 0.717, units="kg/m**3")
prob.set_val("BreakdownModel.gas_properties_pressure", 101325.0, units="Pa")
# Design variables (using ARC002 naming convention)
prob.set_val("G_UMAX_OUT", V, units="V")
prob.set_val("G_e", e, units="m")
prob.set_val("TP_D_OUT", d_t, units="m")
prob.set_val("G_F", f, units="Hz")
prob.set_val("G_Ep", E_p, units="J")
prob.run_model()
results = {
"V": V,
"e": e,
"d_t": d_t,
"f": f,
"E_p": E_p,
"T_0": prob.get_val("T_0", units="K")[0],
"T_max": prob.get_val("T_max", units="K")[0],
"breakdown_margin": prob.get_val("breakdown_margin", units="V")[0],
"coverage_margin": prob.get_val("coverage_margin")[0],
"density_margin": prob.get_val("density_margin", units="J/m**3")[0],
"power_error": prob.get_val("PowerConstraint.power_error")[0],
}
# Check if all values are finite
if not all(
np.isfinite(v) for k, v in results.items() if k not in ["V", "e", "d_t", "f", "E_p"]
):
return None
return results
except Exception:
return None
print("=" * 80)
print("Temperature-Focused Search")
print("=" * 80)
print("\nStrategy: Low frequency + high energy per pulse to maximize T_max/T_0\n")
# For 25 kW target, use lower frequency with higher energy per pulse
# This should give more time for temperature rise
param_sets = [
# (f_kHz, E_p_mJ) pairs that give 25 kW
(10, 2500), # Very low freq, very high energy
(12.5, 2000),
(16.67, 1500),
(20, 1250),
(25, 1000),
(33.33, 750),
(50, 500),
]
# Geometry
V_range = [10e3, 15e3, 20e3, 25e3, 30e3]
e_range = [0.005, 0.0075, 0.01, 0.0125, 0.015]
d_t_range = [0.015, 0.02, 0.025, 0.03]
print(f"Testing {len(param_sets)} frequency/energy combinations")
print(f" x {len(V_range)} voltages")
print(f" x {len(e_range)} gaps")
print(f" x {len(d_t_range)} torch diameters")
print(f"Total: {len(param_sets) * len(V_range) * len(e_range) * len(d_t_range)} evaluations\n")
best_ratio = 0
best_case = None
feasible_cases = []
for idx, (f_kHz, E_p_mJ) in enumerate(param_sets):
f = f_kHz * 1e3
E_p = E_p_mJ * 1e-3
print(f"\n[{idx + 1}/{len(param_sets)}] f={f_kHz:.1f}kHz, E_p={E_p_mJ:.0f}mJ")
print("-" * 80)
local_best_ratio = 0
for V in V_range:
for e in e_range:
for d_t in d_t_range:
results = evaluate_case(V, e, d_t, f, E_p)
if results is None:
continue
ratio = results["T_max"] / (results["T_0"] + 1e-10)
# Track best ratio overall
if ratio > best_ratio:
best_ratio = ratio
best_case = results
# Track local best
if ratio > local_best_ratio:
local_best_ratio = ratio
# Check feasibility
feasible = (
results["breakdown_margin"] > 0
and results["coverage_margin"] > 0
and results["density_margin"] > 0
and results["power_error"] < 0.1
and ratio > 10.0
)
if feasible:
print(
f" [OK] FEASIBLE! V={V / 1e3:.0f}kV, e={e * 1e3:.1f}mm, d_t={d_t * 1e3:.0f}mm"
)
print(
f" T_0={results['T_0']:.0f}K, T_max={results['T_max']:.0f}K, ratio={ratio:.1f}"
)
feasible_cases.append(results)
print(f" Best T_max/T_0 for this (f,E_p): {local_best_ratio:.2f}")
print("\n" + "=" * 80)
print("RESULTS")
print("=" * 80)
if feasible_cases:
print(f"\n[SUCCESS] Found {len(feasible_cases)} feasible cases!\n")
print("First feasible case:")
case = feasible_cases[0]
ratio = case["T_max"] / (case["T_0"] + 1e-10)
print(f" V = {case['V'] / 1e3:.1f} kV")
print(f" e = {case['e'] * 1e3:.1f} mm")
print(f" d_t = {case['d_t'] * 1e3:.1f} mm")
print(f" f = {case['f'] / 1e3:.1f} kHz")
print(f" E_p = {case['E_p'] * 1e3:.1f} mJ")
print(f"\n T_0 = {case['T_0']:.1f} K")
print(f" T_max = {case['T_max']:.1f} K")
print(f" T_max/T_0 = {ratio:.2f}")
print(f" Power error = {case['power_error'] * 100:.1f}%")
print(f" Breakdown margin = {case['breakdown_margin'] / 1e3:.1f} kV")
print(f" Coverage margin = {case['coverage_margin']:.3f}")
else:
print("\n[NO FEASIBLE SOLUTION]")
print(f"Best T_max/T_0 ratio achieved: {best_ratio:.2f} (need > 10.0)\n")
if best_case:
print("Best case (highest T_max/T_0):")
print(f" V = {best_case['V'] / 1e3:.1f} kV")
print(f" e = {best_case['e'] * 1e3:.1f} mm")
print(f" d_t = {best_case['d_t'] * 1e3:.1f} mm")
print(f" f = {best_case['f'] / 1e3:.1f} kHz")
print(f" E_p = {best_case['E_p'] * 1e3:.1f} mJ")
print(f"\n T_0 = {best_case['T_0']:.1f} K")
print(f" T_max = {best_case['T_max']:.1f} K")
print(f" T_max/T_0 = {best_ratio:.2f}")
print(f" Power error = {best_case['power_error'] * 100:.1f}%")
print(f" Breakdown margin = {best_case['breakdown_margin'] / 1e3:.1f} kV")
print(f" Coverage margin = {best_case['coverage_margin']:.3f}")
print(f" Density margin = {best_case['density_margin'] / 1e6:.1f} MJ/m³")
print("\n NOTE: The mobility constraint (T_max/T_0 > 10) may be too strict")
print(" for these operating conditions. Consider reviewing the constraint")
print(" or adjusting model parameters.")
print("\n" + "=" * 80)