VedEf Design Group for 5-parameter operating window exploration.

This is an application-specific OpenMDAO Group for the VedEf operating window example. It integrates all models and constraints for the 5-parameter exploration.

Note: This is an example integration, not a core library component.

import openmdao.api as om

from paroto.models.arc_density import ArcDensityConstraintModel
from paroto.models.arc_mobility import BrownianArcMobilityModel
from paroto.models.breakdown_voltage import RepetitivePulseBreakdownModel
from paroto.models.coverage import ArcCoverageModel
from paroto.models.initial_temperature import PulseInitialTemperatureModel
from paroto.models.thermal_diameter import ThermalDiameterDiffusionModel
from paroto.systems.generator import HighVoltageGeneratorSystem
from paroto.systems.torch.flow import TorchFlowSystem


class VedEfDesignGroup(om.Group):
    """VedEf operating window exploration group with 5 parameters and 10 constraints.

    This group integrates:
    - 5 design parameters: V, e, d_t, f, E_p
    - 10 constraints: Breakdown, HV window, Mechanical, Flow, Min/Max velocity,
      Mobility, Coverage, Power, Arc density

    Parameters (ARC002 naming convention)
    ----------
    G_UMAX_OUT : float (kV)
        Generator voltage (High voltage)
    G_e : float (mm)
        Interelectrode gap distance
    TP_D_OUT : float (mm)
        Torch outlet diameter
    G_F : float (kHz)
        PRF - Pulse Frequency
    G_Ep : float (J)
        Energy per pulse
    TP_QM : float (kg/h)
        Mass flow rate per torch
    """

    def initialize(self):
        """Declare options for model configuration."""
        self.options.declare(
            "target_power",
            default=25000.0,
            desc="Target total power (W)",
            recordable=False,
        )
        self.options.declare(
            "max_energy_density",
            default=1.0e9,
            desc="Maximum arc energy density (J/m³)",
            recordable=False,
        )

    def setup(self):
        """Set up the VedEf design group structure."""
        # Add Thermal Diameter Model (based on diffusion physics)
        self.add_subsystem(
            "ThermalDiameterModel",
            ThermalDiameterDiffusionModel(),
            promotes_inputs=[("pulse_duration", "pulse_duration"), ("gap_distance", "G_e")],
            promotes_outputs=["thermal_diameter"],
        )

        # Add Initial Temperature Model
        self.add_subsystem(
            "InitialTemperatureModel",
            PulseInitialTemperatureModel(),
            promotes_inputs=[
                ("energy_per_pulse", "G_Ep"),
                ("pulse_frequency", "G_F"),
                "thermal_diameter",
                "pulse_duration",
            ],
            promotes_outputs=["T_0", "T_max"],
        )

        # Connect arc_length from ThermalDiameterModel if it outputs it, otherwise set as input
        # For now, set arc_length as a simple geometric relationship to G_e
        self.add_subsystem(
            "ArcGeometryModel",
            om.ExecComp(
                "arc_length = G_e * expansion_factor",
                G_e={"units": "m"},
                arc_length={"units": "m"},
                expansion_factor={"val": 1.1},  # 10% expansion
            ),
            promotes_inputs=["G_e"],
            promotes_outputs=["arc_length"],
        )

        # Add Breakdown Model
        self.add_subsystem(
            "BreakdownModel",
            RepetitivePulseBreakdownModel(
                model_parameters={
                    "A_paschen": 24.4,
                    "B_paschen": 6.73,
                    "field_enhancement_factor": 1.0,
                    "ionization_time": 1e-7,
                    "pulse_amplitude_factor": 4.0,
                    "decay_constant": 1.0,
                }
            ),
            promotes_inputs=[
                ("pulse_frequency", "G_F"),
                ("pulse_duration", "pulse_duration"),
                ("gap_distance", "G_e"),
            ],
            promotes_outputs=["breakdown_voltage"],
        )

        # Connect T_0 to BreakdownModel
        self.connect("T_0", "BreakdownModel.preheat_temperature")

        # Add Coverage Model
        self.add_subsystem(
            "CoverageModel",
            ArcCoverageModel(),
            promotes_inputs=[
                ("gap_distance", "G_e"),
                ("torch_diameter", "TP_D_OUT"),
                ("pulse_frequency", "G_F"),
                "thermal_diameter",
            ],
            promotes_outputs=[
                "coverage_fraction",
                "coverage_margin",
                "coverage_constraint_satisfied",
            ],
        )

        # Add Arc Density Constraint Model
        self.add_subsystem(
            "ArcDensityModel",
            ArcDensityConstraintModel(
                model_parameters={"max_energy_density": self.options["max_energy_density"]}
            ),
            promotes_inputs=[("energy_per_pulse", "G_Ep"), "thermal_diameter", "arc_length"],
            promotes_outputs=["energy_density", "density_margin", "density_constraint_satisfied"],
        )

        # Constraint 1: Breakdown constraint V >= V_breakdown
        self.add_subsystem(
            "BreakdownConstraint",
            om.ExecComp(
                "breakdown_margin = G_UMAX_OUT - breakdown_voltage",
                breakdown_margin={"units": "V"},
                G_UMAX_OUT={"units": "V"},
                breakdown_voltage={"units": "V"},
            ),
            promotes_inputs=["G_UMAX_OUT", "breakdown_voltage"],
            promotes_outputs=["breakdown_margin"],
        )

        # HV Generator System with Operating Window Constraint
        self.add_subsystem(
            "HighVoltageGeneratorSystem",
            HighVoltageGeneratorSystem(
                design_impedance=50.0,
                pulse_duration=1e-6,
                max_power=50000.0,  # 50 kW limit
            ),
            promotes_inputs=[("pulse_frequency", "G_F"), ("hv_voltage", "G_UMAX_OUT")],
            promotes_outputs=["operating_window_power", "hv_operating_window_satisfied"],
        )

        # Constraint 3: Mechanical gap minimum
        self.add_subsystem(
            "MechanicalGapConstraint",
            om.ExecComp(
                "gap_margin = G_e - min_gap",
                gap_margin={"units": "m"},
                G_e={"units": "m"},
                min_gap={"units": "m", "val": 0.001},  # 1 mm minimum
            ),
            promotes_inputs=["G_e"],
        )

        # Flow System (calculates velocity, residence_time, pressure_drop)
        self.add_subsystem(
            "TorchFlowSystem",
            TorchFlowSystem(gas_density=0.717),  # CH4 at 1 bar, 300K
            promotes_inputs=[
                ("mass_flow", "TP_QM"),
                ("torch_diameter", "TP_D_OUT"),
                ("geometry_params", "G_e"),
                ("torch_length", "arc_length"),
            ],
            promotes_outputs=["flow_velocity", "residence_time"],
        )

        # Connect flow velocity to CoverageModel
        self.connect("flow_velocity", "CoverageModel.flow_velocity")

        # Constraint 5: Minimum velocity
        self.add_subsystem(
            "MinVelocityConstraint",
            om.ExecComp(
                "min_velocity_margin = flow_velocity - min_vel",
                min_velocity_margin={"units": "m/s"},
                flow_velocity={"units": "m/s"},
                min_vel={"units": "m/s", "val": 0.1},
            ),
            promotes_inputs=["flow_velocity"],
        )

        # Constraint 6: Maximum velocity (Mach < 0.3)
        self.add_subsystem(
            "MaxVelocityConstraint",
            om.ExecComp(
                "max_velocity_margin = max_mach - flow_velocity / speed_of_sound",
                flow_velocity={"units": "m/s"},
                speed_of_sound={"units": "m/s", "val": 340.0},  # Approx for gas
                max_mach={"val": 0.3},
                max_velocity_margin={"val": 0.3},
            ),
            promotes_inputs=["flow_velocity"],
        )

        # Arc Mobility Model (replaces simple T_max/T_0 constraint)
        self.add_subsystem(
            "ArcMobilityModel",
            BrownianArcMobilityModel(),
            promotes_inputs=[("pulse_frequency", "G_F"), ("gap_distance", "G_e")],
        )

        # Connect temperature from InitialTemperatureModel
        self.connect("T_0", "ArcMobilityModel.preheat_temperature")

        # Constraint 7: Mobility (uses physics-based mobility factor)
        self.add_subsystem(
            "MobilityConstraint",
            om.ExecComp(
                "mobility_margin = T_max / (T_0 + 1e-10) - min_ratio",
                T_max={"units": "K"},
                T_0={"units": "K"},
                min_ratio={"val": 10.0},
                mobility_margin={"val": 0.0},
            ),
            promotes_inputs=["T_max", "T_0"],
        )

        # Constraint 9: Total Power P = f * E_p (using ARC002 naming)
        # First calculate total power: G_PW_OUT = G_F * G_Ep
        self.add_subsystem(
            "PowerCalculation",
            om.ExecComp(
                "G_PW_OUT = G_F * G_Ep",
                G_F={"units": "Hz"},
                G_Ep={"units": "J"},
                G_PW_OUT={"units": "W"},
            ),
            promotes_inputs=["G_F", "G_Ep"],
            promotes_outputs=["G_PW_OUT"],
        )

        # Then calculate power error
        self.add_subsystem(
            "PowerConstraint",
            om.ExecComp(
                "power_error = abs(G_PW_OUT - target_power) / target_power",
                G_PW_OUT={"units": "W"},
                target_power={"units": "W", "val": self.options["target_power"]},
                power_error={"val": 0.0},
            ),
            promotes_inputs=["G_PW_OUT"],
        )

        # Set input defaults (using ARC002 parameter tags)
        self.set_input_defaults("G_UMAX_OUT", val=20000.0, units="V")  # Generator voltage
        self.set_input_defaults("G_e", val=0.01, units="m")  # Interelectrode gap
        self.set_input_defaults("TP_D_OUT", val=0.02, units="m")  # Torch outlet diameter
        self.set_input_defaults("G_F", val=50000.0, units="Hz")  # Pulse frequency (PRF)
        self.set_input_defaults("G_Ep", val=10e-3, units="J")  # Energy per pulse
        self.set_input_defaults("pulse_duration", val=1e-6, units="s")
        # thermal_diameter now computed by thermal_diam component (diffusion-based)
        # arc_length now computed by arc_geom component
        self.set_input_defaults("TP_QM", val=160.0 / 86400.0, units="kg/s")  # Mass flow rate

        # Set inputs for ThermalDiameterModel component
        self.set_input_defaults("ThermalDiameterModel.gas_density", val=0.717, units="kg/m**3")
        self.set_input_defaults(
            "ThermalDiameterModel.thermal_conductivity", val=0.03, units="W/(m*K)"
        )
        self.set_input_defaults(
            "ThermalDiameterModel.gas_heat_capacity", val=2200.0, units="J/(kg*K)"
        )

        # Set gas density for InitialTemperatureModel
        self.set_input_defaults("InitialTemperatureModel.gas_density", val=0.717, units="kg/m**3")

        # Connect residence time from FlowSystem to BreakdownModel
        self.connect("residence_time", "BreakdownModel.residence_time")

        # Add gas properties for BreakdownModel
        self.set_input_defaults("BreakdownModel.gas_properties_pressure", val=101325.0, units="Pa")

        # ArcMobilityModel requires additional inputs
        self.set_input_defaults("ArcMobilityModel.sustainer_pulse_duration", val=1e-6, units="s")
        self.set_input_defaults("ArcMobilityModel.arc_current", val=100.0, units="A")