diff --git a/examples/single_model_evaluation.ex.py b/examples/single_model_evaluation.ex.py index f9e66c2ae..af37205f9 100644 --- a/examples/single_model_evaluation.ex.py +++ b/examples/single_model_evaluation.ex.py @@ -96,7 +96,9 @@ def run_impurities(w_imp_fracs): single_run.models.physics.run() # Evaluate constraint equation 15 (L-H threshold constraint) - con15_value = ConstraintManager.evaluate_constraint(15).normalised_residual + con15_value = ConstraintManager.evaluate_constraint( + 15, single_run.data + ).normalised_residual # Need to copy values p_plasma_rad_mw[i] = data_structure.physics_variables.p_plasma_rad_mw.item() diff --git a/process/core/caller.py b/process/core/caller.py index ae3f29cf7..80e45d899 100644 --- a/process/core/caller.py +++ b/process/core/caller.py @@ -18,6 +18,7 @@ from process.models.tfcoil.base import TFConductorModel if TYPE_CHECKING: + from process.core.model import DataStructure from process.main import Models logger = logging.getLogger(__name__) @@ -26,17 +27,21 @@ class Caller: """Calls physics and engineering models.""" - def __init__(self, models: Models): + def __init__(self, models: Models, data: DataStructure): """Initialise all physics and engineering models. To ensure that, at the start of a run, all physics/engineering variables are fully initialised with consistent values, the models are - called with the initial optimisation paramters, x. + called with the initial optimisation parameters, x. :param models: physics and engineering model objects :type models: Models + :param data: data structure object to be passed on to the constraint + evaluators + :type data: DataStructure """ self.models = models + self.data = data @staticmethod def check_agreement( @@ -92,7 +97,7 @@ def call_models(self, xc: np.ndarray, m: int) -> tuple[float, np.ndarray]: self._call_models_once(xc) # Evaluate objective function and constraints objf = objective_function(data_structure.numerics.minmax) - conf, _, _, _, _ = constraints.constraint_eqns(m, -1) + conf, _, _, _, _ = constraints.constraint_eqns(m, -1, self.data) if objf_prev is None and conf_prev is None: # First run: run again to check idempotence @@ -156,7 +161,7 @@ def call_models_and_write_output(self, xc: np.ndarray, ifail: int): OutputFileManager.open_idempotence_files() self._call_models_once(xc) # Write mfile - finalise(self.models, ifail) + finalise(self.models, self.data, ifail) # Extract data from intermediate idempotence-checking mfile mfile_path = ( @@ -198,7 +203,7 @@ def call_models_and_write_output(self, xc: np.ndarray, ifail: int): # now idempotence checking complete OutputFileManager.close_idempotence_files() # Write final output file and mfile - finalise(self.models, ifail) + finalise(self.models, self.data, ifail) return # Mfiles not yet idempotent: need to re-evaluate models @@ -225,6 +230,7 @@ def call_models_and_write_output(self, xc: np.ndarray, ifail: int): OutputFileManager.close_idempotence_files() finalise( self.models, + self.data, ifail, non_idempotent_msg=non_idempotent_warning + "\n" + non_idempotent_table, ) @@ -380,20 +386,24 @@ def _call_models_once(self, xc: np.ndarray): # FISPACT and LOCA model (not used)- removed -def write_output_files(models: Models, ifail: int, *, runtime: float | None = None): +def write_output_files( + models: Models, data: DataStructure, ifail: int, *, runtime: float | None = None +): """Evaluate models and write output files (OUT.DAT and MFILE.DAT). Parameters ---------- models : Models physics and engineering models + data: DataStructure + data structure object ifail : int solver return code """ n = data_structure.numerics.nvar x = data_structure.numerics.xcm[:n] # Call models, ensuring output mfiles are fully idempotent - caller = Caller(models) + caller = Caller(models, data) if runtime is not None: ovarre( constants.MFILE, diff --git a/process/core/final.py b/process/core/final.py index 52fc0adda..113ea58df 100644 --- a/process/core/final.py +++ b/process/core/final.py @@ -10,7 +10,7 @@ from process.data_structure import numerics -def finalise(models, ifail: int, non_idempotent_msg: str | None = None): +def finalise(models, data, ifail: int, non_idempotent_msg: str | None = None): """Routine to print out the final point in the scan. Writes to OUT.DAT and MFILE.DAT. @@ -19,6 +19,8 @@ def finalise(models, ifail: int, non_idempotent_msg: str | None = None): ---------- models : process.main.Models physics and engineering model objects + data: DataStructure + data structure object to provide data to evaluate the constraints ifail : int error flag non_idempotent_msg : None | str, optional @@ -31,7 +33,7 @@ def finalise(models, ifail: int, non_idempotent_msg: str | None = None): # Output relevant to no optimisation if numerics.ioptimz == -2: - output_evaluation() + output_evaluation(data) # Print non-idempotence warning to OUT.DAT only if non_idempotent_msg: @@ -42,8 +44,13 @@ def finalise(models, ifail: int, non_idempotent_msg: str | None = None): op.write(models, constants.NOUT) -def output_evaluation(): - """Write output for an evaluation run of PROCESS""" +def output_evaluation(data): + """Write output for an evaluation run of PROCESS + Parameters + ---------- + data: DataStructure + data structure object to provide data to evaluate the constraints + """ po.oheadr(constants.NOUT, "Numerics") po.ocmmnt(constants.NOUT, "PROCESS has performed an evaluation run.") po.oblnkl(constants.NOUT) @@ -55,7 +62,7 @@ def output_evaluation(): # Print the residuals of the constraint equations residual_error, value, residual, symbols, units = constraints.constraint_eqns( - numerics.neqns + numerics.nineqns, -1 + numerics.neqns + numerics.nineqns, -1, data ) labels = [ diff --git a/process/core/init.py b/process/core/init.py index 1990f8c7b..34050a5ff 100644 --- a/process/core/init.py +++ b/process/core/init.py @@ -22,7 +22,6 @@ from process.data_structure.ccfe_hcpb_module import init_ccfe_hcpb_module from process.data_structure.constraint_variables import init_constraint_variables from process.data_structure.cost_variables import init_cost_variables -from process.data_structure.cs_fatigue_variables import init_cs_fatigue_variables from process.data_structure.current_drive_variables import init_current_drive_variables from process.data_structure.dcll_variables import init_dcll_module from process.data_structure.divertor_variables import init_divertor_variables @@ -299,7 +298,6 @@ def init_all_module_vars(): init_pulse_variables() init_rebco_variables() init_reinke_variables() - init_cs_fatigue_variables() init_blanket_library() init_dcll_module() init_power_variables() diff --git a/process/core/input.py b/process/core/input.py index 2b0001aac..7bfc2dae9 100644 --- a/process/core/input.py +++ b/process/core/input.py @@ -255,9 +255,7 @@ def __post_init__(self): "bioshld_thk": InputVariable( data_structure.buildings_variables, float, range=(0.25, 25.0) ), - "bkt_life_csf": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.0, 1.0) - ), + "bkt_life_csf": InputVariable("cs_fatigue", float, range=(0.0, 1.0)), "blbmith": InputVariable(data_structure.build_variables, float, range=(0.0, 2.0)), "blbmoth": InputVariable(data_structure.build_variables, float, range=(0.0, 2.0)), "blbpith": InputVariable(data_structure.build_variables, float, range=(0.0, 2.0)), @@ -813,9 +811,7 @@ def __post_init__(self): "f_p_shld_coolant_pump_total_heat": InputVariable( data_structure.heat_transport_variables, float, range=(0.0, 0.2) ), - "fracture_toughness": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.1, 100000000.0) - ), + "fracture_toughness": InputVariable("cs_fatigue", float, range=(0.1, 100000000.0)), "fradpwr": InputVariable( data_structure.constraint_variables, float, range=(0.0, 1.0) ), @@ -1056,9 +1052,7 @@ def __post_init__(self): "mvalim": InputVariable( data_structure.constraint_variables, float, range=(0.0, 1000.0) ), - "n_cycle_min": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.0, 100000000.0) - ), + "n_cycle_min": InputVariable("cs_fatigue", float, range=(0.0, 100000000.0)), "n_tf_coils": InputVariable( data_structure.tfcoil_variables, float, range=(0.0, 100.0) ), @@ -1110,12 +1104,8 @@ def __post_init__(self): "p_he": InputVariable( data_structure.primary_pumping_variables, float, range=(0.0, 100000000.0) ), - "paris_coefficient": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1e-20, 10.0) - ), - "paris_power_law": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1.0, 10.0) - ), + "paris_coefficient": InputVariable("cs_fatigue", float, range=(1e-20, 10.0)), + "paris_power_law": InputVariable("cs_fatigue", float, range=(1.0, 10.0)), "pres_vv_chamber_base": InputVariable( data_structure.vacuum_variables, float, range=(1e-08, 0.001) ), @@ -1281,9 +1271,7 @@ def __post_init__(self): "redun_vacp": InputVariable( data_structure.cost_variables, float, range=(0.0, 100.0) ), - "residual_sig_hoop": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.0, 1000000000.0) - ), + "residual_sig_hoop": InputVariable("cs_fatigue", float, range=(0.0, 1000000000.0)), "rho_tf_bus": InputVariable( data_structure.tfcoil_variables, float, range=(0.0, 1e-05) ), @@ -1339,15 +1327,9 @@ def __post_init__(self): "sec_buildings_w": InputVariable( data_structure.buildings_variables, float, range=(10.0, 1000.0) ), - "sf_fast_fracture": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1.0, 10.0) - ), - "sf_radial_crack": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1.0, 10.0) - ), - "sf_vertical_crack": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1.0, 10.0) - ), + "sf_fast_fracture": InputVariable("cs_fatigue", float, range=(1.0, 10.0)), + "sf_radial_crack": InputVariable("cs_fatigue", float, range=(1.0, 10.0)), + "sf_vertical_crack": InputVariable("cs_fatigue", float, range=(1.0, 10.0)), "shdr": InputVariable(data_structure.ife_variables, float, range=(0.0, 10.0)), "shdzl": InputVariable(data_structure.ife_variables, float, range=(0.0, 10.0)), "shdzu": InputVariable(data_structure.ife_variables, float, range=(0.0, 10.0)), @@ -1413,12 +1395,8 @@ def __post_init__(self): "dx_tf_turn_cable_space_general": InputVariable( data_structure.tfcoil_variables, float, range=(0.0, 0.1) ), - "t_crack_radial": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1e-05, 1.0) - ), - "t_crack_vertical": InputVariable( - data_structure.cs_fatigue_variables, float, range=(1e-05, 1.0) - ), + "t_crack_radial": InputVariable("cs_fatigue", float, range=(1e-05, 1.0)), + "t_crack_vertical": InputVariable("cs_fatigue", float, range=(1e-05, 1.0)), "t_crit_nbti": InputVariable( data_structure.tfcoil_variables, float, range=(0.0, 15.0) ), @@ -1443,12 +1421,8 @@ def __post_init__(self): "t_plant_pulse_plasma_current_ramp_down": InputVariable( data_structure.times_variables, float, range=(0.0, 10000.0) ), - "dr_cs_turn_conduit": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.001, 1.0) - ), - "dz_cs_turn_conduit": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.001, 1.0) - ), + "dr_cs_turn_conduit": InputVariable("cs_fatigue", float, range=(0.001, 1.0)), + "dz_cs_turn_conduit": InputVariable("cs_fatigue", float, range=(0.001, 1.0)), "dx_tf_turn_general": InputVariable( data_structure.tfcoil_variables, float, range=(0.0, 0.1) ), @@ -1751,9 +1725,7 @@ def __post_init__(self): "pflux_fw_neutron_max_mw": InputVariable( data_structure.constraint_variables, float, range=(0.001, 50.0) ), - "walker_coefficient": InputVariable( - data_structure.cs_fatigue_variables, float, range=(0.1, 10.0) - ), + "walker_coefficient": InputVariable("cs_fatigue", float, range=(0.1, 10.0)), "wallpf": InputVariable(data_structure.fwbs_variables, float, range=(1.0, 2.0)), "warm_shop_h": InputVariable( data_structure.buildings_variables, float, range=(1.0, 100.0) diff --git a/process/core/io/plot/scans.py b/process/core/io/plot/scans.py index b3c5d30f0..c9a18c7d6 100644 --- a/process/core/io/plot/scans.py +++ b/process/core/io/plot/scans.py @@ -139,9 +139,9 @@ def plot_scan( 62: "j_cs_flat_top_end", 63: "dr_cs", 64: "f_z_cs_tf_internal", - 65: "csfv.n_cycle_min", - 66: "pfv.f_a_cs_turn_steel", - 67: "csfv.t_crack_vertical", + 65: "n_cycle_min", + 66: "f_a_cs_turn_steel", + 67: "t_crack_vertical", 68: "inlet_temp_liq", 69: "outlet_temp_liq", 70: "blpressure_liq", diff --git a/process/core/io/variable_metadata.py b/process/core/io/variable_metadata.py index 2e4cdbdf4..b0bf93117 100644 --- a/process/core/io/variable_metadata.py +++ b/process/core/io/variable_metadata.py @@ -299,17 +299,17 @@ class VariableMetadata: description="Minimum burn time", units="s", ), - "pfv.f_a_cs_turn_steel": VariableMetadata( + "f_a_cs_turn_steel": VariableMetadata( latex=r"$f_{\mathrm{Steel}}^{\mathrm{CS}}$", description="Steel fraction in CS coil", units="", ), - "csfv.dr_cs_turn_conduit": VariableMetadata( + "dr_cs_turn_conduit": VariableMetadata( latex=r"$Turn_{\mathrm{radial}}^{\mathrm{CS}}[$m$]$", description="Radial turn length", units="m", ), - "csfv.t_crack_vertical": VariableMetadata( + "t_crack_vertical": VariableMetadata( latex=r"$Crack_{\mathrm{vertical}}^{\mathrm{CS}}[$m$]$", description="Vertical crack length", units="m", diff --git a/process/core/model.py b/process/core/model.py index 5b6dc659b..d586f80ac 100644 --- a/process/core/model.py +++ b/process/core/model.py @@ -2,6 +2,7 @@ from dataclasses import dataclass, fields from process.data_structure.cost_2015_variables import Cost2015Data +from process.data_structure.cs_fatigue_variables import CSFatigueData from process.data_structure.water_usage_variables import WaterUseData initialise_later = object() @@ -11,6 +12,7 @@ class DataStructure: water_use: WaterUseData = initialise_later costs_2015: Cost2015Data = initialise_later + cs_fatigue: CSFatigueData = initialise_later def __post_init__(self): for f in fields(self): diff --git a/process/core/scan.py b/process/core/scan.py index 77631c4b0..815ca97eb 100644 --- a/process/core/scan.py +++ b/process/core/scan.py @@ -16,7 +16,6 @@ build_variables, constraint_variables, cost_variables, - cs_fatigue_variables, current_drive_variables, divertor_variables, fwbs_variables, @@ -195,17 +194,20 @@ def _missing_(cls, var): class Scan: """Perform a parameter scan using the Fortran scan module.""" - def __init__(self, models, solver): + def __init__(self, models, solver, data): """Immediately run the run_scan() method. :param models: physics and engineering model objects :type models: process.main.Models :param solver: which solver to use, as specified in solver.py :type solver: str + :param data: data structure object + :type data: DataStructure """ self.models = models self.solver = solver - self.solver_handler = SolverHandler(models, solver) + self.data = data + self.solver_handler = SolverHandler(models, solver, data) self.run_scan() def run_scan(self): @@ -223,7 +225,10 @@ def run_scan(self): start_time = time.time() ifail = self.doopt() write_output_files( - models=self.models, ifail=ifail, runtime=time.time() - start_time + models=self.models, + data=self.data, + ifail=ifail, + runtime=time.time() - start_time, ) show_errors(constants.NOUT) return @@ -521,7 +526,7 @@ def post_optimise(self, ifail: int): ) con1, con2, err, _, lab = constraints.constraint_eqns( - numerics.neqns + numerics.nineqns, -1 + numerics.neqns + numerics.nineqns, -1, self.data ) # Write equality constraints to mfile @@ -597,7 +602,7 @@ def post_optimise(self, ifail: int): for i in range(numerics.neqns, numerics.neqns + numerics.nineqns): name = numerics.lablcc[numerics.icc[i] - 1] constraint = constraints.ConstraintManager.evaluate_constraint( - int(numerics.icc[i]) + int(numerics.icc[i]), self.data ) inequality_constraint_table.append([ @@ -825,7 +830,10 @@ def scan_1d(self): ifail = self.doopt() scan_1d_ifail_dict[iscan] = ifail write_output_files( - models=self.models, ifail=ifail, runtime=time.time() - start_time + models=self.models, + data=self.data, + ifail=ifail, + runtime=time.time() - start_time, ) show_errors(constants.NOUT) @@ -882,7 +890,10 @@ def scan_2d(self): start_time = time.time() ifail = self.doopt() write_output_files( - models=self.models, ifail=ifail, runtime=time.time() - start_time + models=self.models, + data=self.data, + ifail=ifail, + runtime=time.time() - start_time, ) show_errors(constants.NOUT) @@ -1171,11 +1182,11 @@ def scan_select(self, nwp, swp, iscn): case 64: pfcoil_variables.ohhghf = swp[iscn - 1] case 65: - cs_fatigue_variables.n_cycle_min = swp[iscn - 1] + self.data.cs_fatigue.n_cycle_min = swp[iscn - 1] case 66: pfcoil_variables.oh_steel_frac = swp[iscn - 1] case 67: - cs_fatigue_variables.t_crack_vertical = swp[iscn - 1] + self.data.cs_fatigue.t_crack_vertical = swp[iscn - 1] case 68: fwbs_variables.inlet_temp_liq = swp[iscn - 1] case 69: diff --git a/process/core/solver/constraints.py b/process/core/solver/constraints.py index 6f757a4a2..8b5927326 100644 --- a/process/core/solver/constraints.py +++ b/process/core/solver/constraints.py @@ -8,6 +8,7 @@ from process import data_structure from process.core import constants from process.core.exceptions import ProcessError, ProcessValueError +from process.core.model import DataStructure from process.models.physics.physics import BetaComponentLimits from process.models.tfcoil.base import TFConductorModel @@ -126,13 +127,16 @@ def get_constraint(cls, name: Hashable): return cls._constraint_registry.get(name) @classmethod - def evaluate_constraint(cls, name: Hashable): + def evaluate_constraint(cls, name: Hashable, data: DataStructure): """Evalutes a constraint with a given name. Parameters ---------- name : Hashable the name of the constraint + data: DataStructure + data structure object for providing constraint data to evaluate + constraint equations Returns ------- @@ -145,7 +149,7 @@ def evaluate_constraint(cls, name: Hashable): error_msg = f"Constraint '{name}' cannot be found." raise ProcessError(error_msg) - return registration.constraint_equation(registration) + return registration.constraint_equation(registration, data) def leq(value: float, bound: float, registration: ConstraintRegistration): @@ -188,7 +192,7 @@ def eq(value: float, bound: float, registration: ConstraintRegistration): @ConstraintManager.register_constraint(1, "", "=") -def constraint_equation_1(constraint_registration): +def constraint_equation_1(constraint_registration, _data): """Relationship between beta, temperature (keV) and density beta_total_vol_avg: total plasma beta @@ -222,7 +226,7 @@ def constraint_equation_1(constraint_registration): @ConstraintManager.register_constraint(2, "MW/m3", "=") -def constraint_equation_2(constraint_registration): +def constraint_equation_2(constraint_registration, _data): """ i_rad_loss: switch for radiation loss term usage in power balance (see User Guide): @@ -283,7 +287,7 @@ def constraint_equation_2(constraint_registration): @ConstraintManager.register_constraint(3, "MW/m3", "=") -def constraint_equation_3(constraint_registration): +def constraint_equation_3(constraint_registration, _data): """Global power balance equation for ions i_plasma_ignited: switch for ignition assumption - 0 do not assume plasma ignition; @@ -327,7 +331,7 @@ def constraint_equation_3(constraint_registration): @ConstraintManager.register_constraint(4, "MW/m3", "=") -def constraint_equation_4(constraint_registration): +def constraint_equation_4(constraint_registration, _data): """Global power balance equation for electrons i_rad_loss: switch for radiation loss term usage in power balance @@ -381,7 +385,7 @@ def constraint_equation_4(constraint_registration): @ConstraintManager.register_constraint(5, "/m3", "<=") -def constraint_equation_5(constraint_registration): +def constraint_equation_5(constraint_registration, _data): """Equation for density upper limit fdene: density limit scale @@ -424,7 +428,7 @@ def constraint_equation_5(constraint_registration): @ConstraintManager.register_constraint(6, "", "<=") -def constraint_equation_6(constraint_registration): +def constraint_equation_6(constraint_registration, _data): """Equation for epsilon beta-poloidal upper limit beta_poloidal_eps_max: maximum (eps*beta_poloidal) @@ -442,7 +446,7 @@ def constraint_equation_6(constraint_registration): @ConstraintManager.register_constraint(7, "/m3", "=") -def constraint_equation_7(constraint_registration): +def constraint_equation_7(constraint_registration, _data): """Equation for hot beam ion density i_plasma_ignited: switch for ignition assumption: @@ -466,7 +470,7 @@ def constraint_equation_7(constraint_registration): @ConstraintManager.register_constraint(8, "MW/m2", "<=") -def constraint_equation_8(constraint_registration): +def constraint_equation_8(constraint_registration, _data): """Equation for neutron wall load upper limit pflux_fw_neutron_max_mw: allowable wall-load (MW/m2) @@ -480,7 +484,7 @@ def constraint_equation_8(constraint_registration): @ConstraintManager.register_constraint(9, "MW", "<=") -def constraint_equation_9(constraint_registration): +def constraint_equation_9(constraint_registration, _data): """Equation for fusion power upper limit p_fusion_total_max_mw: maximum fusion power (MW) @@ -494,7 +498,7 @@ def constraint_equation_9(constraint_registration): @ConstraintManager.register_constraint(11, "m", "=") -def constraint_equation_11(constraint_registration): +def constraint_equation_11(constraint_registration, _data): """Equation for radial build rbld: sum of thicknesses to the major radius (m) @@ -508,7 +512,7 @@ def constraint_equation_11(constraint_registration): @ConstraintManager.register_constraint(12, "V.sec", ">=") -def constraint_equation_12(constraint_registration): +def constraint_equation_12(constraint_registration, _data): """Equation for volt-second capability lower limit vs_plasma_total_required: total V-s needed (Wb) @@ -524,7 +528,7 @@ def constraint_equation_12(constraint_registration): @ConstraintManager.register_constraint(13, "sec", ">=") -def constraint_equation_13(constraint_registration): +def constraint_equation_13(constraint_registration, _data): """Equation for burn time lower limit t_plant_pulse_burn: burn time (s) (calculated if i_pulsed_plant=1) @@ -538,7 +542,7 @@ def constraint_equation_13(constraint_registration): @ConstraintManager.register_constraint(14, "", "=") -def constraint_equation_14(constraint_registration): +def constraint_equation_14(constraint_registration, _data): """Equation to fix number of NBI decay lengths to plasma centre n_beam_decay_lengths_core: neutral beam e-decay lengths to plasma centre @@ -552,7 +556,7 @@ def constraint_equation_14(constraint_registration): @ConstraintManager.register_constraint(15, "MW", ">=") -def constraint_equation_15(constraint_registration): +def constraint_equation_15(constraint_registration, _data): """Equation for L-H power threshold limit to enforce H-mode f_h_mode_margin: a margin on the constraint @@ -576,7 +580,7 @@ def constraint_equation_15(constraint_registration): @ConstraintManager.register_constraint(16, "MW", ">=") -def constraint_equation_16(constraint_registration): +def constraint_equation_16(constraint_registration, _data): """Equation for net electric power lower limit p_plant_electric_net_mw: net electric power (MW) @@ -590,7 +594,7 @@ def constraint_equation_16(constraint_registration): @ConstraintManager.register_constraint(17, "MW/m3", "<=") -def constraint_equation_17(constraint_registration): +def constraint_equation_17(constraint_registration, _data): """Equation for radiation power upper limit f_p_alpha_plasma_deposited: fraction of alpha power deposited in plasma @@ -624,7 +628,7 @@ def constraint_equation_17(constraint_registration): @ConstraintManager.register_constraint(18, "MW/m2", "<=") -def constraint_equation_18(constraint_registration): +def constraint_equation_18(constraint_registration, _data): """Equation for divertor heat load upper limit pflux_div_heat_load_max_mw: heat load limit (MW/m2) @@ -638,7 +642,7 @@ def constraint_equation_18(constraint_registration): @ConstraintManager.register_constraint(19, "MVA", "<=") -def constraint_equation_19(constraint_registration): +def constraint_equation_19(constraint_registration, _data): """Equation for MVA (power) upper limit: resistive TF coil set p_cp_resistive_mw: peak resistive TF coil inboard leg power (total) (MW) @@ -656,7 +660,7 @@ def constraint_equation_19(constraint_registration): @ConstraintManager.register_constraint(20, "m", "<=") -def constraint_equation_20(constraint_registration): +def constraint_equation_20(constraint_registration, _data): """Equation for neutral beam tangency radius upper limit radius_beam_tangency_max: maximum tangency radius for centreline of beam (m) @@ -670,7 +674,7 @@ def constraint_equation_20(constraint_registration): @ConstraintManager.register_constraint(21, "", ">=") -def constraint_equation_21(constraint_registration): +def constraint_equation_21(constraint_registration, _data): """Equation for minor radius lower limit rminor: plasma minor radius (m) @@ -684,7 +688,7 @@ def constraint_equation_21(constraint_registration): @ConstraintManager.register_constraint(22, "MW", ">=") -def constraint_equation_22(constraint_registration): +def constraint_equation_22(constraint_registration, _data): """Equation for L-H power threshold limit to enforce L-mode f_l_mode_margin: a margin on the constraint @@ -708,7 +712,7 @@ def constraint_equation_22(constraint_registration): @ConstraintManager.register_constraint(23, "m", "<=") -def constraint_equation_23(constraint_registration): +def constraint_equation_23(constraint_registration, _data): """Equation for conducting shell radius / rminor upper limit rminor: plasma minor radius (m) @@ -735,7 +739,7 @@ def constraint_equation_23(constraint_registration): @ConstraintManager.register_constraint(24, "", "<=") -def constraint_equation_24(constraint_registration): +def constraint_equation_24(constraint_registration, _data): """Equation for beta upper limit i_beta_component: switch for beta limit scaling (constraint equation 24): @@ -798,7 +802,7 @@ def constraint_equation_24(constraint_registration): @ConstraintManager.register_constraint(25, "T", "<=") -def constraint_equation_25(constraint_registration): +def constraint_equation_25(constraint_registration, _data): """Equation for peak toroidal field upper limit b_tf_inboard_max: maximum peak toroidal field (T) @@ -812,7 +816,7 @@ def constraint_equation_25(constraint_registration): @ConstraintManager.register_constraint(26, "A/m2", "<=") -def constraint_equation_26(constraint_registration): +def constraint_equation_26(constraint_registration, _data): """Equation for Central Solenoid current density upper limit at EOF fjohc: margin for central solenoid current at end-of-flattop @@ -830,7 +834,7 @@ def constraint_equation_26(constraint_registration): @ConstraintManager.register_constraint(27, "A/m2", "<=") -def constraint_equation_27(constraint_registration): +def constraint_equation_27(constraint_registration, _data): """Equation for Central Solenoid current density upper limit at BOP fjohc0: margin for central solenoid current at beginning of pulse @@ -848,7 +852,7 @@ def constraint_equation_27(constraint_registration): @ConstraintManager.register_constraint(28, "", ">=") -def constraint_equation_28(constraint_registration): +def constraint_equation_28(constraint_registration, _data): """Equation for fusion gain (big Q) lower limit big_q_plasma: Fusion gain; P_fusion / (P_injection + P_ohmic) @@ -872,7 +876,7 @@ def constraint_equation_28(constraint_registration): @ConstraintManager.register_constraint(29, "m", "=") -def constraint_equation_29(constraint_registration): +def constraint_equation_29(constraint_registration, _data): """Equation for inboard major radius: This is a consistency equation rmajor: plasma major radius (m) (iteration variable 3) @@ -890,7 +894,7 @@ def constraint_equation_29(constraint_registration): @ConstraintManager.register_constraint(30, "MW", "<=") -def constraint_equation_30(constraint_registration): +def constraint_equation_30(constraint_registration, _data): """Equation for injection power upper limit p_hcd_injected_total_mw: total auxiliary injected power (MW) @@ -904,7 +908,7 @@ def constraint_equation_30(constraint_registration): @ConstraintManager.register_constraint(31, "Pa", "<=") -def constraint_equation_31(constraint_registration): +def constraint_equation_31(constraint_registration, _data): """Equation for TF coil case stress upper limit (SCTF) sig_tf_case_max: Allowable maximum shear stress in TF coil case (Tresca criterion) (Pa) @@ -918,7 +922,7 @@ def constraint_equation_31(constraint_registration): @ConstraintManager.register_constraint(32, "Pa", "<=") -def constraint_equation_32(constraint_registration): +def constraint_equation_32(constraint_registration, _data): """Equation for TF coil conduit stress upper limit (SCTF) sig_tf_wp_max: Allowable maximum shear stress in TF coil conduit (Tresca criterion) (Pa) @@ -932,7 +936,7 @@ def constraint_equation_32(constraint_registration): @ConstraintManager.register_constraint(33, "A/m2", "<=") -def constraint_equation_33(constraint_registration): +def constraint_equation_33(constraint_registration, _data): """Equation for TF coil operating/critical J upper limit (SCTF) args : output structure : residual error; constraint value; @@ -955,7 +959,7 @@ def constraint_equation_33(constraint_registration): @ConstraintManager.register_constraint(34, "V", "<=") -def constraint_equation_34(constraint_registration): +def constraint_equation_34(constraint_registration, _data): """Equation for TF coil dump voltage upper limit (SCTF) v_tf_coil_dump_quench_max_kv: max voltage across TF coil during quench (kV) @@ -969,7 +973,7 @@ def constraint_equation_34(constraint_registration): @ConstraintManager.register_constraint(35, "A/m2", "<=") -def constraint_equation_35(constraint_registration): +def constraint_equation_35(constraint_registration, _data): """Equation for TF coil J_wp/J_prot upper limit (SCTF) j_tf_wp_quench_heat_max: allowable TF coil winding pack current density, for dump temperature @@ -984,7 +988,7 @@ def constraint_equation_35(constraint_registration): @ConstraintManager.register_constraint(36, "K", ">=") -def constraint_equation_36(constraint_registration): +def constraint_equation_36(constraint_registration, _data): """Equation for TF coil s/c temperature margin lower limit (SCTF) temp_tf_superconductor_margin: TF coil temperature margin (K) @@ -998,7 +1002,7 @@ def constraint_equation_36(constraint_registration): @ConstraintManager.register_constraint(37, "1E20 A/Wm2", "<=") -def constraint_equation_37(constraint_registration): +def constraint_equation_37(constraint_registration, _data): """Equation for current drive gamma upper limit eta_cd_norm_hcd_primary_max: maximum current drive gamma @@ -1012,7 +1016,7 @@ def constraint_equation_37(constraint_registration): @ConstraintManager.register_constraint(39, "K", "<=") -def constraint_equation_39(constraint_registration): +def constraint_equation_39(constraint_registration, _data): """Equation for first wall temperature upper limit temp_fw_max: maximum temperature of first wall material (K) (i_thermal_electric_conversion>1) @@ -1031,7 +1035,7 @@ def constraint_equation_39(constraint_registration): @ConstraintManager.register_constraint(40, "MW", ">=") -def constraint_equation_40(constraint_registration): +def constraint_equation_40(constraint_registration, _data): """Equation for auxiliary power lower limit p_hcd_injected_total_mw: total auxiliary injected power (MW) @@ -1045,7 +1049,7 @@ def constraint_equation_40(constraint_registration): @ConstraintManager.register_constraint(41, "sec", ">=") -def constraint_equation_41(constraint_registration): +def constraint_equation_41(constraint_registration, _data): """Equation for plasma current ramp-up time lower limit t_plant_pulse_plasma_current_ramp_up: plasma current ramp-up time for current initiation (s) @@ -1059,7 +1063,7 @@ def constraint_equation_41(constraint_registration): @ConstraintManager.register_constraint(42, "sec", ">=") -def constraint_equation_42(constraint_registration): +def constraint_equation_42(constraint_registration, _data): """Equation for cycle time lower limit t_plant_pulse_total: full cycle time (s) @@ -1078,7 +1082,7 @@ def constraint_equation_42(constraint_registration): @ConstraintManager.register_constraint(43, "deg C", "=") -def constraint_equation_43(constraint_registration): +def constraint_equation_43(constraint_registration, _data): """Equation for average centrepost temperature: This is a consistency equation (TART) temp_cp_average: average temp of TF coil inboard leg conductor (C)e @@ -1101,7 +1105,7 @@ def constraint_equation_43(constraint_registration): @ConstraintManager.register_constraint(44, "deg C", "<=") -def constraint_equation_44(constraint_registration): +def constraint_equation_44(constraint_registration, _data): """Equation for centrepost temperature upper limit (TART) temp_cp_max: maximum peak centrepost temperature (K) @@ -1126,7 +1130,7 @@ def constraint_equation_44(constraint_registration): @ConstraintManager.register_constraint(45, "", ">=") -def constraint_manager_45(constraint_registration): +def constraint_manager_45(constraint_registration, _data): """Equation for edge safety factor lower limit (TART) q95 : safety factor 'near' plasma edge @@ -1147,7 +1151,7 @@ def constraint_manager_45(constraint_registration): @ConstraintManager.register_constraint(46, "", "<=") -def constraint_equation_46(constraint_registration): +def constraint_equation_46(constraint_registration, _data): """Equation for Ip/Irod upper limit (TART) eps: inverse aspect ratio @@ -1174,7 +1178,7 @@ def constraint_equation_46(constraint_registration): @ConstraintManager.register_constraint(48, "", "<=") -def constraint_equation_48(constraint_registration): +def constraint_equation_48(constraint_registration, _data): """Equation for poloidal beta upper limit beta_poloidal_max: maximum poloidal beta @@ -1188,7 +1192,7 @@ def constraint_equation_48(constraint_registration): @ConstraintManager.register_constraint(50, "Hz", "<=") -def constraint_equation_50(constraint_registration): +def constraint_equation_50(constraint_registration, _data): """IFE option: Equation for repetition rate upper limit""" return leq( data_structure.ife_variables.reprat, @@ -1198,7 +1202,7 @@ def constraint_equation_50(constraint_registration): @ConstraintManager.register_constraint(51, "V.s", "=") -def constraint_equation_51(constraint_registration): +def constraint_equation_51(constraint_registration, _data): """Equation to enforce startup flux = available startup flux vs_plasma_res_ramp: resistive losses in startup V-s (Wb) @@ -1216,7 +1220,7 @@ def constraint_equation_51(constraint_registration): @ConstraintManager.register_constraint(52, "", ">=") -def constraint_equation_52(constraint_registration): +def constraint_equation_52(constraint_registration, _data): """Equation for tritium breeding ratio lower limit The tritium breeding ratio is only calculated when using the IFE model. @@ -1237,7 +1241,7 @@ def constraint_equation_52(constraint_registration): @ConstraintManager.register_constraint(53, "neutron/m2", "<=") -def constraint_equation_53(constraint_registration): +def constraint_equation_53(constraint_registration, _data): """Equation for fast neutron fluence on TF coil upper limit nflutfmax: max fast neutron fluence on TF coil (n/m2) @@ -1251,7 +1255,7 @@ def constraint_equation_53(constraint_registration): @ConstraintManager.register_constraint(54, "MW/m3", "<=") -def constraint_equation_54(constraint_registration): +def constraint_equation_54(constraint_registration, _data): """Equation for peak TF coil nuclear heating upper limit ptfnucmax: maximum nuclear heating in TF coil (MW/m3) @@ -1265,7 +1269,7 @@ def constraint_equation_54(constraint_registration): @ConstraintManager.register_constraint(56, "MW/m", "<=") -def constraint_equation_56(constraint_registration): +def constraint_equation_56(constraint_registration, _data): """Equation for power through separatrix / major radius upper limit pseprmax: maximum ratio of power crossing the separatrix to plasma major radius (Psep/R) (MW/m) @@ -1283,7 +1287,7 @@ def constraint_equation_56(constraint_registration): @ConstraintManager.register_constraint(59, "", "<=") -def constraint_equation_59(constraint_registration): +def constraint_equation_59(constraint_registration, _data): """Equation for neutral beam shine-through fraction upper limit f_p_beam_shine_through_max: maximum neutral beam shine-through fraction @@ -1297,7 +1301,7 @@ def constraint_equation_59(constraint_registration): @ConstraintManager.register_constraint(60, "K", ">=") -def constraint_equation_60(constraint_registration): +def constraint_equation_60(constraint_registration, _data): """Equation for Central Solenoid s/c temperature margin lower limit temp_cs_superconductor_margin: Central solenoid temperature margin (K) @@ -1311,7 +1315,7 @@ def constraint_equation_60(constraint_registration): @ConstraintManager.register_constraint(61, "", ">=") -def constraint_equation_61(constraint_registration): +def constraint_equation_61(constraint_registration, _data): """Equation for availability lower limit f_t_plant_available: Total plant availability fraction @@ -1325,7 +1329,7 @@ def constraint_equation_61(constraint_registration): @ConstraintManager.register_constraint(62, "", ">=") -def constraint_equation_62(constraint_registration): +def constraint_equation_62(constraint_registration, _data): """Lower limit on f_alpha_energy_confinement the ratio of alpha particle to energy confinement times t_alpha_confinement: alpha particle confinement time (s) @@ -1341,7 +1345,7 @@ def constraint_equation_62(constraint_registration): @ConstraintManager.register_constraint(63, "", "<=") -def constraint_equation_63(constraint_registration): +def constraint_equation_63(constraint_registration, _data): """Upper limit on n_iter_vacuum_pumps (i_vacuum_pumping = simple) tfno: number of TF coils (default = 50 for stellarators) @@ -1355,7 +1359,7 @@ def constraint_equation_63(constraint_registration): @ConstraintManager.register_constraint(64, "", "<=") -def constraint_equation_64(constraint_registration): +def constraint_equation_64(constraint_registration, _data): """Upper limit on Zeff zeff_max: maximum value for Zeff @@ -1369,7 +1373,7 @@ def constraint_equation_64(constraint_registration): @ConstraintManager.register_constraint(65, "Pa", "<=") -def constraint_equation_65(constraint_registration): +def constraint_equation_65(constraint_registration, _data): """Upper limit on stress of the vacuum vessel that occurs when the TF coil quenches. max_vv_stress: Maximum permitted stress of the VV (Pa) @@ -1383,7 +1387,7 @@ def constraint_equation_65(constraint_registration): @ConstraintManager.register_constraint(66, "MW", "<=") -def constrain_equation_66(constraint_registration): +def constrain_equation_66(constraint_registration, _data): """Upper limit on rate of change of energy in poloidal field maxpoloidalpower: Maximum permitted absolute rate of change of stored energy in poloidal field (MW) @@ -1397,7 +1401,7 @@ def constrain_equation_66(constraint_registration): @ConstraintManager.register_constraint(67, "MW/m2", "<=") -def constraint_equation_67(constraint_registration): +def constraint_equation_67(constraint_registration, _data): """Simple upper limit on radiation wall load pflux_fw_rad_max: Maximum permitted radiation wall load (MW/m^2) @@ -1411,7 +1415,7 @@ def constraint_equation_67(constraint_registration): @ConstraintManager.register_constraint(68, "MWT/m", "<=") -def constraint_equation_68(constraint_registration): +def constraint_equation_68(constraint_registration, _data): """Upper limit on Psep scaling (PsepB/qAR) psepbqarmax: maximum permitted value of ratio of Psep*Bt/qAR (MWT/m) @@ -1458,7 +1462,7 @@ def constraint_equation_68(constraint_registration): @ConstraintManager.register_constraint(72, "Pa", "<=") -def constraint_equation_72(constraint_registration): +def constraint_equation_72(constraint_registration, _data): """Upper limit on central Solenoid Tresca yield stress In the case if the bucked and wedged option ( i_tf_bucking >= 2 ) the constrained @@ -1499,7 +1503,7 @@ def constraint_equation_72(constraint_registration): @ConstraintManager.register_constraint(73, "MW", ">=") -def constraint_equation_73(constraint_registration): +def constraint_equation_73(constraint_registration, _data): """Lower limit to ensure separatrix power is greater than the L-H power + auxiliary power Related to constraint 15 @@ -1518,7 +1522,7 @@ def constraint_equation_73(constraint_registration): @ConstraintManager.register_constraint(74, "K", "<=") -def constraint_equation_74(constraint_registration): +def constraint_equation_74(constraint_registration, _data): """Upper limit to ensure TF coil quench temperature < temp_croco_quench_max ONLY used for croco HTS coil @@ -1533,7 +1537,7 @@ def constraint_equation_74(constraint_registration): @ConstraintManager.register_constraint(75, "A/m2", "<=") -def constraint_equation_75(constraint_registration): +def constraint_equation_75(constraint_registration, _data): """Upper limit to ensure that TF coil current / copper area < Maximum value ONLY used for croco HTS coil @@ -1548,7 +1552,7 @@ def constraint_equation_75(constraint_registration): @ConstraintManager.register_constraint(76, "m-3", "<=") -def constraint_equation_76(constraint_registration): +def constraint_equation_76(constraint_registration, _data): """Upper limit for Eich critical separatrix density model: Added for issue 558 Eich critical separatrix density model @@ -1586,7 +1590,7 @@ def constraint_equation_76(constraint_registration): @ConstraintManager.register_constraint(77, "A/turn", "<=") -def constraint_equation_77(constraint_registration): +def constraint_equation_77(constraint_registration, _data): """Equation for maximum TF current per turn upper limit c_tf_turn_max : allowable TF coil current per turn [A/turn] @@ -1600,7 +1604,7 @@ def constraint_equation_77(constraint_registration): @ConstraintManager.register_constraint(78, "", ">=") -def constraint_equation_78(constraint_registration): +def constraint_equation_78(constraint_registration, _data): """Equation for Reinke criterion, divertor impurity fraction lower limit fzmin : input : minimum impurity fraction from Reinke model @@ -1614,7 +1618,7 @@ def constraint_equation_78(constraint_registration): @ConstraintManager.register_constraint(79, "A/turn", "<=") -def constraint_equation_79(constraint_registration): +def constraint_equation_79(constraint_registration, _data): """Equation for maximum CS field b_cs_limit_max: Central solenoid max field limit [T] @@ -1633,7 +1637,7 @@ def constraint_equation_79(constraint_registration): @ConstraintManager.register_constraint(80, "MW", ">=") -def constraint_equation_80(constraint_registration): +def constraint_equation_80(constraint_registration, _data): """Equation for p_plasma_separatrix_mw lower limit args : output structure : residual error; constraint value; residual error in physical units; @@ -1651,7 +1655,7 @@ def constraint_equation_80(constraint_registration): @ConstraintManager.register_constraint(81, "m-3", ">=") -def constraint_equation_81(constraint_registration): +def constraint_equation_81(constraint_registration, _data): """Lower limit to ensure central density is larger that the pedestal one args : output structure : residual error; constraint value; @@ -1669,7 +1673,7 @@ def constraint_equation_81(constraint_registration): @ConstraintManager.register_constraint(82, "m", ">=") -def constraint_equation_82(constraint_registration): +def constraint_equation_82(constraint_registration, _data): """Equation for toroidal consistency of stellarator build toroidalgap: minimal gap between two stellarator coils @@ -1683,7 +1687,7 @@ def constraint_equation_82(constraint_registration): @ConstraintManager.register_constraint(83, "m", ">=") -def constraint_equation_83(constraint_registration): +def constraint_equation_83(constraint_registration, _data): """Equation for radial consistency of stellarator build available_radial_space: avaible space in radial direction as given by each s.-configuration @@ -1697,7 +1701,7 @@ def constraint_equation_83(constraint_registration): @ConstraintManager.register_constraint(84, "", ">=") -def constraint_equation_84(constraint_registration): +def constraint_equation_84(constraint_registration, _data): """Equation for the lower limit of beta beta_vol_avg_min: Lower limit for beta @@ -1711,7 +1715,7 @@ def constraint_equation_84(constraint_registration): @ConstraintManager.register_constraint(85, "years", "=") -def constraint_equation_85(constraint_registration): +def constraint_equation_85(constraint_registration, _data): """Equality constraint for the centerpost (CP) lifetime Depending on the chosen option i_cp_lifetime: @@ -1744,7 +1748,7 @@ def constraint_equation_85(constraint_registration): @ConstraintManager.register_constraint(86, "m", "<=") -def constraint_equation_86(constraint_registration): +def constraint_equation_86(constraint_registration, _data): """Upper limit on the turn edge length in the TF winding pack dx_tf_turn_general: TF coil turn edge length including turn insulation [m] @@ -1758,7 +1762,7 @@ def constraint_equation_86(constraint_registration): @ConstraintManager.register_constraint(87, "MW", "<=") -def constraint_equation_87(constraint_registration): +def constraint_equation_87(constraint_registration, _data): """Equation for TF coil cryogenic power upper limit p_cryo_plant_electric_mw: cryogenic plant power (MW) @@ -1772,7 +1776,7 @@ def constraint_equation_87(constraint_registration): @ConstraintManager.register_constraint(88, "", "<=") -def constraint_equation_88(constraint_registration): +def constraint_equation_88(constraint_registration, _data): """Equation for TF coil vertical strain upper limit (absolute value) str_wp_max: Allowable maximum TF coil vertical strain @@ -1786,7 +1790,7 @@ def constraint_equation_88(constraint_registration): @ConstraintManager.register_constraint(89, "A/m2", "<=") -def constraint_equation_89(constraint_registration): +def constraint_equation_89(constraint_registration, _data): """Upper limit to ensure that the Central Solenoid [OH] coil current / copper area < Maximum value copperaoh_m2: CS coil current at EOF / copper area [A/m2] @@ -1800,7 +1804,7 @@ def constraint_equation_89(constraint_registration): @ConstraintManager.register_constraint(90, "", ">=") -def constraint_equation_90(constraint_registration): +def constraint_equation_90(constraint_registration, data): """Lower limit for CS coil stress load cycles n_cycle: Allowable number of cycles for CS @@ -1808,21 +1812,19 @@ def constraint_equation_90(constraint_registration): """ if ( data_structure.cost_variables.ibkt_life == 1 - and data_structure.cs_fatigue_variables.bkt_life_csf == 1 + and data.cs_fatigue.bkt_life_csf == 1 ): - data_structure.cs_fatigue_variables.n_cycle_min = ( - data_structure.cost_variables.bktcycles - ) + data.cs_fatigue.n_cycle_min = data_structure.cost_variables.bktcycles return geq( - data_structure.cs_fatigue_variables.n_cycle, - data_structure.cs_fatigue_variables.n_cycle_min, + data.cs_fatigue.n_cycle, + data.cs_fatigue.n_cycle_min, constraint_registration, ) @ConstraintManager.register_constraint(91, "MW", ">=") -def constraint_equation_91(constraint_registration): +def constraint_equation_91(constraint_registration, _data): """Lower limit to ensure ECRH te is greater than required te for ignition at lower values for n and B. Or if the design point is ECRH heatable (if i_plasma_ignited==0) stellarators only (but in principle usable also for tokamaks). @@ -1847,7 +1849,7 @@ def constraint_equation_91(constraint_registration): @ConstraintManager.register_constraint(92, "", "=") -def constraint_equation_92(constraint_registration): +def constraint_equation_92(constraint_registration, _data): """Equation for checking is D/T ratio is consistent, and sums to 1. f_plasma_fuel_deuterium: fraction of deuterium ions @@ -1863,7 +1865,7 @@ def constraint_equation_92(constraint_registration): ) -def constraint_eqns(m: int, ieqn: int): +def constraint_eqns(m: int, ieqn: int, data: DataStructure): """Evaluates the constraints given the current state of PROCESS. Parameters @@ -1873,6 +1875,9 @@ def constraint_eqns(m: int, ieqn: int): ieqn : Evaluates the 'ieqn'th constraint equation (index starts at 1) or all equations if <= 0 + data: + DataStructure object to provide data to evaluate the constraints + NOTE: this is only for getting data, not setting it """ @@ -1887,7 +1892,7 @@ def constraint_eqns(m: int, ieqn: int): for i in range(i1, i2): constraint_id = data_structure.numerics.icc[i] - result = ConstraintManager.evaluate_constraint(constraint_id) + result = ConstraintManager.evaluate_constraint(constraint_id, data) tmp_cc, tmp_con, tmp_err = ( result.normalised_residual, diff --git a/process/core/solver/evaluators.py b/process/core/solver/evaluators.py index 3a177b5f5..493c21028 100644 --- a/process/core/solver/evaluators.py +++ b/process/core/solver/evaluators.py @@ -16,15 +16,18 @@ class Evaluators: """Calls models to evaluate function and gradient functions.""" - def __init__(self, models, _x): + def __init__(self, models, data, _x): """Instantiate Caller with model objects. :param models: physics and engineering model objects :type models: process.main.Models + :param data: data structure object for providing constraint + data to the Caller + :type data: DataStructure :param x: optimisation parameters :type x: np.ndarray """ - self.caller = Caller(models) + self.caller = Caller(models, data) def fcnvmc1(self, _n, m, xv, ifail): """Function evaluator for VMCON. diff --git a/process/core/solver/solver_handler.py b/process/core/solver/solver_handler.py index 41ff0dbfa..148ff8c09 100644 --- a/process/core/solver/solver_handler.py +++ b/process/core/solver/solver_handler.py @@ -18,11 +18,15 @@ class SolverHandler: physics and engineering model objects solver_name : str which solver to use, as specified in solver.py + data: DataStructure + data structure object for providing objective/constraint data to + the solver """ - def __init__(self, models, solver_name): + def __init__(self, models, solver_name, data): self.models = models self.solver_name = solver_name + self.data = data def run(self): """Run solver and retry if it fails in certain ways.""" @@ -44,7 +48,7 @@ def run(self): # Evaluators() calculates the objective and constraint functions and # their gradients for a given vector x - evaluators = Evaluators(self.models, x) + evaluators = Evaluators(self.models, self.data, x) # Configure solver for problem self.solver = get_solver(self.solver_name) diff --git a/process/data_structure/cs_fatigue_variables.py b/process/data_structure/cs_fatigue_variables.py index daa826286..b7f726003 100644 --- a/process/data_structure/cs_fatigue_variables.py +++ b/process/data_structure/cs_fatigue_variables.py @@ -1,91 +1,52 @@ -residual_sig_hoop: float = None -"""residual hoop stress in strucutal material (Pa)""" +from dataclasses import dataclass -n_cycle: float = None -"""Allowable number of cycles for CS stress model""" -n_cycle_min: float = None -"""Minimum allowable number of cycles for CS stress model""" +@dataclass +class CSFatigueData: + residual_sig_hoop: float = 2.4e8 + """residual hoop stress in strucutal material (Pa)""" -t_crack_radial: float = None -"""Initial depth of crack in thickness of conduit (m)""" + n_cycle: float = 0.0 + """Allowable number of cycles for CS stress model""" -t_crack_vertical: float = None -"""Inital vertical crack size (m)""" + n_cycle_min: float = 2.0e4 + """Minimum allowable number of cycles for CS stress model""" -dr_cs_turn_conduit: float = None -"""Thickness of CS conductor conduit (m)""" + t_crack_radial: float = 6.0e-3 + """Initial depth of crack in thickness of conduit (m)""" -dz_cs_turn_conduit: float = None -"""Vertical thickness of CS conductor conduit (m)""" + t_crack_vertical: float = 0.89e-3 + """Inital vertical crack size (m)""" -bkt_life_csf: float = None -"""Switch to pass bkt_life cycles to n_cycle_min""" + dr_cs_turn_conduit: float = 0.07 + """Thickness of CS conductor conduit (m)""" -sf_vertical_crack: float = None -"""Safety factor for vertical crack size (-)""" + dz_cs_turn_conduit: float = 0.022 + """Vertical thickness of CS conductor conduit (m)""" -sf_radial_crack: float = None -"""Safety factor for radial crack size (-)""" + bkt_life_csf: float = 0.0 + """Switch to pass bkt_life cycles to n_cycle_min""" -sf_fast_fracture: float = None -"""safety factor for stress intensity factor (-)""" + sf_vertical_crack: float = 2.0 + """Safety factor for vertical crack size (-)""" -paris_coefficient: float = None -"""Paris equation material coefficient (-)""" + sf_radial_crack: float = 2.0 + """Safety factor for radial crack size (-)""" -paris_power_law: float = None -"""Paris equation material power law (-)""" + sf_fast_fracture: float = 1.5 + """safety factor for stress intensity factor (-)""" -walker_coefficient: float = None -"""walker coefficent (-)""" + paris_coefficient: float = 65.0e-14 + """Paris equation material coefficient (-)""" -fracture_toughness: float = None -"""fracture toughness (MPa m^1/2)""" + paris_power_law: float = 3.5 + """Paris equation material power law (-)""" + walker_coefficient: float = 0.436 + """walker coefficent (-)""" -def init_cs_fatigue_variables(): - global residual_sig_hoop - residual_sig_hoop = 2.4e8 + fracture_toughness: float = 2.0e2 + """fracture toughness (MPa m^1/2)""" - global t_crack_radial - t_crack_radial = 6.0e-3 - global t_crack_vertical - t_crack_vertical = 0.89e-3 - - global n_cycle - n_cycle = 0.0 - - global n_cycle_min - n_cycle_min = 2.0e4 - - global dz_cs_turn_conduit - dz_cs_turn_conduit = 0.022 - - global dr_cs_turn_conduit - dr_cs_turn_conduit = 0.07 - - global bkt_life_csf - bkt_life_csf = 0.0 - - global sf_vertical_crack - sf_vertical_crack = 2.0 - - global sf_radial_crack - sf_radial_crack = 2.0 - - global sf_fast_fracture - sf_fast_fracture = 1.5 - - global paris_coefficient - paris_coefficient = 65.0e-14 - - global paris_power_law - paris_power_law = 3.5 - - global walker_coefficient - walker_coefficient = 0.436 - - global fracture_toughness - fracture_toughness = 2.0e2 +CREATE_DICTS_FROM_DATACLASS = CSFatigueData diff --git a/process/main.py b/process/main.py index 17952822c..a66404f69 100644 --- a/process/main.py +++ b/process/main.py @@ -71,7 +71,7 @@ from process.models.divertor import Divertor from process.models.fw import FirstWall from process.models.ife import IFE -from process.models.pfcoil import PFCoil +from process.models.pfcoil import CSCoil, PFCoil from process.models.physics.bootstrap_current import PlasmaBootstrapCurrent from process.models.physics.confinement_time import PlasmaConfinementTime from process.models.physics.current_drive import ( @@ -427,7 +427,7 @@ def run_scan(self): f"Invalid ioptimz value: {data_structure.numerics.ioptimz}. Please " "select either 1 (optimise) or -2 (no optimisation)." ) - self.scan = Scan(self.models, self.solver) + self.scan = Scan(self.models, self.solver, self.data) def show_errors(self): """Report all informational/error messages encountered.""" @@ -600,7 +600,8 @@ def __init__(self, data: DataStructure): self._costs_1990 = Costs() self._costs_2015 = Costs2015() self.cs_fatigue = CsFatigue() - self.pfcoil = PFCoil(cs_fatigue=self.cs_fatigue) + self.cs_coil = CSCoil(cs_fatigue=self.cs_fatigue) + self.pfcoil = PFCoil(cs_fatigue=self.cs_fatigue, cs_coil=self.cs_coil) self.power = Power() self.cryostat = Cryostat() self.build = Build() @@ -703,7 +704,13 @@ def costs(self, value: CostsProtocol): def models(self) -> tuple[Model, ...]: # At the moment, this property just returns models that implement the Model interface. # Eventually every Model will comply and then this method can be used as the caller/outputter! - return (self.water_use, self._costs_2015) + return ( + self.water_use, + self._costs_2015, + self.cs_fatigue, + self.cs_coil, + self.pfcoil, + ) def setup_data_structure(self): # This Models class should be replaced with a dataclass so we can diff --git a/process/models/cs_fatigue.py b/process/models/cs_fatigue.py index ad62c4f33..4e6e98d2a 100644 --- a/process/models/cs_fatigue.py +++ b/process/models/cs_fatigue.py @@ -1,14 +1,20 @@ import numpy as np from numba import njit -from process import data_structure from process.core import constants +from process.core.model import Model -class CsFatigue: +class CsFatigue(Model): def __init__(self): self.outfile = constants.NOUT + def output(self): + """This model doesn't have any output""" + + def run(self): + """This model doesn't need to be run""" + def ncycle( self, max_hoop_stress, @@ -36,8 +42,8 @@ def ncycle( # X. Sarasola et al, IEEE Transactions on Applied Superconductivity, # vol. 30, no. 4, pp. 1-5, June 2020, Art no. 4200705 - n = -data_structure.cs_fatigue_variables.paris_power_law * ( - data_structure.cs_fatigue_variables.walker_coefficient - 1.0e0 + n = -self.data.cs_fatigue.paris_power_law * ( + self.data.cs_fatigue.walker_coefficient - 1.0e0 ) # Set units to MPa @@ -59,7 +65,7 @@ def ncycle( # Calculated constant for a given stress ratio using Walker equation # https://en.wikipedia.org/wiki/Crack_growth_equation#Walker_equation - cr = data_structure.cs_fatigue_variables.paris_coefficient / (1.0e0 - r) ** n + cr = self.data.cs_fatigue.paris_coefficient / (1.0e0 - r) ** n # select given increase in crack area delta = 1.0e-4 @@ -72,20 +78,12 @@ def ncycle( # Default CS steel undergoes fast fracture when SIF > 200 MPa, under a saftey factor 1.5 we use 133MPa pi_2_arr = np.array([np.pi / 2.0e0, 0]) while ( - ( - a - <= dz_cs_turn_conduit - / data_structure.cs_fatigue_variables.sf_vertical_crack - ) - and ( - c - <= dr_cs_turn_conduit - / data_structure.cs_fatigue_variables.sf_radial_crack - ) + (a <= dz_cs_turn_conduit / self.data.cs_fatigue.sf_vertical_crack) + and (c <= dr_cs_turn_conduit / self.data.cs_fatigue.sf_radial_crack) and ( k_max - <= data_structure.cs_fatigue_variables.fracture_toughness - / data_structure.cs_fatigue_variables.sf_fast_fracture + <= self.data.cs_fatigue.fracture_toughness + / self.data.cs_fatigue.sf_fast_fracture ) ): # find SIF max from SIF_a and SIF_c @@ -100,21 +98,11 @@ def ncycle( k_max = max(k_a, k_c) # run euler_method and find number of cycles needed to give crack increase - delta_n = delta / ( - cr * (k_max**data_structure.cs_fatigue_variables.paris_power_law) - ) + delta_n = delta / (cr * (k_max**self.data.cs_fatigue.paris_power_law)) # update a and c, N (+= doesnt work for fortran (?) reasons) - a = ( - a - + delta - * (k_a / k_max) ** data_structure.cs_fatigue_variables.paris_power_law - ) - c = ( - c - + delta - * (k_c / k_max) ** data_structure.cs_fatigue_variables.paris_power_law - ) + a = a + delta * (k_a / k_max) ** self.data.cs_fatigue.paris_power_law + c = c + delta * (k_c / k_max) ** self.data.cs_fatigue.paris_power_law n_pulse = n_pulse + delta_n # two pulses - ramp to Vsmax and ramp down per cycle diff --git a/process/models/pfcoil.py b/process/models/pfcoil.py index 5bce84744..3441c6c5d 100644 --- a/process/models/pfcoil.py +++ b/process/models/pfcoil.py @@ -13,7 +13,6 @@ from process.core.model import Model from process.data_structure import build_variables as bv from process.data_structure import constraint_variables as ctv -from process.data_structure import cs_fatigue_variables as csfv from process.data_structure import fwbs_variables as fwbsv from process.data_structure import ( pfcoil_variables, @@ -33,13 +32,13 @@ class PFCoil(Model): """Calculate poloidal field coil system parameters.""" - def __init__(self, cs_fatigue): + def __init__(self, cs_fatigue, cs_coil): """Initialise Fortran module variables.""" self.outfile = constants.NOUT # output file unit self.mfile = constants.MFILE # mfile file unit pfcoil_variables.init_pfcoil_module() self.cs_fatigue = cs_fatigue - self.cs_coil = CSCoil(cs_fatigue) + self.cs_coil = cs_coil def run(self): """Run the PF coil model.""" @@ -2285,7 +2284,7 @@ def outpf(self): self.outfile, "Residual hoop stress in CS Steel (Pa)", "(residual_sig_hoop)", - csfv.residual_sig_hoop, + self.data.cs_fatigue.residual_sig_hoop, ) op.ovarre( self.outfile, @@ -2297,13 +2296,13 @@ def outpf(self): self.outfile, "Initial vertical crack size (m)", "(t_crack_vertical)", - csfv.t_crack_vertical, + self.data.cs_fatigue.t_crack_vertical, ) op.ovarre( self.outfile, "Initial radial crack size (m)", "(t_crack_radial)", - csfv.t_crack_radial, + self.data.cs_fatigue.t_crack_radial, ) op.ovarre( self.outfile, @@ -2333,26 +2332,26 @@ def outpf(self): self.outfile, "CS structural vertical thickness (m)", "(dz_cs_turn_conduit)", - csfv.dz_cs_turn_conduit, + self.data.cs_fatigue.dz_cs_turn_conduit, ) op.ovarre( self.outfile, "CS structural radial thickness (m)", "(dr_cs_turn_conduit)", - csfv.dr_cs_turn_conduit, + self.data.cs_fatigue.dr_cs_turn_conduit, ) op.ovarre( self.outfile, "Allowable number of cycles till CS fracture", "(n_cycle)", - csfv.n_cycle, + self.data.cs_fatigue.n_cycle, "OP ", ) op.ovarre( self.outfile, "Minimum number of cycles required till CS fracture", "(n_cycle_min)", - csfv.n_cycle_min, + self.data.cs_fatigue.n_cycle_min, "OP ", ) # Check whether CS coil is hitting any limits @@ -2937,7 +2936,7 @@ def waveform(self): pfcoil_variables.f_c_pf_cs_peak_time_array[ic, 5] = 0.0e0 -class CSCoil: +class CSCoil(Model): """Calculate central solenoid coil system parameters.""" def __init__(self, cs_fatigue): @@ -2946,6 +2945,12 @@ def __init__(self, cs_fatigue): self.mfile = constants.MFILE # mfile file unit self.cs_fatigue = cs_fatigue + def output(self): + """This model doesn't have any output""" + + def run(self): + """This model doesn't need to be run""" + def calculate_cs_geometry( self, z_tf_inside_half: float, @@ -3236,8 +3241,8 @@ def ohcalc(self): pfcoil_variables.dz_cs_turn, pfcoil_variables.dr_cs_turn, pfcoil_variables.radius_cs_turn_cable_space, - csfv.dr_cs_turn_conduit, - csfv.dz_cs_turn_conduit, + self.data.cs_fatigue.dr_cs_turn_conduit, + self.data.cs_fatigue.dz_cs_turn_conduit, ) = self.calculate_cs_turn_geometry_eu_demo( a_cs_turn=pfcoil_variables.a_cs_turn, f_dr_dz_cs_turn=pfcoil_variables.f_dr_dz_cs_turn, @@ -3349,12 +3354,15 @@ def ohcalc(self): # Calculation of CS fatigue # this is only valid for pulsed reactor design if pv.f_c_plasma_inductive > 0.0e-4: - csfv.n_cycle, csfv.t_crack_radial = self.cs_fatigue.ncycle( + ( + self.data.cs_fatigue.n_cycle, + self.data.cs_fatigue.t_crack_radial, + ) = self.cs_fatigue.ncycle( pfcoil_variables.sig_hoop, - csfv.residual_sig_hoop, - csfv.t_crack_vertical, - csfv.dz_cs_turn_conduit, - csfv.dr_cs_turn_conduit, + self.data.cs_fatigue.residual_sig_hoop, + self.data.cs_fatigue.t_crack_vertical, + self.data.cs_fatigue.dz_cs_turn_conduit, + self.data.cs_fatigue.dr_cs_turn_conduit, ) # Now steel area fraction is iteration variable and constraint @@ -3660,14 +3668,14 @@ def output_cs_structure(self): self.outfile, "Radial thickness of steel conduit to cable space [m]", "(dr_cs_turn_conduit)", - csfv.dr_cs_turn_conduit, + self.data.cs_fatigue.dr_cs_turn_conduit, "OP ", ) op.ovarre( self.outfile, "Vertical thickness of steel conduit to cable space [m]", "(dz_cs_turn_conduit)", - csfv.dz_cs_turn_conduit, + self.data.cs_fatigue.dz_cs_turn_conduit, "OP ", ) op.ovarre( diff --git a/tests/conftest.py b/tests/conftest.py index 67b83f430..72b1f2b31 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -238,7 +238,10 @@ def _plot_show_and_close_class(request): @pytest.fixture def process_models(): - return Models(DataStructure()) + models = Models(DataStructure()) + for model in models.models: + model.data = models.data + return models @pytest.fixture() diff --git a/tests/integration/test_pfcoil_int.py b/tests/integration/test_pfcoil_int.py index 0286241db..564e0538a 100644 --- a/tests/integration/test_pfcoil_int.py +++ b/tests/integration/test_pfcoil_int.py @@ -14,16 +14,13 @@ from numpy.testing import assert_array_almost_equal from process.core import constants -from process.core.init import init_all_module_vars from process.data_structure import build_variables as bv from process.data_structure import fwbs_variables as fwbsv from process.data_structure import pfcoil_variables, superconducting_tf_coil_variables from process.data_structure import physics_variables as pv from process.data_structure import tfcoil_variables as tfv from process.data_structure import times_variables as tv -from process.models.cs_fatigue import CsFatigue from process.models.pfcoil import ( - CSCoil, PFCoil, fixb, mtrx, @@ -33,25 +30,24 @@ @pytest.fixture -def pfcoil(): +def pfcoil(process_models): """Fixture to create PFCoil object. :return: a PFCoil instance :rtype: process.pfcoil.PFCoil """ - init_all_module_vars() - return PFCoil(cs_fatigue=CsFatigue()) + return process_models.pfcoil @pytest.fixture -def cs_coil(): - """Fixture to create CSCoil object. +def cs_coil(process_models): + """Fixture to get the CSCoil from process_models. :return: a CSCoil instance :rtype: process.pfcoil.CSCoil """ - return CSCoil(cs_fatigue=CsFatigue()) + return process_models.cs_coil def test_pfcoil(monkeypatch, pfcoil): diff --git a/tests/unit/test_constraints.py b/tests/unit/test_constraints.py index de24a38e8..d4c3aa886 100644 --- a/tests/unit/test_constraints.py +++ b/tests/unit/test_constraints.py @@ -4,6 +4,7 @@ from process.core.exceptions import ProcessValueError from process.core.init import init_all_module_vars +from process.core.model import DataStructure from process.core.solver.constraints import ConstraintManager @@ -24,4 +25,6 @@ def test_constraint_functions(constraint_registration): # default flags (i_pulsed_plant=0 or itart=0). with contextlib.suppress(ZeroDivisionError, ProcessValueError): # call the constraint equation and check no error occurs - constraint_registration.constraint_equation(constraint_registration) + constraint_registration.constraint_equation( + constraint_registration, DataStructure() + ) diff --git a/tests/unit/test_cs_fatigue.py b/tests/unit/test_cs_fatigue.py index a2441f314..fbf9fc0d8 100644 --- a/tests/unit/test_cs_fatigue.py +++ b/tests/unit/test_cs_fatigue.py @@ -1,38 +1,27 @@ -from typing import Any, NamedTuple +from typing import NamedTuple import pytest -from process.models.cs_fatigue import CsFatigue - @pytest.fixture -def cs_fatigue_python(): +def cs_fatigue_python(process_models): """Fixture to create a CsFatigue object. :return: an instance of CsFatigue :rtype: process.cs_fatigue.CsFatigue """ - return CsFatigue() + return process_models.cs_fatigue class NcycleParam(NamedTuple): - max_hoop_stress: Any = None - - residual_stress: Any = None - - t_crack_vertical: Any = None - - dz_cs_turn_conduit: Any = None - - dr_cs_turn_conduit: Any = None - - n_cycle: Any = None - - t_crack_radial: Any = None - - expected_n_cycle: Any = None - - expected_t_crack_radial: Any = None + max_hoop_stress: float + residual_stress: float + t_crack_vertical: float + dz_cs_turn_conduit: float + dr_cs_turn_conduit: float + t_crack_radial: float + expected_n_cycle: float + expected_t_crack_radial: float @pytest.mark.parametrize( diff --git a/tests/unit/test_pfcoil.py b/tests/unit/test_pfcoil.py index f78e0aaed..8884dc6f3 100644 --- a/tests/unit/test_pfcoil.py +++ b/tests/unit/test_pfcoil.py @@ -10,18 +10,18 @@ from process.data_structure import pfcoil_variables from process.data_structure import tfcoil_variables as tfv from process.models.cs_fatigue import CsFatigue -from process.models.pfcoil import CSCoil, PFCoil, calculate_b_field_at_point, rsid +from process.models.pfcoil import CSCoil, calculate_b_field_at_point, rsid @pytest.fixture -def pfcoil(): +def pfcoil(process_models): """Fixture to create a PFCoil object. :return: an instance of PFCoil :rtype: process.pfcoil.PFCoil """ - return PFCoil(cs_fatigue=CsFatigue()) + return process_models.pfcoil @pytest.fixture