diff --git a/documentation/source/physics-models/plasma_beta/plasma_beta.md b/documentation/source/physics-models/plasma_beta/plasma_beta.md index 20e81f5a14..634adc1850 100644 --- a/documentation/source/physics-models/plasma_beta/plasma_beta.md +++ b/documentation/source/physics-models/plasma_beta/plasma_beta.md @@ -390,7 +390,7 @@ $$ ### Beta consistency -This constraint can be activated by stating `icc = 1` in the input file. +`beta_total_vol_avg` is calculated directly from the current plasma state and is no longer an iteration variable. Constraint 1 may be activated as an optional global consistency check or safeguard against potential non-idempotent behaviour by stating `icc = 1` in the input file. Ensures the relationship between $\beta$, density, temperature and total magnetic field is withheld by checking the fixed input or iteration variable $\texttt{beta}$ is consistent in value with the rest of the physics parameters @@ -404,8 +404,6 @@ $$ \langle n_{\text{e}}T_{\text{e}} \rangle_{\text{V}} = \overbrace{\langle n_{\text{e}} \rangle_{\text{V}}}^{\texttt{nd_plasma_electrons_vol_avg}} \times \overbrace{\langle T_{\text{e}} \rangle_{\text{n}}}^{\texttt{temp_plasma_electron_density_weighted_kev}} $$ -**It is highly recommended to always have this constraint on as it is a global consistency checker** - ---------------- ### Poloidal beta and inverse aspect upper limit diff --git a/documentation/source/physics-models/profiles/plasma_profiles.md b/documentation/source/physics-models/profiles/plasma_profiles.md index 362da0f76e..4697d4f4df 100644 --- a/documentation/source/physics-models/profiles/plasma_profiles.md +++ b/documentation/source/physics-models/profiles/plasma_profiles.md @@ -405,7 +405,7 @@ $$ With the coefficients used to turn the temperature from $\text{keV}$ back to Joules. -We calculate the volume averaged thermal plasma pressure using the volume averaged density and the density weighted volume averaged temperature for the electrons and ions. The density weighted temperature is used as $\langle nT \rangle \neq \langle n \rangle \langle T \rangle$. This can be seen also in the [equality constraint for thermal $\beta$](../plasma_beta/plasma_beta.md#beta-consistency) +We calculate the volume averaged thermal plasma pressure using the volume averaged density and the density weighted volume averaged temperature for the electrons and ions. The density weighted temperature is used as $\langle nT \rangle \neq \langle n \rangle \langle T \rangle$. $$ \langle p_{\text{thermal}} \rangle = \overbrace{\langle n_{\text{e}} \rangle_{\text{V}}}^{\texttt{nd_plasma_electrons_vol_avg}} \times \overbrace{\langle T_{\text{e}} \rangle_{\text{n}}}^{\texttt{temp_plasma_electron_density_weighted_kev}}\\ diff --git a/documentation/source/solver/solver-guide.md b/documentation/source/solver/solver-guide.md index be891cce1a..ec43d3bdb3 100644 --- a/documentation/source/solver/solver-guide.md +++ b/documentation/source/solver/solver-guide.md @@ -1,163 +1,247 @@ +# Solver PROCESS docs -# Constraint Equations +## Constraint Equations -Any computer program naturally contains many equations. The built-in equation -solvers within PROCESS act on a special class, known as constraint equations, -all of which are formulated in the source file `process/core/solver/constraints.py`. These -can be split into two types: - -**Equality constraints (consistency equations)** that enforce consistency between the physics and -engineering parameters in the various models. +PROCESS models a fusion power plant at a **self-consistent operating point**. It contains a large set of coupled physics and engineering models that describe different aspects of the plant. For a design to constitute a PROCESS solution, quantities calculated by different parts of the code must satisfy the required relationships between them. -**Inequality constraints (limit equations)** that enforce a value be greater/less than or equal to some bound (limit). +PROCESS represents these relationships using **constraint equations**, which are defined in `process/core/solver/constraints.py`. -A PROCESS input file will, for example, define which constraint equations are being used as follows: +Constraint equations fall into two classes: -``` -... -n_equality_constraints = 3 +- **Equality constraints (consistency equations)** define relationships that must be satisfied for the calculated plant to be self-consistent. +- **Inequality constraints (limit equations)** define physics or engineering limits that determine whether a solution is feasible. -* Equalities -icc = 1 * Beta -icc = 2 * Global power balance -icc = 11 * Radial build +PROCESS uses user-defined **iteration variables** as the degrees of freedom available to the solver when solving the selected constraint problem. -* Inequalities -icc = 9 * Fusion power upper limit -icc = 5 * Density upper limit -icc = 24 * Beta upper limit -icc = 15 * LH power threshold limit -... -``` +## PROCESS Solutions + +At the heart of a PROCESS solution is a plasma operating in a self-consistent equilibrium state. Its state cannot be chosen arbitrarily: plasma temperature, density, heating, fusion reactions, radiation, transport losses and fuelling are all coupled. + +The principal plasma state is described by the electron temperature and density, + +$$ +T_e,\qquad n_e. +$$ + +To determine these two quantities, two independent governing equations are required. + +### Governing plasma equations + +In PROCESS these are: + +1. **Global plasma power balance - constraint 2** +2. **Fuel-ion equilibrium - constraint 93** + +Together these equations close the plasma system and define the equilibrium plasma state solved by PROCESS. + +For a given set of design parameters, a plasma solution is therefore a pair + +$$ +(T_e,n_e) +$$ + +for which both governing equations are satisfied simultaneously. + +#### Global plasma power balance + +For the plasma to remain in a steady operating state, the power supplied to the plasma must balance the power lost from it. + +Schematically, + +$$ +P_{\mathrm{heating}} = P_{\mathrm{loss}}. +$$ + +The heating and loss terms depend on the plasma state and on the selected PROCESS physics models. The power balance therefore provides one relationship between plasma temperature and density. + +PROCESS represents this condition using **constraint 2**. + +#### Fuel-ion equilibrium + +The plasma must also satisfy particle balance for the fuel ions. + +Fusion reactions consume fuel ions. To maintain the plasma state, this consumption must be consistent with the calculated fuelling rate and fuel burn-up. + +Schematically, + +$$ +\text{fuel supplied} = \text{fuel required to sustain the fusion reaction rate}. +$$ + +PROCESS represents this condition using **constraint 93**. + +Together, global plasma power balance and fuel-ion equilibrium form a coupled system in $T_e$ and $n_e$, whose simultaneous solution determines the self-consistent plasma state. + +### Other consistency equations + +The governing plasma equations are not necessarily the only equality constraints required for a complete PROCESS calculation. -Here each `icc=n` statement tells PROCESS to activate a constraint with the name `n`. A list of the constraints and -their corresponding names can be found [here](../../source/reference/process/data_structure/numerics/#process.data_structure.numerics.lablcc). +Other physics and engineering models may introduce additional consistency relationships, for example for the machine radial build. -The `n_equality_constraints = 3` statement is telling PROCESS to treat the first `3` equations as equality constraints, and the rest as inequality constraints. Therefore, it is imperative that all equality constraints are stated before any inequality constraints. +A **PROCESS solution** is therefore a point at which the governing plasma equations and all other required equality constraints are satisfied. +### Feasible solutions -In both types of equations, an optimiser/solver uses the normalised residuals $c_i$ of the constraints (and sometimes its gradient, depending on the solver/optimiser) to guide the solution towards one that satisfies all of the constraints. +A self-consistent PROCESS solution is not necessarily feasible. -## Consistency Equations +PROCESS therefore also defines **inequality constraints**, or **limit equations**, which represent physics and engineering limits on the design. -Consistency equations are equalities that ensure that the machine produced by PROCESS is -self-consistent. This means, therefore, that many of these constraint equations should -always be used, namely equations 1, 2, 10 and 11. Equation 7 should also be activated -if neutral beam injection is used. The other consistency equations can be activated if -required. A typical consistency equation ensures that two functions $g$ and $h$ are equal: +Examples include limits on plasma density, beta, fusion power, neutron wall load, divertor power loading, magnet stresses and burn time. + +This gives an important distinction: $$ -g(x, y, z, \cdots) = h(x, y, z, \cdots) +\text{PROCESS solution} += +\text{required equality constraints satisfied} $$ +whereas + $$ -c_i = 1 - \frac{g}{h} +\text{feasible PROCESS solution} += +\text{solution that also satisfies the required inequality constraints} $$ -The optimiser/solver will attempt to find a solution that produces $c_i = 0$ for all equality constraints. +A point may therefore be a valid solution of the PROCESS governing equations while still being infeasible because one or more design limits are violated. + +## Specifying constraint equations + +Constraint equations are selected using the `icc` array. For example: + +```text +n_equality_constraints = 3 + +* Equalities +icc = 2 * Global plasma power balance +icc = 11 * Radial build +icc = 93 * Fuel-ion equilibrium + +* Inequalities +icc = 5 * Density upper limit +icc = 9 * Fusion power upper limit +icc = 15 * L-H power threshold limit +icc = 24 * Beta upper limit +``` + +Each `icc = n` statement activates constraint `n`. -## Limit Equations +`n_equality_constraints` specifies how many entries at the beginning of the `icc` array are treated as equality constraints. All equality constraints must therefore be listed before any inequality constraints. -The limit equations are inequalities that ensure that various physics or engineering -limits are not exceeded. +A list of the available constraints and their corresponding names can be found [here](../../source/reference/process/data_structure/numerics/#process.data_structure.numerics.lablcc). -As with the consistency equations, the general form of the limit equations is +For a typical equality constraint, $$ -c_i = 1 - \frac{h_{max}}{h} +g(\mathbf{x}) = h(\mathbf{x}), $$ -where $h_{max}$ is the maximum allowed value of the quantity $h$, or +PROCESS forms the normalised residual $$ -c_i = 1 - f.\frac{h}{h_{min}} +c_i = 1 - \frac{g}{h}, $$ -where $h_{min}$ is the minimum allowed value of the quantity $h$. +which is satisfied when -For example, to set the net electric power to a certain value, the following -should be carried out: +$$ +c_i = 0. +$$ -1. Activate `constraint 16` (net electric power lower limit) by including it in the `icc` array -2. Set `p_plant_electric_net_required_mw` to the required net electric power. +For inequality constraints, the residual is formulated so that the permitted region satisfies -The optimiser/solver will attempt to find a solution that produces $c_i \geq 0$ for all inequality constraints. +$$ +c_i \geq 0. +$$ ## Iteration Variables -... +Iteration variables are the quantities that the solver is allowed to vary in order to satisfy the active constraints. + +They are selected using the `ixc` array, for example: + +```text +ixc = 4 * temp_plasma_electron_vol_avg_kev +ixc = 6 * nd_plasma_electrons_vol_avg +``` + +The equations are coupled, so iteration variables should not generally be thought of as corresponding one-to-one with individual constraints. Changing one iteration variable may affect several constraint residuals. ## Figure of Merit -In optimisation mode, PROCESS finds the self-consistent set of iteration -variable values that maximises or minimises a certain function of them, -known as the figure of merit. +In optimisation mode there may be many feasible solutions. PROCESS selects between them using a **figure of merit**, or objective function. -Several possible figures of merit are available, all of which are in the -source file `evaluators.f90`. +The switch `i_figure_merit` selects the objective. A positive value indicates minimisation, while a negative value indicates maximisation. -Switch `i_figure_merit` is used to control which figure of merit is to be used. If the -figure of merit is to be minimised, `i_figure_merit` should be **positive**, and if a -maximised figure of merit is desired, `i_figure_merit` should be **negative**. +The optimisation problem can therefore be viewed as finding the feasible PROCESS solution that gives the best value of the selected figure of merit. ## Convergence -... +PROCESS solves a nonlinear, coupled constrained problem iteratively. Numerical convergence indicates that the solver has satisfied its convergence criteria. The constraint residuals should then be inspected to determine whether the required consistency equations and limits have been satisfied. -## Optimisation mode +For equality constraints, -Switch `i_process_run_mode` should be set to 1 for optimisation mode. +$$ +c_i \approx 0. +$$ -If `i_process_run_mode = 0`, a non-optimisation pass is performed first. Occasionally this provides a feasible set of initial conditions that aids convergence of the optimiser, but it is recommended to use `i_process_run_mode = 1`. +Numerical convergence should not be confused with feasibility: a converged calculation may still violate one or more inequality constraints. + +## Optimisation mode -Enable all the relevant consistency equations, and it is advisable to enable the corresponding iterations variables. A number of limit equations (inequality constraints) can also be activated. In optimisation mode, the number of iteration variables is unlimited. +Optimisation mode is selected using -It may still be difficult, if not impossible, to reconcile the fusion power and the net electric power with the required values. This may well be due to the power conversion efficiency values being used. +```text +i_process_run_mode = 1 +``` -If scan of a given variable are to be made over a large range of values, it is often a good idea to start the scan in the middle of the desired range, and to split the in two - one going downwards from the initial value, and the other upwards. This ensures that the whole range of the scan produces well-converged machines (assuming a "good" initial point), without sharp changes in gradient in the parameter values. +In this mode, PROCESS varies the active iteration variables in order to: -It should be remembered that the value of the scan variable is set in the array `sweep`, and this overrules any value set for the variable elsewhere in the input file. +1. satisfy the equality constraints; +2. satisfy the inequality constraints; and +3. optimise the selected figure of merit. -The output from an optimisation run contains an indication as to which iteration variables lie at their limit values. +The number of iteration variables may exceed the number of equality constraints because the additional degrees of freedom allow the optimiser to explore the feasible design space. ## Evaluation mode -Evaluation mode is used to evaluate models for a given set of input parameters (a "point") whilst ensuring that the models are self-consistent. It can also be used to perform benchmark comparison, whereby the machine size, output power etc. are known and one only wishes to find the calculated stresses, beta values and fusion powers, for example. +Evaluation mode is selected using -Running `PROCESS` in evaluation mode requires few changes to be made to the input file from the optimisation case. The main differences between optimisation mode and evaluation mode are: +```text +i_process_run_mode = -2 +``` -1. Evaluation mode does not apply lower or upper bounds to the iteration variables. +It is intended for evaluating a specified design while maintaining model consistency, rather than searching for an optimum. -1. Inequality constraints are not enforced. +In evaluation mode: -2. The number of iteration variables must be equal to the number of equality constraints being solved. +1. equality constraints are solved; +2. inequality constraints may be reported but are not enforced; +3. the number of iteration variables must equal the number of equality constraints; +4. no figure of merit is optimised; and +5. iteration-variable optimisation bounds are not applied. -3. An objective function is not available, as no optimisation is taking place. +For example, for a specified design point the plasma state can be obtained by solving the two governing plasma equations for $T_e$ and $n_e$: -As before, the user must decide which constraint equations and iteration variables to activate. For example, an extract from an input file might look like: -``` -* Evaluation problem: evaluate models consistently by solving equality constraints only -i_process_run_mode = -2 * evaluation mode +```text +i_process_run_mode = -2 -*---------------Constraint Equations---------------* -* Define number of equality constraints n_equality_constraints = 2 * Equalities -icc = 1 * Beta -icc = 2 * Global power balance +icc = 2 * Global plasma power balance +icc = 93 * Fuel-ion equilibrium * Inequalities -* Not enforced, but values reported -icc = 5 * Density upper limit -icc = 8 * Neutron wall load upper limit -icc = 9 * Fusion power upper limit -... - -*---------------Iteration Variables----------------* -* Used to solve equality constraints -ixc = 4 * temp_plasma_electron_vol_avg_kev -ixc = 6 * nd_plasma_electrons_vol_avg -... +* Reported but not enforced +icc = 5 * Density upper limit +icc = 8 * Neutron wall load upper limit +icc = 9 * Fusion power upper limit + +* Iteration variables +ixc = 4 * temp_plasma_electron_vol_avg_kev +ixc = 6 * nd_plasma_electrons_vol_avg ``` -The beta and plasma power balance equality constraints are recommended to ensure model consistency, and the `temp_plasma_electron_vol_avg_kev` and `nd_plasma_electrons_vol_avg` iteration variables used to solve them. + +This is useful for benchmark calculations or for evaluating a prescribed machine design without performing an optimisation. \ No newline at end of file diff --git a/examples/data/large_tokamak_eval_IN.DAT b/examples/data/large_tokamak_eval_IN.DAT index 96fe08dcc4..067d4c25fe 100644 --- a/examples/data/large_tokamak_eval_IN.DAT +++ b/examples/data/large_tokamak_eval_IN.DAT @@ -6,8 +6,8 @@ i_process_run_mode = -2 n_equality_constraints = 2 * Equalities -icc = 1 * Beta icc = 2 * Global power balance +icc = 93 * Fuel-ion equilibrium * Inequalities icc = 5 * Density upper limit @@ -101,7 +101,7 @@ p_fusion_total_max_mw = 3000 * maximum fusion power (MW) (`constraint equation p_div_bt_q_aspect_rmajor_max_mw = 10.0 * maximum ratio of Psep*Bt/qAR (MWT/m) (`constraint equation 68`) t_burn_min = 7200.0 * minimum burn time (s) (KE - no longer itv;; see issue #706) pflux_fw_neutron_max_mw = 2.0 * allowable neutron wall-load (MW/m2) (`constraint equation 8`) -f_t_alpha_energy_confinement_min = 5.0 * Lower limit on f_t_alpha_energy_confinementnt the ratio of alpha particle to energy confinement +f_t_alpha_energy_confinement_min = 5.0 * Lower limit on f_t_alpha_energy_confinement the ratio of alpha particle to energy confinement *-------------------Constraints--------------------* @@ -269,11 +269,11 @@ zref(10) = 1.0 *----------------Physics Variables-----------------* - +molflow_plasma_fuelling = 5.0e+21 alphan = 1.00 * density profile index alphat = 1.45 * temperature profile index +i_alphaj = 1 * Current profile index selection (0 = user input) aspect = 3.0 * aspect ratio (`iteration variable 1`) -beta_total_vol_avg = 0.03230408815355488 * total plasma beta (`iteration variable 5`) (calculated if stellarator) b_plasma_toroidal_on_axis = 5.318322174644904 * toroidal field on axis (T) (`iteration variable 2`) nd_plasma_electrons_vol_avg = 7.796223900029837e+19 * electron density (/m3) (`iteration variable 6`) beta_norm_max = 3.0 * Troyon-like coefficient for beta scaling calculated @@ -294,9 +294,9 @@ plasma_res_factor = 0.7 * plasma resistivity pre-factor radius_plasma_pedestal_density_norm = 0.94 * r/a of density pedestal (`i_plasma_pedestal==1`) radius_plasma_pedestal_temp_norm = 0.94 * r/a of temperature pedestal (`i_plasma_pedestal==1`) tbeta = 2.0 * temperature profile index beta (`i_plasma_pedestal==1) -temp_plasma_separatrix_electron_kev = 5.5 * electron temperature of pedestal (keV) (`i_plasma_pedestal==1') +temp_plasma_pedestal_electron_kev = 5.5 * electron temperature of pedestal (keV) (`i_plasma_pedestal==1') temp_plasma_separatrix_electron_kev = 0.1 * electron temperature at separatrix (keV) (`i_plasma_pedestal==1`) calculated if reinke -i_confinement_time = 34 * switch for energy confinement time scaling law (see description in `tauscl`) +i_confinement_time = 34 * switch for energy confinement time scaling law (see description in `labels_confinement_scalings`) i_plasma_geometry = 0 * switch for plasma cross-sectional shape calculation; kappa = 1.85 * plasma separatrix elongation (calculated if `i_plasma_geometry = 1-5; 7 or 9-10`) q95 = 3.7339078193128556 * Safety factor 'near' plasma edge (`iteration variable 18`) equal to q95 @@ -307,6 +307,7 @@ i_single_null = 1 * switch for single null / double null plasma; f_sync_reflect = 0.6 * synchrotron wall reflectivity factor temp_plasma_electron_vol_avg_kev = 12.221383528378944 * volume averaged electron temperature (keV) (`iteration variable 4`) triang = 0.5 * plasma separatrix triangularity (calculated if `i_plasma_geometry = 1; 3-5 or 7`) +i_ind_plasma_internal_norm = 1 * Normalised plasma internal induction selection switch (0 = user input) *----------------------Power-----------------------* @@ -375,13 +376,15 @@ t_tf_superconductor_quench = 17.97282589344206 * fast discharge time for TF co n_tf_coils = 16 * Number of TF coils (default = 50 for stellarators); Number of TF coils outer legs for ST tftmp = 4.75 * peak helium coolant temperature in TF coils and PF coils (K) dr_tf_nose_case = 0.2816873221155309 * inboard TF coil case outer (non-plasma side) thickness (m) (`iteration variable 57`) -dr_tf_wp_with_insulation = 0.5153787768966674 * radial thickness of winding pack (m) (`iteration variable 140`) (issue #514) +dr_tf_wp_with_insulation = 0.5153787768966674 * radial thickness of winding pack (m) (issue #514) dx_tf_turn_steel = 0.008012110032981922 * TF coil conduit case thickness (m) (`iteration variable 58`) dx_tf_wp_insulation = 0.008 * Thickness of the ground insulation layer surrounding (m) temp_cs_superconductor_margin_min = 1.5 * minimum allowable temperature margin ; CS (K) tmargmin = 1.5 * minimum allowable temperature margin ; TFC AND CS (K) v_tf_coil_dump_quench_max_kv = 10.0 * max voltage across TF coil during quench (kV) (`iteration variable 52`) f_a_tf_turn_cable_space_extra_void = 0.3 * coolant fraction of TFC 'cable' (`i_tf_sup=1`); or of TFC leg (`i_tf_sup=0`) +flu_tf_neutron_fast_max = 1e22 * max fluence in the TF coil + *-----------------Times Variables------------------* diff --git a/examples/data/large_tokamak_varied_min_net_electric_IN.DAT b/examples/data/large_tokamak_varied_min_net_electric_IN.DAT index f60f49cdb1..b913f7347d 100644 --- a/examples/data/large_tokamak_varied_min_net_electric_IN.DAT +++ b/examples/data/large_tokamak_varied_min_net_electric_IN.DAT @@ -17,17 +17,11 @@ i_figure_merit = 1 * Error tolerance for VMCON epsvmc = 1e-7 -n_equality_constraints = 3 +n_equality_constraints = 1 * Constraint Equations - Consistency Equations * ************************************************ -* Beta consistency * -*------------------* -icc = 1 -ixc = 5 * beta_total_vol_avg -beta_total_vol_avg = 0.03 - * Global power balance * *----------------------* icc = 2 @@ -36,6 +30,10 @@ icc = 2 *--------------------------* icc = 11 +* Fuel-ion equilibrium * +*----------------------* +icc = 93 + * Constraint Equations - Limit Equations * ****************************************** diff --git a/examples/data/large_tokamak_varyrun_IN.DAT b/examples/data/large_tokamak_varyrun_IN.DAT index 1b76e61add..cf0dbb269a 100644 --- a/examples/data/large_tokamak_varyrun_IN.DAT +++ b/examples/data/large_tokamak_varyrun_IN.DAT @@ -17,16 +17,11 @@ i_figure_merit = 1 * Error tolerance for VMCON epsvmc = 1e-7 -n_equality_constraints = 3 +n_equality_constraints = 2 * Constraint Equations - Consistency Equations * ************************************************ -* Beta consistency * -*------------------* -icc = 1 -ixc = 5 * beta_total_vol_avg -beta_total_vol_avg = 0.03 * Global power balance * *----------------------* @@ -36,6 +31,10 @@ icc = 2 *--------------------------* icc = 11 +*Fuel ion equilibrium * +*----------------------* +icc = 93 + * Constraint Equations - Limit Equations * ****************************************** diff --git a/process/core/input.py b/process/core/input.py index 31b6c6f6cd..475a3c0c4e 100644 --- a/process/core/input.py +++ b/process/core/input.py @@ -164,7 +164,6 @@ def bounds(self) -> tuple[NumberType | None, NumberType | None]: "alphat": InputVariable("physics", float, range=(0.0, 10.0)), "aspect": InputVariable("physics", float, range=(1.001, 40.0)), "beamfus0": InputVariable("physics", float, range=(0.01, 10.0)), - "beta_total_vol_avg": InputVariable("physics", float, range=(0.0, 1.0)), "beta_vol_avg_max": InputVariable("physics", float, range=(0.0, 1.0)), "beta_vol_avg_min": InputVariable("physics", float, range=(0.0, 1.0)), "betbm0": InputVariable("physics", float, range=(0.0, 10.0)), @@ -1186,6 +1185,7 @@ def bounds(self) -> tuple[NumberType | None, NumberType | None]: "f_len_sol_power_decay_inboard_outboard": InputVariable( "physics", float, range=(0.01, 2.0) ), + "molflow_plasma_fuelling": InputVariable("physics", float, range=(1.0e20, 1.0e23)), } diff --git a/process/core/io/plot/summary.py b/process/core/io/plot/summary.py index 6d1379d3d5..6c3ad4181e 100644 --- a/process/core/io/plot/summary.py +++ b/process/core/io/plot/summary.py @@ -2982,7 +2982,7 @@ def plot_main_plasma_information( f" - Average mass of all plasma ions: {mfile.get('m_ions_total_amu', scan=scan):.3f} amu\n" f"Fuel mass: {mfile.get('m_plasma_fuel_ions', scan=scan) * 1000:.4f} g\n" f" - Average mass of all fuel ions: {mfile.get('m_fuel_amu', scan=scan):.3f} amu\n\n" - f"Fueling rate: {mfile.get('molflow_plasma_fuelling_required', scan=scan):.3e} nucleus-pairs/s\n" + f"Fueling rate: {mfile.get('molflow_plasma_fuelling', scan=scan):.3e} nucleus-pairs/s\n" f"Fuel burn-up rate: {mfile.get('rndfuel', scan=scan):.3e} reactions/s \n" f"Burn-up fraction: {mfile.get('burnup', scan=scan):.4f} \n" ) diff --git a/process/core/solver/constraints.py b/process/core/solver/constraints.py index c5d786d0f8..08792a461e 100644 --- a/process/core/solver/constraints.py +++ b/process/core/solver/constraints.py @@ -1981,6 +1981,20 @@ def constraint_equation_92(constraint_registration, data): ) +@ConstraintManager.register_constraint(93, "", "=") +def constraint_equation_93(constraint_registration, data): + """Fuel-ion equilibrium constraint. + + Requires the plasma fuelling rate to equal the fuelling rate required + for fuel-ion equilibrium. + """ + return eq( + data.physics.molflow_plasma_fuelling, + data.physics.molflow_plasma_fuelling_equilibrium, + constraint_registration, + ) + + def constraint_eqns(m: int, ieqn: int, data: DataStructure): """Evaluates the constraints given the current state of PROCESS. diff --git a/process/core/solver/iteration_variables.py b/process/core/solver/iteration_variables.py index fd468c5537..745aad46e8 100644 --- a/process/core/solver/iteration_variables.py +++ b/process/core/solver/iteration_variables.py @@ -46,7 +46,6 @@ class IterationVariable: 2: IterationVariable("b_plasma_toroidal_on_axis", "physics", 0.010, 30.00), 3: IterationVariable("rmajor", "physics", 0.1, 50.00), 4: IterationVariable("temp_plasma_electron_vol_avg_kev", "physics", 5.0, 150.0), - 5: IterationVariable("beta_total_vol_avg", "physics", 0.001, 1.0), 6: IterationVariable("nd_plasma_electrons_vol_avg", "physics", 2.0e19, 1.0e21), 7: IterationVariable("f_nd_beam_electron", "physics", 1.0e-6, 1.0), 10: IterationVariable("hfact", "physics", 0.1, 3.0), @@ -239,6 +238,12 @@ class IterationVariable: 175: IterationVariable("kappa", "physics", 0.00, 10.00), 176: IterationVariable("f_st_coil_aspect", "stellarator", 0.70, 1.30), 177: IterationVariable("f_a_tf_turn_cable_space_extra_void", "tfcoil", 0.01, 1.0), + 178: IterationVariable( + "molflow_plasma_fuelling", + "physics", + 1.0e20, + 1.0e23, + ), } diff --git a/process/data_structure/numerics.py b/process/data_structure/numerics.py index 85785576a4..b34af3cf6b 100644 --- a/process/data_structure/numerics.py +++ b/process/data_structure/numerics.py @@ -292,14 +292,15 @@ class NumericsData: "CS achievable stress load cycles lower limit ", "ECRH ignitability ", # Stellarator constraint "Fuel composition consistency ", + "Fuel ion equilibrium ", ] ) """Labels describing constraint equations (corresponding itvs)