device.dat):
+ - 0 use tokamak model;
+ - 1 use stellarator model
+ fbeta_max: f-value for beta limit
+ beta_max: allowable beta
+ beta: total plasma beta (calculated if ipedestal =3)
+ beta_fast_alpha: fast alpha beta component
+ beta_beam: neutral beam beta component
+ bt: toroidal field
+ btot: total field
+ """
+ # Include all beta components: relevant for both tokamaks and stellarators
+ if (
+ fortran.physics_variables.i_beta_component == 0
+ or fortran.stellarator_variables.istell != 0
+ ):
+ cc = (
+ fortran.physics_variables.beta / fortran.physics_variables.beta_max
+ - 1.0 * fortran.constraint_variables.fbeta_max
+ )
+ con = fortran.physics_variables.beta_max
+ err = (
+ fortran.physics_variables.beta_max
+ - fortran.physics_variables.beta / fortran.constraint_variables.fbeta_max
+ )
+ # Here, the beta limit applies to only the thermal component, not the fast alpha or neutral beam parts
+ elif fortran.physics_variables.i_beta_component == 1:
+ cc = (
+ (
+ fortran.physics_variables.beta
+ - fortran.physics_variables.beta_fast_alpha
+ - fortran.physics_variables.beta_beam
+ )
+ / fortran.physics_variables.beta_max
+ - 1.0 * fortran.constraint_variables.fbeta_max
+ )
+ con = fortran.physics_variables.beta_max
+ err = (
+ fortran.physics_variables.beta_max
+ - (
+ fortran.physics_variables.beta
+ - fortran.physics_variables.beta_fast_alpha
+ - fortran.physics_variables.beta_beam
+ )
+ / fortran.constraint_variables.fbeta_max
+ )
+ # Beta limit applies to thermal + neutral beam: components of the total beta, i.e. excludes alphas
+ elif fortran.physics_variables.i_beta_component == 2:
+ cc = (
+ (fortran.physics_variables.beta - fortran.physics_variables.beta_fast_alpha)
+ / fortran.physics_variables.beta_max
+ - 1.0 * fortran.constraint_variables.fbeta_max
+ )
+ con = fortran.physics_variables.beta_max * (1.0 - cc)
+ err = (
+ fortran.physics_variables.beta - fortran.physics_variables.beta_fast_alpha
+ ) * cc
+ # Beta limit applies to toroidal beta
+ elif fortran.physics_variables.i_beta_component == 3:
+ cc = (
+ (
+ fortran.physics_variables.beta
+ * (fortran.physics_variables.btot / fortran.physics_variables.bt) ** 2
+ )
+ / fortran.physics_variables.beta_max
+ - 1.0 * fortran.constraint_variables.fbeta_max
+ )
+ con = fortran.physics_variables.beta_max
+ err = (
+ fortran.physics_variables.beta_max
+ - (
+ fortran.physics_variables.beta
+ * (fortran.physics_variables.btot / fortran.physics_variables.bt) ** 2
+ )
+ / fortran.constraint_variables.fbeta_max
+ )
+
+ return ConstraintResult(cc, con, err)
+
+
+@ConstraintManager.register_constraint(25, "T", "<=")
+def constraint_equation_25():
+ """Equation for peak toroidal field upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fpeakb: f-value for maximum toroidal field
+ bmxlim: maximum peak toroidal field (T)
+ b_tf_inboard_peak: mean peak field at TF coil (T)
+ """
+ cc = (
+ fortran.tfcoil_variables.b_tf_inboard_peak / fortran.constraint_variables.bmxlim
+ - 1.0 * fortran.constraint_variables.fpeakb
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.bmxlim * (1.0 - cc),
+ fortran.tfcoil_variables.b_tf_inboard_peak * cc,
+ )
+
+
+@ConstraintManager.register_constraint(26, "A/m2", "<=")
+def constraint_equation_26():
+ """Equation for Central Solenoid current density upper limit at EOF
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fjohcreal: f-value for central solenoid current at end-of-flattop
+ j_cs_critical_flat_top_end: allowable central solenoid current density at end of flat-top (A/m2)
+ j_cs_flat_top_end: central solenoid overall current density at end of flat-top (A/m2)
+ """
+ return ConstraintResult(
+ fortran.pfcoil_variables.j_cs_flat_top_end
+ / fortran.pfcoil_variables.j_cs_critical_flat_top_end
+ - 1.0 * fortran.constraint_variables.fjohc,
+ fortran.pfcoil_variables.j_cs_critical_flat_top_end,
+ fortran.pfcoil_variables.j_cs_critical_flat_top_end
+ - fortran.pfcoil_variables.j_cs_flat_top_end
+ / fortran.constraint_variables.fjohc,
+ )
+
+
+@ConstraintManager.register_constraint(27, "A/m2", "<=")
+def constraint_equation_27():
+ """Equation for Central Solenoid current density upper limit at BOP
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fjohc0: f-value for central solenoid current at beginning of pulse
+ j_cs_critical_pulse_start: allowable central solenoid current density at beginning of pulse (A/m2)
+ j_cs_pulse_start: central solenoid overall current density at beginning of pulse (A/m2)
+ """
+ return ConstraintResult(
+ fortran.pfcoil_variables.j_cs_pulse_start
+ / fortran.pfcoil_variables.j_cs_critical_pulse_start
+ - 1.0 * fortran.constraint_variables.fjohc0,
+ fortran.pfcoil_variables.j_cs_critical_pulse_start,
+ fortran.pfcoil_variables.j_cs_critical_pulse_start
+ - fortran.pfcoil_variables.j_cs_pulse_start
+ / fortran.constraint_variables.fjohc0,
+ )
+
+
+@ConstraintManager.register_constraint(28, "", ">=")
+def constraint_equation_28():
+ """Equation for fusion gain (big Q) lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fqval: pf-value for Q
+ bigq: Fusion gain; P_fusion / (P_injection + P_ohmic)
+ bigqmin: minimum fusion gain Q
+ i_plasma_ignited : input integer : switch for ignition assumption:
+ - 0 do not assume plasma ignition;
+ - 1 assume ignited (but include auxiliary power in costs)
+ Obviously, ignite must be zero if current drive is required.
+ If i_plasma_ignited=1, any auxiliary power is assumed to be used only
+ during plasma start-up, and is excluded from all steady-state
+ power balance calculations.
+ """
+ if fortran.physics_variables.i_plasma_ignited != 0:
+ raise ProcessValueError("Do not use constraint 28 if i_plasma_ignited=1")
+
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fqval
+ * fortran.current_drive_variables.bigq
+ / fortran.constraint_variables.bigqmin
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.bigqmin * (1.0 - cc),
+ fortran.constraint_variables.bigqmin * cc,
+ )
+
+
+@ConstraintManager.register_constraint(29, "m", "=")
+def constraint_equation_29():
+ """Equation for inboard major radius: This is a consistency equation
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ rmajor: plasma major radius (m) (iteration variable 3)
+ rminor: plasma minor radius (m)
+ rinboard: plasma inboard radius (m)
+ """
+ cc = (
+ 1.0
+ - (fortran.physics_variables.rmajor - fortran.physics_variables.rminor)
+ / fortran.build_variables.rinboard
+ )
+ return ConstraintResult(
+ cc,
+ fortran.build_variables.rinboard * (1.0 - cc),
+ fortran.build_variables.rinboard * cc,
+ )
+
+
+@ConstraintManager.register_constraint(30, "MW", "<=")
+def constraint_equation_30():
+ """Equation for injection power upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ p_hcd_injected_total_mw: total auxiliary injected power (MW)
+ fpinj: f-value for injection power
+ p_hcd_injected_max: Maximum allowable value for injected power (MW)
+ """
+ return ConstraintResult(
+ fortran.current_drive_variables.p_hcd_injected_total_mw
+ / fortran.current_drive_variables.p_hcd_injected_max
+ - 1.0 * fortran.constraint_variables.fpinj,
+ fortran.current_drive_variables.p_hcd_injected_max,
+ fortran.current_drive_variables.p_hcd_injected_max
+ - fortran.current_drive_variables.p_hcd_injected_total_mw
+ / fortran.constraint_variables.fpinj,
+ )
+
+
+@ConstraintManager.register_constraint(31, "Pa", "<=")
+def constraint_equation_31():
+ """Equation for TF coil case stress upper limit (SCTF)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fstrcase: f-value for TF coil case stress
+ sig_tf_case_max: Allowable maximum shear stress in TF coil case (Tresca criterion) (Pa)
+ sig_tf_case: Constrained stress in TF coil case (Pa)
+ """
+ return ConstraintResult(
+ fortran.tfcoil_variables.sig_tf_case / fortran.tfcoil_variables.sig_tf_case_max
+ - 1.0 * fortran.constraint_variables.fstrcase,
+ fortran.tfcoil_variables.sig_tf_case_max,
+ fortran.tfcoil_variables.sig_tf_case_max
+ - fortran.tfcoil_variables.sig_tf_case / fortran.constraint_variables.fstrcase,
+ )
+
+
+@ConstraintManager.register_constraint(32, "Pa", "<=")
+def constraint_equation_32():
+ """Equation for TF coil conduit stress upper limit (SCTF)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fstrcond: f-value for TF coil conduit stress
+ sig_tf_wp_max: Allowable maximum shear stress in TF coil conduit (Tresca criterion) (Pa)
+ sig_tf_wp: Constrained stress in TF conductor conduit (Pa)
+ """
+ return ConstraintResult(
+ fortran.tfcoil_variables.sig_tf_wp / fortran.tfcoil_variables.sig_tf_wp_max
+ - 1.0 * fortran.constraint_variables.fstrcond,
+ fortran.tfcoil_variables.sig_tf_wp_max,
+ fortran.tfcoil_variables.sig_tf_wp_max
+ - fortran.tfcoil_variables.sig_tf_wp / fortran.constraint_variables.fstrcond,
+ )
+
+
+@ConstraintManager.register_constraint(33, "A/m2", "<=")
+def constraint_equation_33():
+ """Equation for TF coil operating/critical J upper limit (SCTF)
+ author: P B Lloyd, CCFE, Culham Science Centre
+ args : output structure : residual error; constraint value;
+
+ fiooic: f-value for TF coil operating current / critical
+ jwdgcrt: critical current density for winding pack (A/m2)
+ j_tf_wp: winding pack current density (A/m2)
+ """
+ if fortran.constraint_variables.fiooic > 0.7:
+ WarningManager.create_warning(
+ "fiooic shouldn't be above 0.7 for engineering reliability"
+ )
+
+ cc = (
+ fortran.tfcoil_variables.j_tf_wp / fortran.tfcoil_variables.jwdgcrt
+ - 1.0 * fortran.constraint_variables.fiooic
+ )
+ return ConstraintResult(
+ cc,
+ fortran.tfcoil_variables.jwdgcrt * (1.0 - cc),
+ fortran.tfcoil_variables.j_tf_wp * cc,
+ )
+
+
+@ConstraintManager.register_constraint(34, "V", "<=")
+def constraint_equation_34():
+ """Equation for TF coil dump voltage upper limit (SCTF)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fvdump: f-value for dump voltage
+ vdalw: max voltage across TF coil during quench (kV)
+ vtfskv: voltage across a TF coil during quench (kV)
+ """
+ return ConstraintResult(
+ fortran.tfcoil_variables.vtfskv / fortran.tfcoil_variables.vdalw
+ - 1.0 * fortran.constraint_variables.fvdump,
+ fortran.tfcoil_variables.vdalw,
+ fortran.tfcoil_variables.vdalw - fortran.tfcoil_variables.vtfskv,
+ )
+
+
+@ConstraintManager.register_constraint(35, "A/m2", "<=")
+def constraint_equation_35():
+ """Equation for TF coil J_wp/J_prot upper limit (SCTF)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fjprot: f-value for TF coil winding pack current density
+ jwdgpro: allowable TF coil winding pack current density, for dump temperature
+ rise protection (A/m2)
+ j_tf_wp: winding pack current density (A/m2)
+ """
+ return ConstraintResult(
+ fortran.tfcoil_variables.j_tf_wp / fortran.tfcoil_variables.jwdgpro
+ - 1.0 * fortran.constraint_variables.fjprot,
+ fortran.tfcoil_variables.jwdgpro,
+ fortran.tfcoil_variables.j_tf_wp - fortran.tfcoil_variables.jwdgpro,
+ )
+
+
+@ConstraintManager.register_constraint(36, "K", ">=")
+def constraint_equation_36():
+ """Equation for TF coil s/c temperature margin lower limit (SCTF)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ ftmargtf: f-value for TF coil temperature margin
+ tmargtf: TF coil temperature margin (K)
+ tmargmin_tf: minimum allowable temperature margin : TF coils (K)
+ """
+ return ConstraintResult(
+ 1.0
+ - fortran.constraint_variables.ftmargtf
+ * fortran.tfcoil_variables.tmargtf
+ / fortran.tfcoil_variables.tmargmin_tf,
+ fortran.tfcoil_variables.tmargmin_tf,
+ fortran.tfcoil_variables.tmargmin_tf - fortran.tfcoil_variables.tmargtf,
+ )
+
+
+@ConstraintManager.register_constraint(37, "1E20 A/Wm2", "<=")
+def constraint_equation_37():
+ """Equation for current drive gamma upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fgamcd: f-value for current drive gamma
+ gammax: maximum current drive gamma
+ eta_cd_norm_hcd_primary: normalised current drive efficiency (1.0e20 A/W-m2)
+ """
+ cc = (
+ fortran.current_drive_variables.eta_cd_norm_hcd_primary
+ / fortran.constraint_variables.gammax
+ - 1.0 * fortran.constraint_variables.fgamcd
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.gammax * (1.0 - cc),
+ fortran.current_drive_variables.eta_cd_norm_hcd_primary * cc,
+ )
+
+
+@ConstraintManager.register_constraint(39, "K", "<=")
+def constraint_equation_39():
+ """Equation for first wall temperature upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ ftpeak: f-value for first wall peak temperature
+ temp_fw_max: maximum temperature of first wall material (K) (i_thermal_electric_conversion>1)
+ temp_fw_peak: peak first wall temperature (K)
+ """
+ if fortran.fwbs_variables.temp_fw_peak < 1.0:
+ raise ProcessValueError(
+ "temp_fw_peak = 0 implies i_pulsed_plant=0; do not use constraint 39 if i_pulsed_plant=0"
+ )
+ cc = (
+ fortran.fwbs_variables.temp_fw_peak / fortran.fwbs_variables.temp_fw_max
+ - 1.0 * fortran.constraint_variables.ftpeak
+ )
+ return ConstraintResult(
+ cc,
+ fortran.fwbs_variables.temp_fw_max * (1.0 - cc),
+ fortran.fwbs_variables.temp_fw_peak * cc,
+ )
+
+
+@ConstraintManager.register_constraint(40, "MW", ">=")
+def constraint_equation_40():
+ """Equation for auxiliary power lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fauxmn: f-value for minimum auxiliary power
+ p_hcd_injected_total_mw: total auxiliary injected power (MW)
+ auxmin: minimum auxiliary power (MW)
+ """
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fauxmn
+ * fortran.current_drive_variables.p_hcd_injected_total_mw
+ / fortran.constraint_variables.auxmin
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.auxmin * (1.0 - cc),
+ fortran.constraint_variables.auxmin * cc,
+ )
+
+
+@ConstraintManager.register_constraint(41, "sec", ">=")
+def constraint_equation_41():
+ """Equation for plasma current ramp-up time lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ ft_current_ramp_up: f-value for plasma current ramp-up time
+ t_current_ramp_up: plasma current ramp-up time for current initiation (s)
+ t_current_ramp_up_min: minimum plasma current ramp-up time (s)
+ """
+ cc = (
+ 1.0
+ - fortran.constraint_variables.ft_current_ramp_up
+ * fortran.times_variables.t_current_ramp_up
+ / fortran.constraint_variables.t_current_ramp_up_min
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.t_current_ramp_up_min * (1.0 - cc),
+ fortran.constraint_variables.t_current_ramp_up_min * cc,
+ )
+
+
+@ConstraintManager.register_constraint(42, "sec", ">=")
+def constraint_equation_42():
+ """Equation for cycle time lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ ftcycl: f-value for cycle time
+ t_cycle: full cycle time (s)
+ tcycmn: minimum cycle time (s)
+ """
+ if fortran.constraint_variables.tcycmn < 1.0:
+ raise ProcessValueError(
+ "tcycmn = 0 implies that i_pulsed_plant=0; do not use constraint 42 if i_pulsed_plant=0"
+ )
+
+ cc = (
+ 1.0
+ - fortran.constraint_variables.ftcycl
+ * fortran.times_variables.t_cycle
+ / fortran.constraint_variables.tcycmn
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.tcycmn * (1.0 - cc),
+ fortran.constraint_variables.tcycmn * cc,
+ )
+
+
+@ConstraintManager.register_constraint(43, "deg C", "=")
+def constraint_equation_43():
+ """Equation for average centrepost temperature: This is a consistency equation (TART)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ temp_cp_average: average temp of TF coil inboard leg conductor (C)e
+ tcpav2: centrepost average temperature (C) (for consistency)
+ itart: switch for spherical tokamak (ST) models:
+ - 0 use conventional aspect ratio models;
+ - 1 use spherical tokamak models
+ """
+ if fortran.physics_variables.itart == 0:
+ raise ProcessValueError("Do not use constraint 43 if itart=0")
+
+ if fortran.tfcoil_variables.i_tf_sup == 0:
+ temp_cp_average = fortran.tfcoil_variables.temp_cp_average - 273.15
+ tcpav2 = fortran.tfcoil_variables.tcpav2 - 273.15
+ else:
+ temp_cp_average = fortran.tfcoil_variables.temp_cp_average
+ tcpav2 = fortran.tfcoil_variables.tcpav2
+
+ cc = 1.0 - temp_cp_average / tcpav2
+
+ return ConstraintResult(cc, tcpav2 * (1.0 - cc), tcpav2 * cc)
+
+
+@ConstraintManager.register_constraint(44, "deg C", "<=")
+def constraint_equation_44():
+ """Equation for centrepost temperature upper limit (TART)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fptemp: f-value for peak centrepost temperature
+ ptempalw: maximum peak centrepost temperature (K)
+ tcpmax: peak centrepost temperature (K)
+ itart: switch for spherical tokamak (ST) models:
+ - 0: use conventional aspect ratio models;
+ - 1: use spherical tokamak models
+ """
+ if fortran.physics_variables.itart == 0:
+ raise ProcessValueError("Do not use constraint 44 if itart=0")
+
+ if fortran.tfcoil_variables.i_tf_sup == 0: # ! Copper case
+ ptempalw = fortran.tfcoil_variables.ptempalw - 273.15
+ tcpmax = fortran.tfcoil_variables.tcpmax - 273.15
+ else:
+ ptempalw = fortran.tfcoil_variables.ptempalw
+ tcpmax = fortran.tfcoil_variables.tcpmax
+
+ cc = tcpmax / ptempalw - 1.0 * fortran.constraint_variables.fptemp
+ return ConstraintResult(cc, ptempalw * (1.0 - cc), tcpmax * cc)
+
+
+@ConstraintManager.register_constraint(45, "", ">=")
+def constraint_manager_45():
+ """Equation for edge safety factor lower limit (TART)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fq: f-value for edge safety factor
+ q95 : safety factor 'near' plasma edge
+ (unless i_plasma_current = 2 (ST current scaling), in which case q = mean edge safety factor qbar)
+ q95_min: lower limit for edge safety factor
+ itart : input integer : switch for spherical tokamak (ST) models:
+ - 0 use conventional aspect ratio models;
+ - 1 use spherical tokamak models"""
+ if fortran.physics_variables.itart == 0:
+ raise ProcessValueError("Do not use constraint 45 if itart=0")
+
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fq
+ * fortran.physics_variables.q95
+ / fortran.physics_variables.q95_min
+ )
+ return ConstraintResult(
+ cc,
+ fortran.physics_variables.q95_min * (1.0 - cc),
+ fortran.physics_variables.q95_min * cc,
+ )
+
+
+@ConstraintManager.register_constraint(46, "", "<=")
+def constraint_equation_46():
+ """Equation for Ip/Irod upper limit (TART)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ eps: inverse aspect ratio
+ fipir: f-value for Ip/Irod upper limit
+ c_tf_total: total (summed) current in TF coils (A)
+ plasma_current: plasma current (A)
+ itart: switch for spherical tokamak (ST) models:
+ - 0: use conventional aspect ratio models;
+ - 1: use spherical tokamak models
+ """
+ if fortran.physics_variables.itart == 0:
+ raise ProcessValueError("Do not use constraint 46 if itart=0")
+
+ # maximum ratio of plasma current to centrepost current
+ cratmx = 1.0 + 4.91 * (fortran.physics_variables.eps - 0.62)
+ cc = (
+ fortran.physics_variables.plasma_current / fortran.tfcoil_variables.c_tf_total
+ ) / cratmx - 1.0 * fortran.constraint_variables.fipir
+
+ return ConstraintResult(
+ cc,
+ cratmx * (1.0 - cc),
+ fortran.physics_variables.plasma_current
+ / fortran.tfcoil_variables.c_tf_total
+ * cc,
+ )
+
+
+@ConstraintManager.register_constraint(48, "", "<=")
+def constraint_equation_48():
+ """Equation for poloidal beta upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fbeta_poloidal: rf-value for poloidal beta
+ beta_poloidal_max: maximum poloidal beta
+ beta_poloidal: poloidal beta
+ """
+ cc = (
+ fortran.physics_variables.beta_poloidal
+ / fortran.constraint_variables.beta_poloidal_max
+ - 1.0 * fortran.constraint_variables.fbeta_poloidal
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.beta_poloidal_max * (1.0 - cc),
+ fortran.physics_variables.beta_poloidal * cc,
+ )
+
+
+@ConstraintManager.register_constraint(50, "Hz", "<=")
+def constraint_equation_50():
+ """IFE option: Equation for repetition rate upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+ author: S I Muldrew, CCFE, Culham Science Centre
+ """
+ cc = (
+ fortran.ife_variables.reprat / fortran.ife_variables.rrmax
+ - 1.0 * fortran.ife_variables.frrmax
+ )
+ return ConstraintResult(
+ cc, fortran.ife_variables.rrmax * (1.0 - cc), fortran.ife_variables.reprat * cc
+ )
+
+
+@ConstraintManager.register_constraint(51, "V.s", "=")
+def constraint_equation_51():
+ """Equation to enforce startup flux = available startup flux
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ vs_plasma_res_ramp: resistive losses in startup V-s (Wb)
+ vs_plasma_ind_ramp: internal and external plasma inductance V-s (Wb))
+ vs_cs_pf_total_ramp: total flux swing for startup (Wb)
+ """
+ cc = 1.0 - fortran.pfcoil_variables.fvs_cs_pf_total_ramp * abs(
+ (
+ fortran.physics_variables.vs_plasma_res_ramp
+ + fortran.physics_variables.vs_plasma_ind_ramp
+ )
+ / fortran.pfcoil_variables.vs_cs_pf_total_ramp
+ )
+ return ConstraintResult(
+ cc,
+ fortran.pfcoil_variables.vs_cs_pf_total_ramp * (1.0 - cc),
+ fortran.pfcoil_variables.vs_cs_pf_total_ramp * cc,
+ )
+
+
+@ConstraintManager.register_constraint(52, "", ">=")
+def constraint_equation_52():
+ """Equation for tritium breeding ratio lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ ftbr: f-value for minimum tritium breeding ratio
+ tbr: tritium breeding ratio (i_blanket_type=2,3 (KIT HCPB/HCLL))
+ tbrmin: minimum tritium breeding ratio (If i_blanket_type=1, tbrmin=minimum 5-year time-averaged tritium breeding ratio)
+ """
+ cc = (
+ 1.0
+ - fortran.constraint_variables.ftbr
+ * fortran.fwbs_variables.tbr
+ / fortran.constraint_variables.tbrmin
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.tbrmin * (1.0 - cc),
+ fortran.constraint_variables.tbrmin * cc,
+ )
+
+
+@ConstraintManager.register_constraint(53, "neutron/m2", "<=")
+def constraint_equation_53():
+ """Equation for fast neutron fluence on TF coil upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fflutf: f-value for maximum TF coil nuclear heating
+ nflutfmax: max fast neutron fluence on TF coil (n/m2)
+ nflutf: peak fast neutron fluence on TF coil superconductor (n/m2)
+ """
+ cc = (
+ fortran.fwbs_variables.nflutf / fortran.constraint_variables.nflutfmax
+ - 1.0 * fortran.constraint_variables.fflutf
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.nflutfmax * (1.0 - cc),
+ fortran.fwbs_variables.nflutf * cc,
+ )
+
+
+@ConstraintManager.register_constraint(54, "MW/m3", "<=")
+def constraint_equation_54():
+ """Equation for peak TF coil nuclear heating upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fptfnuc: f-value for maximum TF coil nuclear heating
+ ptfnucmax: maximum nuclear heating in TF coil (MW/m3)
+ ptfnucpm3: nuclear heating in the TF coil (MW/m3) (blktmodel>0)
+ """
+ cc = (
+ fortran.fwbs_variables.ptfnucpm3 / fortran.constraint_variables.ptfnucmax
+ - 1.0 * fortran.constraint_variables.fptfnuc
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.ptfnucmax * (1.0 - cc),
+ fortran.fwbs_variables.ptfnucpm3 * cc,
+ )
+
+
+@ConstraintManager.register_constraint(56, "MW/m", "<=")
+def constraint_equation_56():
+ """Equation for power through separatrix / major radius upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fpsepr: f-value for maximum Psep/R limit
+ pseprmax: maximum ratio of power crossing the separatrix to plasma major radius (Psep/R) (MW/m)
+ p_plasma_separatrix_mw: power to be conducted to the divertor region (MW)
+ rmajor: plasma major radius (m)
+ """
+ cc = (
+ (
+ fortran.physics_variables.p_plasma_separatrix_mw
+ / fortran.physics_variables.rmajor
+ )
+ / fortran.constraint_variables.pseprmax
+ - 1.0 * fortran.constraint_variables.fpsepr
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.pseprmax * (1.0 - cc),
+ (
+ fortran.physics_variables.p_plasma_separatrix_mw
+ / fortran.physics_variables.rmajor
+ )
+ * cc,
+ )
+
+
+@ConstraintManager.register_constraint(59, "", "<=")
+def constraint_equation_59():
+ """Equation for neutral beam shine-through fraction upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fnbshinef: f-value for maximum neutral beam shine-through fraction
+ f_p_beam_shine_through_max: maximum neutral beam shine-through fraction
+ f_p_beam_shine_through: neutral beam shine-through fraction
+ """
+ cc = (
+ fortran.current_drive_variables.f_p_beam_shine_through
+ / fortran.constraint_variables.f_p_beam_shine_through_max
+ - 1.0 * fortran.constraint_variables.fnbshinef
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.f_p_beam_shine_through_max * (1.0 - cc),
+ fortran.current_drive_variables.f_p_beam_shine_through * cc,
+ )
+
+
+@ConstraintManager.register_constraint(60, "K", ">=")
+def constraint_equation_60():
+ """Equation for Central Solenoid s/c temperature margin lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ ftmargoh: f-value for central solenoid temperature margin
+ temp_cs_margin: Central solenoid temperature margin (K)
+ tmargmin_cs: Minimum allowable temperature margin : CS (K)
+ """
+ return ConstraintResult(
+ 1.0
+ - fortran.constraint_variables.ftmargoh
+ * fortran.pfcoil_variables.temp_cs_margin
+ / fortran.tfcoil_variables.tmargmin_cs,
+ fortran.tfcoil_variables.tmargmin_cs,
+ fortran.tfcoil_variables.tmargmin_cs - fortran.pfcoil_variables.temp_cs_margin,
+ )
+
+
+@ConstraintManager.register_constraint(61, "", ">=")
+def constraint_equation_61():
+ """Equation for availability lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ favail: F-value for minimum availability
+ cfactr: Total plant availability fraction
+ avail_min: Minimum availability
+ """
+ cc = (
+ 1.0
+ - fortran.cost_variables.favail
+ * fortran.cost_variables.cfactr
+ / fortran.cost_variables.avail_min
+ )
+ return ConstraintResult(
+ cc,
+ fortran.cost_variables.avail_min * (1.0 - cc),
+ fortran.cost_variables.cfactr * cc,
+ )
+
+
+@ConstraintManager.register_constraint(62, "", ">=")
+def constraint_equation_62():
+ """Lower limit on f_alpha_energy_confinement the ratio of alpha particle to energy confinement times
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ falpha_energy_confinement: f-value for lower limit on f_alpha_energy_confinement the ratio of alpha particle to energy confinement
+ t_alpha_confinement: alpha particle confinement time (s)
+ t_energy_confinement: global thermal energy confinement time (sec)
+ f_alpha_energy_confinement_min: Lower limit on f_alpha_energy_confinement the ratio of alpha particle to energy confinement times
+ """
+ cc = (
+ 1.0
+ - fortran.constraint_variables.falpha_energy_confinement
+ * (
+ fortran.physics_variables.t_alpha_confinement
+ / fortran.physics_variables.t_energy_confinement
+ )
+ / fortran.constraint_variables.f_alpha_energy_confinement_min
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.f_alpha_energy_confinement_min,
+ (
+ fortran.physics_variables.t_alpha_confinement
+ / fortran.physics_variables.t_energy_confinement
+ )
+ * cc,
+ )
+
+
+@ConstraintManager.register_constraint(63, "", "<=")
+def constraint_equation_63():
+ """Upper limit on niterpump (vacuum_model = simple)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fniterpump: f-value for constraint that number of pumps < tfno
+ tfno: number of TF coils (default = 50 for stellarators)
+ niterpump: number of high vacuum pumps (real number), each with the throughput
+ """
+ cc = (
+ fortran.vacuum_variables.niterpump / fortran.tfcoil_variables.n_tf_coils
+ - 1.0 * fortran.constraint_variables.fniterpump
+ )
+ return ConstraintResult(
+ cc,
+ fortran.tfcoil_variables.n_tf_coils,
+ fortran.tfcoil_variables.n_tf_coils * cc,
+ )
+
+
+@ConstraintManager.register_constraint(64, "", "<=")
+def constraint_equation_64():
+ """Upper limit on Zeff
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fzeffmax: f-value for maximum zeff
+ zeffmax: maximum value for Zeff
+ zeff: plasma effective charge
+ """
+ cc = (
+ fortran.physics_variables.zeff / fortran.constraint_variables.fzeffmax
+ - 1.0 * fortran.constraint_variables.ffzeffmax
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.fzeffmax,
+ fortran.constraint_variables.fzeffmax * cc,
+ )
+
+
+@ConstraintManager.register_constraint(65, "Pa", "<=")
+def constraint_equation_65():
+ """Upper limit on stress of the vacuum vessel that occurs when the TF coil quenches.
+ author: Timothy Nunn, UKAEA
+
+ fmaxvvstress: f-value for constraint on maximum VV stress
+ max_vv_stress: Maximum permitted stress of the VV (Pa)
+ vv_stress_quench: Stress of the VV (Pa)
+ """
+ cc = (
+ fortran.sctfcoil_module.vv_stress_quench
+ / fortran.tfcoil_variables.max_vv_stress
+ - 1.0 * fortran.constraint_variables.fmaxvvstress
+ )
+ return ConstraintResult(
+ cc,
+ fortran.tfcoil_variables.max_vv_stress,
+ fortran.tfcoil_variables.max_vv_stress * cc,
+ )
+
+
+@ConstraintManager.register_constraint(66, "MW", "<=")
+def constrain_equation_66():
+ """Upper limit on rate of change of energy in poloidal field
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fpoloidalpower: f-value for constraint on rate of change of energy in poloidal field
+ maxpoloidalpower: Maximum permitted absolute rate of change of stored energy in poloidal field (MW)
+ peakpoloidalpower: Peak absolute rate of change of stored energy in poloidal field (MW) (11/01/16)
+ """
+ cc = (
+ fortran.pf_power_variables.peakpoloidalpower
+ / fortran.pf_power_variables.maxpoloidalpower
+ - 1.0 * fortran.constraint_variables.fpoloidalpower
+ )
+ return ConstraintResult(
+ cc,
+ fortran.pf_power_variables.maxpoloidalpower,
+ fortran.pf_power_variables.maxpoloidalpower * cc,
+ )
+
+
+@ConstraintManager.register_constraint(67, "MW/m2", "<=")
+def constraint_equation_67():
+ """Simple upper limit on radiation wall load
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fradwall: f-value for upper limit on radiation wall load
+ pflux_fw_rad_max: Maximum permitted radiation wall load (MW/m^2)
+ pflux_fw_rad_max_mw: Peak radiation wall load (MW/m^2)
+ """
+ cc = (
+ fortran.constraint_variables.pflux_fw_rad_max_mw
+ / fortran.constraint_variables.pflux_fw_rad_max
+ - 1.0 * fortran.constraint_variables.fradwall
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.pflux_fw_rad_max,
+ fortran.constraint_variables.pflux_fw_rad_max * cc,
+ )
+
+
+@ConstraintManager.register_constraint(68, "MWT/m", "<=")
+def constraint_equation_68():
+ """Upper limit on Psep scaling (PsepB/qAR)
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fpsepbqar: f-value for upper limit on psepbqar, maximum Psep*Bt/qAR limit
+ psepbqarmax: maximum permitted value of ratio of Psep*Bt/qAR (MWT/m)
+ p_plasma_separatrix_mw: Power to conducted to the divertor region (MW)
+ bt: toroidal field on axis (T) (iteration variable 2)
+ q95: safety factor q at 95% flux surface
+ aspect: aspect ratio (iteration variable 1)
+ rmajor: plasma major radius (m) (iteration variable 3)
+ i_q95_fixed: Switch that allows for fixing q95 only in this constraint.
+ q95_fixed: fixed safety factor q at 95% flux surface
+ """
+ if fortran.constraint_variables.i_q95_fixed == 1:
+ cc = (
+ (
+ (
+ fortran.physics_variables.p_plasma_separatrix_mw
+ * fortran.physics_variables.bt
+ )
+ / (
+ fortran.constraint_variables.q95_fixed
+ * fortran.physics_variables.aspect
+ * fortran.physics_variables.rmajor
+ )
+ )
+ / fortran.constraint_variables.psepbqarmax
+ - 1.0 * fortran.constraint_variables.fpsepbqar
+ )
+ err = (
+ fortran.physics_variables.p_plasma_separatrix_mw
+ * fortran.physics_variables.bt
+ ) / (
+ fortran.constraint_variables.q95_fixed
+ * fortran.physics_variables.aspect
+ * fortran.physics_variables.rmajor
+ ) - fortran.constraint_variables.psepbqarmax
+ else:
+ cc = (
+ (
+ (
+ fortran.physics_variables.p_plasma_separatrix_mw
+ * fortran.physics_variables.bt
+ )
+ / (
+ fortran.physics_variables.q95
+ * fortran.physics_variables.aspect
+ * fortran.physics_variables.rmajor
+ )
+ )
+ / fortran.constraint_variables.psepbqarmax
+ - 1.0 * fortran.constraint_variables.fpsepbqar
+ )
+ err = (
+ fortran.physics_variables.p_plasma_separatrix_mw
+ * fortran.physics_variables.bt
+ ) / (
+ fortran.physics_variables.q95
+ * fortran.physics_variables.aspect
+ * fortran.physics_variables.rmajor
+ ) - fortran.constraint_variables.psepbqarmax
+
+ return ConstraintResult(cc, fortran.constraint_variables.psepbqarmax, err)
+
+
+@ConstraintManager.register_constraint(72, "Pa", "<=")
+def constraint_equation_72():
+ """Upper limit on central Solenoid Tresca yield stress
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ In the case if the bucked and wedged option ( i_tf_bucking >= 2 ) the constrained
+ stress is the largest the largest stress of the
+ - CS stress at maximum current (conservative as the TF inward pressure is not taken
+ into account)
+ - CS stress at flux swing (no current in CS) from the TF inward pressure
+ This allow to cover the 2 worst stress scenario in the bucked and wedged design
+ Otherwise (free standing TF), the stress limits are only set by the CS stress at max current
+ Reverse the sign so it works as an inequality constraint (tmp_cc > 0)
+ This will have no effect if it is used as an equality constraint because it will be squared.
+
+ foh_stress: f-value for Tresca yield criterion in Central Solenoid
+ alstroh: allowable hoop stress in Central Solenoid structural material (Pa)
+ s_shear_cs_peak: Maximum shear stress coils/central solenoid (Pa)
+ sig_tf_cs_bucked: Maximum shear stress in CS case at flux swing (no current in CS)
+ can be significant for the bucked and weged design
+ i_tf_bucking: switch for TF structure design
+ """
+ # bucked and wedged desing
+ if (
+ fortran.tfcoil_variables.i_tf_bucking >= 2
+ and fortran.build_variables.i_tf_inside_cs == 0
+ ):
+ cc = (
+ max(
+ fortran.pfcoil_variables.s_shear_cs_peak,
+ fortran.tfcoil_variables.sig_tf_cs_bucked,
+ )
+ / fortran.pfcoil_variables.alstroh
+ - 1.0 * fortran.constraint_variables.foh_stress
+ )
+ err = fortran.pfcoil_variables.alstroh - max(
+ fortran.pfcoil_variables.s_shear_cs_peak,
+ fortran.tfcoil_variables.sig_tf_cs_bucked,
+ )
+ # Free standing CS
+ else:
+ cc = (
+ fortran.pfcoil_variables.s_shear_cs_peak / fortran.pfcoil_variables.alstroh
+ - 1.0 * fortran.constraint_variables.foh_stress
+ )
+ err = (
+ fortran.pfcoil_variables.alstroh - fortran.pfcoil_variables.s_shear_cs_peak
+ )
+
+ return ConstraintResult(cc, fortran.pfcoil_variables.alstroh, err)
+
+
+@ConstraintManager.register_constraint(73, "MW", ">=")
+def constraint_equation_73():
+ """Lower limit to ensure separatrix power is greater than the L-H power + auxiliary power
+ Related to constraint 15
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fplhsep: F-value for Psep >= Plh + Paux : for consistency of two values of separatrix power
+ p_l_h_threshold_mw: L-H mode power threshold (MW)
+ p_plasma_separatrix_mw: power to be conducted to the divertor region (MW)
+ p_hcd_injected_total_mw : inout real : total auxiliary injected power (MW)
+ """
+ cc = (
+ 1.0
+ - fortran.physics_variables.fplhsep
+ * fortran.physics_variables.p_plasma_separatrix_mw
+ / (
+ fortran.physics_variables.p_l_h_threshold_mw
+ + fortran.current_drive_variables.p_hcd_injected_total_mw
+ )
+ )
+ return ConstraintResult(
+ cc,
+ fortran.physics_variables.p_plasma_separatrix_mw,
+ fortran.physics_variables.p_plasma_separatrix_mw * cc,
+ )
+
+
+@ConstraintManager.register_constraint(74, "K", "<=")
+def constraint_equation_74():
+ """Upper limit to ensure TF coil quench temperature < tmax_croco
+ ONLY used for croco HTS coil
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fcqt: f-value: TF coil quench temparature remains below tmax_croco
+ croco_quench_temperature: CroCo strand: Actual temp reached during a quench (K)
+ tmax_croco: CroCo strand: maximum permitted temp during a quench (K)
+ """
+ cc = (
+ fortran.tfcoil_variables.croco_quench_temperature
+ / fortran.tfcoil_variables.tmax_croco
+ - 1.0 * fortran.constraint_variables.fcqt
+ )
+ return ConstraintResult(
+ cc,
+ fortran.tfcoil_variables.croco_quench_temperature,
+ fortran.tfcoil_variables.croco_quench_temperature * cc,
+ )
+
+
+@ConstraintManager.register_constraint(75, "A/m2", "<=")
+def constraint_equation_75():
+ """Upper limit to ensure that TF coil current / copper area < Maximum value
+ ONLY used for croco HTS coil
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ copperA_m2: TF coil current / copper area (A/m2)
+ copperA_m2_max: Maximum TF coil current / copper area (A/m2)
+ f_coppera_m2: f-value for TF coil current / copper area < copperA_m2_max
+ """
+ cc = (
+ fortran.rebco_variables.coppera_m2 / fortran.rebco_variables.coppera_m2_max
+ - 1.0 * fortran.rebco_variables.f_coppera_m2
+ )
+ return ConstraintResult(
+ cc, fortran.rebco_variables.coppera_m2, fortran.rebco_variables.coppera_m2 * cc
+ )
+
+
+@ConstraintManager.register_constraint(76, "m-3", "<=")
+def constraint_equation_76():
+ """Upper limit for Eich critical separatrix density model: Added for issue 558
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ Eich critical separatrix density model
+ Added for issue 558 with ref to http://iopscience.iop.org/article/10.1088/1741-4326/aaa340/pdf
+
+ alpha_crit: critical ballooning parameter value
+ nesep_crit: critical electron density at separatrix [m-3]
+ kappa: plasma separatrix elongation (calculated if i_plasma_geometry = 1-5, 7 or 9)
+ triang: plasma separatrix triangularity (calculated if i_plasma_geometry = 1, 3-5 or 7)
+ aspect: aspect ratio (iteration variable 1)
+ p_plasma_separatrix_mw: power to conducted to the divertor region (MW)
+ dlimit(7)array : density limit (/m3) as calculated using various models
+ fnesep: f-value for Eich critical separatrix density
+ """
+ # TODO: why on earth are these variables being set here!? Should they be local?
+ fortran.physics_variables.alpha_crit = (fortran.physics_variables.kappa**1.2) * (
+ 1.0 + 1.5 * fortran.physics_variables.triang
+ )
+ fortran.physics_variables.nesep_crit = (
+ 5.9
+ * fortran.physics_variables.alpha_crit
+ * (fortran.physics_variables.aspect ** (-2.0 / 7.0))
+ * (((1.0 + (fortran.physics_variables.kappa**2.0)) / 2.0) ** (-6.0 / 7.0))
+ * ((fortran.physics_variables.p_plasma_separatrix_mw * 1.0e6) ** (-11.0 / 70.0))
+ * fortran.physics_variables.dlimit[6]
+ )
+
+ cc = (
+ fortran.physics_variables.nesep / fortran.physics_variables.nesep_crit
+ - 1.0 * fortran.constraint_variables.fnesep
+ )
+ return ConstraintResult(
+ cc, fortran.physics_variables.nesep, fortran.physics_variables.nesep * cc
+ )
+
+
+@ConstraintManager.register_constraint(77, "A/turn", "<=")
+def constraint_equation_77():
+ """Equation for maximum TF current per turn upper limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fcpttf: f-value for TF coil current per turn
+ cpttf_max : allowable TF coil current per turn [A/turn]
+ cpttf : TF coil current per turn [A/turn]
+ """
+ cc = (
+ fortran.tfcoil_variables.cpttf / fortran.tfcoil_variables.cpttf_max
+ - 1.0 * fortran.constraint_variables.fcpttf
+ )
+ return ConstraintResult(
+ cc, fortran.tfcoil_variables.cpttf_max, fortran.tfcoil_variables.cpttf_max * cc
+ )
+
+
+@ConstraintManager.register_constraint(78, "", ">=")
+def constraint_equation_78():
+ """Equation for Reinke criterion, divertor impurity fraction lower limit
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ freinke : input : f-value for Reinke criterion (itv 147)
+ fzmin : input : minimum impurity fraction from Reinke model
+ fzactual : input : actual impurity fraction
+ """
+ cc = (
+ 1.0
+ - fortran.constraint_variables.freinke
+ * fortran.reinke_variables.fzactual
+ / fortran.reinke_variables.fzmin
+ )
+ return ConstraintResult(
+ cc,
+ fortran.reinke_variables.fzmin * (1.0 - cc),
+ fortran.reinke_variables.fzmin * cc,
+ )
+
+
+@ConstraintManager.register_constraint(79, "A/turn", "<=")
+def constraint_equation_79():
+ """Equation for maximum CS field
+ author: P B Lloyd, CCFE, Culham Science Centre
+
+ fb_cs_limit_max: F-value for CS mmax field (cons. 79, itvar 149)
+ b_cs_limit_max: Central solenoid max field limit [T]
+ b_cs_peak_pulse_start: maximum field in central solenoid at beginning of pulse (T)
+ b_cs_peak_flat_top_end: maximum field in central solenoid at end of flat-top (EoF) (T)
+ (Note: original code has "b_cs_peak_flat_top_end/b_cs_peak_pulse_start | peak CS field [T]".)
+ """
+ cc = (
+ max(
+ fortran.pfcoil_variables.b_cs_peak_flat_top_end,
+ fortran.pfcoil_variables.b_cs_peak_pulse_start,
+ )
+ / fortran.pfcoil_variables.b_cs_limit_max
+ - 1.0 * fortran.pfcoil_variables.fb_cs_limit_max
+ )
+ return ConstraintResult(
+ cc,
+ fortran.pfcoil_variables.b_cs_limit_max,
+ max(
+ fortran.pfcoil_variables.b_cs_peak_flat_top_end,
+ fortran.pfcoil_variables.b_cs_peak_pulse_start,
+ )
+ * cc,
+ )
+
+
+@ConstraintManager.register_constraint(80, "MW", ">=")
+def constraint_equation_80():
+ """Equation for p_plasma_separatrix_mw lower limit
+ author: J Morris, Culham Science Centre
+ args : output structure : residual error; constraint value; residual error in physical units;
+ output string; units string
+ Lower limit p_plasma_separatrix_mw
+ #=# physics
+ #=#=# fp_plasma_separatrix_min_mw, p_plasma_separatrix_mw
+ Logic change during pre-factoring: err, symbol, units will be assigned only if present.
+ fp_plasma_separatrix_min_mw : input : F-value for lower limit on p_plasma_separatrix_mw (cons. 80, itvar 153)
+ p_plasma_separatrix_min_mw : input : Minimum power crossing separatrix p_plasma_separatrix_mw [MW]
+ p_plasma_separatrix_mw : input : Power crossing separatrix [MW]
+ """
+ cc = (
+ 1.0
+ - fortran.physics_variables.fp_plasma_separatrix_min_mw
+ * fortran.physics_variables.p_plasma_separatrix_mw
+ / fortran.constraint_variables.p_plasma_separatrix_min_mw
+ )
+ return ConstraintResult(
+ cc,
+ fortran.constraint_variables.p_plasma_separatrix_min_mw,
+ fortran.constraint_variables.p_plasma_separatrix_min_mw * cc,
+ )
+
+
+@ConstraintManager.register_constraint(81, "m-3", ">=")
+def constraint_equation_81():
+ """Lower limit to ensure central density is larger that the pedestal one
+ author: S Kahn, Culham Science Centre
+ args : output structure : residual error; constraint value;
+ residual error in physical units; output string; units string
+ Lower limit ne0 > neped
+ !#=# physics
+ !#=#=# ne0, neped
+ Logic change during pre-factoring: err, symbol, units will be
+ assigned only if present.
+ fne0 : input : F-value for constraint on ne0 > neped
+ ne0 : input : Central electron density [m-3]
+ neped : input : Electron density at pedestal [m-3]
+ """
+ cc = (
+ 1.0
+ - fortran.physics_variables.fne0
+ * fortran.physics_variables.ne0
+ / fortran.physics_variables.neped
+ )
+ return ConstraintResult(
+ cc, fortran.physics_variables.fne0, fortran.physics_variables.fne0 * cc
+ )
+
+
+@ConstraintManager.register_constraint(82, "m", ">=")
+def constraint_equation_82():
+ """Equation for toroidal consistency of stellarator build
+ author: J Lion, IPP Greifswald
+
+ ftoroidalgap: f-value for constraint toroidalgap > dx_tf_inboard_out_toroidal
+ toroidalgap: minimal gap between two stellarator coils
+ dx_tf_inboard_out_toroidal: total toroidal width of a tf coil
+ """
+ return ConstraintResult(
+ 1.0
+ - fortran.tfcoil_variables.ftoroidalgap
+ * fortran.tfcoil_variables.toroidalgap
+ / fortran.tfcoil_variables.dx_tf_inboard_out_toroidal,
+ fortran.tfcoil_variables.toroidalgap,
+ fortran.tfcoil_variables.toroidalgap
+ - fortran.tfcoil_variables.dx_tf_inboard_out_toroidal
+ / fortran.tfcoil_variables.ftoroidalgap,
+ )
+
+
+@ConstraintManager.register_constraint(83, "m", ">=")
+def constraint_equation_83():
+ """Equation for radial consistency of stellarator build
+ author: J Lion, IPP Greifswald
+
+ f_avspace: f-value for constraint available_radial_space > required_radial_space
+ available_radial_space: avaible space in radial direction as given by each s.-configuration
+ required_radial_space: required space in radial direction
+ """
+ cc = (
+ 1.0
+ - fortran.build_variables.f_avspace
+ * fortran.build_variables.available_radial_space
+ / fortran.build_variables.required_radial_space
+ )
+ return ConstraintResult(
+ cc,
+ fortran.build_variables.available_radial_space * (1.0 - cc),
+ fortran.build_variables.required_radial_space * cc,
+ )
+
+
+@ConstraintManager.register_constraint(84, "", ">=")
+def constraint_equation_84():
+ """Equation for the lower limit of beta
+ author: J Lion, IPP Greifswald
+
+ fbeta_min: f-value for constraint beta-beta_fast_alpha > beta_min
+ beta_min: Lower limit for beta
+ beta: plasma beta
+ """
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fbeta_min
+ * fortran.physics_variables.beta
+ / fortran.physics_variables.beta_min
+ )
+ return ConstraintResult(
+ cc,
+ fortran.physics_variables.beta_min * (1.0 - cc),
+ fortran.physics_variables.beta * cc,
+ )
+
+
+@ConstraintManager.register_constraint(85, "years", "=")
+def constraint_equation_85():
+ """Equality constraint for the centerpost (CP) lifetime
+ Author : S Kahn
+
+ Depending on the chosen option i_cp_lifetime:
+ - 0 : The CP full power year lifelime is set by the user (cplife_input)
+ - 1 : The CP lifelime is equal to the divertor one
+ - 2 : The CP lifetime is equal to the breeding blankets one
+ - 3 : The CP lifetime is equal to the plant one
+
+ cplife: calculated CP full power year lifetime (years)
+ life_blkt_fpy: calculated first wall/blanket power year lifetime (years)
+ divlife: calculated divertor power year lifetime (years)
+ i_cp_lifetime: switch chosing which plant element the CP
+ the CP lifetime must equate
+ """
+ # The CP lifetime is equal to the the divertor one
+ if fortran.cost_variables.i_cp_lifetime == 0:
+ cc = 1.0 - fortran.cost_variables.cplife / fortran.cost_variables.cplife_input
+
+ elif fortran.cost_variables.i_cp_lifetime == 1:
+ cc = 1.0 - fortran.cost_variables.cplife / fortran.cost_variables.divlife
+
+ # The CP lifetime is equal to the tritium breeding blankets / FW one
+ elif fortran.cost_variables.i_cp_lifetime == 2:
+ cc = 1.0 - fortran.cost_variables.cplife / fortran.fwbs_variables.life_blkt_fpy
+
+ elif fortran.cost_variables.i_cp_lifetime == 3:
+ cc = 1.0 - fortran.cost_variables.cplife / fortran.cost_variables.tlife
+
+ return ConstraintResult(
+ cc,
+ fortran.cost_variables.divlife * (1.0 - cc),
+ fortran.cost_variables.divlife * cc,
+ )
+
+
+@ConstraintManager.register_constraint(86, "m", "<=")
+def constraint_equation_86():
+ """Upper limit on the turn edge length in the TF winding pack
+ Author : S Kahn
+
+ t_turn_tf: TF coil turn edge length including turn insulation [m]
+ f_t_turn_tf: f-value for TF turn edge length constraint
+ t_turn_tf_max: TF turn edge length including turn insulation upper limit [m]
+ """
+ cc = (
+ fortran.tfcoil_variables.t_turn_tf / fortran.tfcoil_variables.t_turn_tf_max
+ - 1.0 * fortran.tfcoil_variables.f_t_turn_tf
+ )
+ return ConstraintResult(
+ cc,
+ fortran.tfcoil_variables.t_turn_tf_max * (1.0 - cc),
+ fortran.tfcoil_variables.t_turn_tf_max * cc,
+ )
+
+
+@ConstraintManager.register_constraint(87, "MW", "<=")
+def constraint_equation_87():
+ """Equation for TF coil cryogenic power upper limit
+ author: S. Kahn, CCFE, Culham Science Centre
+
+ crypmw: cryogenic plant power (MW)
+ f_crypmw: f-value for maximum cryogenic plant power
+ crypmw_max: Maximum cryogenic plant power (MW)
+ """
+ cc = (
+ fortran.heat_transport_variables.crypmw
+ / fortran.heat_transport_variables.crypmw_max
+ - 1.0 * fortran.heat_transport_variables.f_crypmw
+ )
+ return ConstraintResult(
+ cc,
+ fortran.heat_transport_variables.crypmw_max * (1.0 - cc),
+ fortran.heat_transport_variables.crypmw * cc,
+ )
+
+
+@ConstraintManager.register_constraint(88, "", "<=")
+def constraint_equation_88():
+ """Equation for TF coil vertical strain upper limit (absolute value)
+ author: CPS Swanson, PPPL, USA
+
+ fstr_wp: f-value for TF coil strain
+ str_wp_max: Allowable maximum TF coil vertical strain
+ str_wp: Constrained TF coil vertical strain
+ """
+ return ConstraintResult(
+ abs(fortran.tfcoil_variables.str_wp) / fortran.tfcoil_variables.str_wp_max
+ - 1.0 * fortran.constraint_variables.fstr_wp,
+ fortran.tfcoil_variables.str_wp_max,
+ fortran.tfcoil_variables.str_wp_max
+ - abs(fortran.tfcoil_variables.str_wp) / fortran.constraint_variables.fstr_wp,
+ )
+
+
+@ConstraintManager.register_constraint(89, "A/m2", "<=")
+def constraint_equation_89():
+ """Upper limit to ensure that the Central Solenoid [OH] coil current / copper area < Maximum value
+ author: G Turkington, CCFE, Culham Science Centre
+
+ copperaoh_m2: CS coil current at EOF / copper area [A/m2]
+ copperaoh_m2_max: maximum coil current / copper area [A/m2]
+ f_copperaoh_m2: f-value for CS coil current / copper area
+ """
+ cc = (
+ fortran.rebco_variables.copperaoh_m2 / fortran.rebco_variables.copperaoh_m2_max
+ - 1.0 * fortran.rebco_variables.f_copperaoh_m2
+ )
+ return ConstraintResult(
+ cc,
+ fortran.rebco_variables.copperaoh_m2,
+ fortran.rebco_variables.copperaoh_m2 * cc,
+ )
+
+
+@ConstraintManager.register_constraint(90, "", ">=")
+def constraint_equation_90():
+ """Lower limit for CS coil stress load cycles
+ author: A. Pearce, G Turkington CCFE, Culham Science Centre
+
+ fncycle: f-value for constraint n_cycle > n_cycle_min
+ n_cycle: Allowable number of cycles for CS
+ n_cycle_min: Minimum required cycles for CS
+ """
+ if (
+ fortran.cost_variables.ibkt_life == 1
+ and fortran.cs_fatigue_variables.bkt_life_csf == 1
+ ):
+ fortran.cs_fatigue_variables.n_cycle_min = fortran.cost_variables.bktcycles
+
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fncycle
+ * fortran.cs_fatigue_variables.n_cycle
+ / fortran.cs_fatigue_variables.n_cycle_min
+ )
+ return ConstraintResult(
+ cc,
+ fortran.cs_fatigue_variables.n_cycle_min * (1.0 - cc),
+ fortran.cs_fatigue_variables.n_cycle * cc,
+ )
+
+
+@ConstraintManager.register_constraint(91, "MW", ">=")
+def constraint_equation_91():
+ """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).
+ author: J Lion, IPP Greifswald
+
+ fecrh_ignition: f-value for constraint powerht_local > powerscaling
+ max_gyrotron_frequency: Max. av. gyrotron frequency
+ te0_ecrh_achievable: Max. achievable electron temperature at ignition point
+ """
+ # Achievable ECRH te needs to be larger than needed te for igntion
+ if fortran.physics_variables.i_plasma_ignited == 0:
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fecrh_ignition
+ * (
+ fortran.stellarator_variables.powerht_constraint
+ + fortran.current_drive_variables.p_hcd_primary_extra_heat_mw
+ )
+ / fortran.stellarator_variables.powerscaling_constraint
+ )
+ else:
+ cc = (
+ 1.0
+ - fortran.constraint_variables.fecrh_ignition
+ * fortran.stellarator_variables.powerht_constraint
+ / fortran.stellarator_variables.powerscaling_constraint
+ )
+
+ return ConstraintResult(
+ cc,
+ fortran.stellarator_variables.powerscaling_constraint * (1.0 - cc),
+ fortran.stellarator_variables.powerht_constraint * cc,
+ )
+
+
+@ConstraintManager.register_constraint(92, "", "=")
+def constraint_equation_92():
+ """Equation for checking is D/T ratio is consistent, and sums to 1.
+ author: G Turkington, UKAEA
+
+ f_deuterium: fraction of deuterium ions
+ f_tritium: fraction of tritium ions
+ f_helium3: fraction of helium-3 ions
+ """
+ f_deuterium = 1.0 - (
+ fortran.physics_variables.f_tritium + fortran.physics_variables.f_helium3
+ )
+ cc = 1.0 - (
+ f_deuterium
+ + fortran.physics_variables.f_tritium
+ + fortran.physics_variables.f_helium3
+ )
+
+ return ConstraintResult(cc, 1.0, cc)
+
+
+def constraint_eqns(m: int, ieqn: int):
+ """Evaluates the constraints given the current state of PROCESS.
+
+ :param m: The number of constraints to evaluate
+ :param ieqn: Evaluates the 'ieqn'th constraint equation (index starts at 1)
+ or all equations if <= 0
+ """
+
+ if ieqn > 0:
+ i1 = ieqn - 1
+ i2 = ieqn
+ else:
+ i1 = 0
+ i2 = m
+
+ cc, con, err, symbol, units = [], [], [], [], []
+
+ for i in range(i1, i2):
+ constraint_id = fortran.numerics.icc[i].item()
+ try:
+ constraint = ConstraintManager.get_constraint(constraint_id)
+
+ if constraint is None:
+ tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units = getattr(
+ fortran.constraints, f"constraint_eqn_{constraint_id:03d}"
+ )()
+ else:
+ result = constraint.result()
+ tmp_cc, tmp_con, tmp_err = result.cc, result.con, result.err
+ tmp_symbol, tmp_units = constraint.symbol, constraint.units
+
+ except AttributeError as e:
+ error_msg = f"Constraint equation {i + 1} cannot be found"
+ raise ProcessError(error_msg) from e
+
+ if np.isnan(tmp_cc) or np.isinf(tmp_cc) or abs(tmp_cc) > 9.99e99:
+ error_msg = (
+ f"Constraint equation {constraint_id} returned an invalid residual"
+ )
+
+ raise ProcessValueError(error_msg, cc=tmp_cc)
+
+ # Reverse the sign so it works as an inequality constraint (cc(i) > 0)
+ # This will have no effect if it is used as an equality constraint because it will be squared.
+ cc.append(-tmp_cc)
+ con.append(tmp_con)
+ err.append(tmp_err)
+ symbol.append(tmp_symbol)
+ units.append(tmp_units)
+
+ return cc, con, err, symbol, units
+
+
+def init_constraint_variables():
+ fortran.constraint_variables.auxmin = 0.1
+ fortran.constraint_variables.beta_poloidal_max = 0.19
+ fortran.constraint_variables.bigqmin = 10.0
+ fortran.constraint_variables.bmxlim = 12.0
+ fortran.constraint_variables.fauxmn = 1.0
+ fortran.constraint_variables.fbeta_poloidal_eps = 1.0
+ fortran.constraint_variables.fbeta_poloidal = 1.0
+ fortran.constraint_variables.fbeta_max = 1.0
+ fortran.constraint_variables.fbeta_min = 1.0
+ fortran.constraint_variables.fcpttf = 1.0
+ fortran.constraint_variables.fr_conducting_wall = 1.0
+ fortran.constraint_variables.fdene = 1.0
+ fortran.constraint_variables.fdtmp = 1.0
+ fortran.constraint_variables.fflutf = 1.0
+ fortran.constraint_variables.ffuspow = 1.0
+ fortran.constraint_variables.fgamcd = 1.0
+ fortran.constraint_variables.fpflux_div_heat_load_mw = 1.0
+ fortran.constraint_variables.fiooic = 0.5
+ fortran.constraint_variables.fipir = 1.0
+ fortran.constraint_variables.q95_fixed = 3.0
+ fortran.constraint_variables.fjohc = 1.0
+ fortran.constraint_variables.fjohc0 = 1.0
+ fortran.constraint_variables.fjprot = 1.0
+ fortran.constraint_variables.fl_h_threshold = 1.0
+ fortran.constraint_variables.fmva = 1.0
+ fortran.constraint_variables.fnbshinef = 1.0
+ fortran.constraint_variables.fncycle = 1.0
+ fortran.constraint_variables.fnesep = 1.0
+ fortran.constraint_variables.foh_stress = 1.0
+ fortran.constraint_variables.fpeakb = 1.0
+ fortran.constraint_variables.fpinj = 1.0
+ fortran.constraint_variables.fpnetel = 1.0
+ fortran.constraint_variables.fportsz = 1.0
+ fortran.constraint_variables.fpsepbqar = 1.0
+ fortran.constraint_variables.fpsepr = 1.0
+ fortran.constraint_variables.fptemp = 1.0
+ fortran.constraint_variables.fptfnuc = 1.0
+ fortran.constraint_variables.fq = 1.0
+ fortran.constraint_variables.fqval = 1.0
+ fortran.constraint_variables.fradpwr = 0.99
+ fortran.constraint_variables.fradwall = 1.0
+ fortran.constraint_variables.freinke = 1.0
+ fortran.constraint_variables.frminor = 1.0
+ fortran.constraint_variables.fstrcase = 1.0
+ fortran.constraint_variables.fstrcond = 1.0
+ fortran.constraint_variables.fstr_wp = 1.0
+ fortran.constraint_variables.fmaxvvstress = 1.0
+ fortran.constraint_variables.ftbr = 1.0
+ fortran.constraint_variables.ft_burn = 1.0
+ fortran.constraint_variables.ftcycl = 1.0
+ fortran.constraint_variables.ftmargoh = 1.0
+ fortran.constraint_variables.ftmargtf = 1.0
+ fortran.constraint_variables.ft_current_ramp_up = 1.0
+ fortran.constraint_variables.ftpeak = 1.0
+ fortran.constraint_variables.fvdump = 1.0
+ fortran.constraint_variables.fvs = 1.0
+ fortran.constraint_variables.fvvhe = 1.0
+ fortran.constraint_variables.fwalld = 1.0
+ fortran.constraint_variables.fzeffmax = 1.0
+ fortran.constraint_variables.gammax = 2.0
+ fortran.constraint_variables.i_q95_fixed = 0
+ fortran.constraint_variables.pflux_fw_rad_max = 1.0
+ fortran.constraint_variables.mvalim = 40.0
+ fortran.constraint_variables.f_p_beam_shine_through_max = 1e-3
+ fortran.constraint_variables.nflutfmax = 1.0e23
+ fortran.constraint_variables.p_plasma_separatrix_min_mw = 150.0
+ fortran.constraint_variables.f_fw_rad_max = 3.33
+ fortran.constraint_variables.pflux_fw_rad_max_mw = 0.0
+ fortran.constraint_variables.pnetelin = 1.0e3
+ fortran.constraint_variables.powfmax = 1.5e3
+ fortran.constraint_variables.psepbqarmax = 9.5
+ fortran.constraint_variables.pseprmax = 25.0
+ fortran.constraint_variables.ptfnucmax = 1e-3
+ fortran.constraint_variables.tbrmin = 1.1
+ fortran.constraint_variables.t_burn_min = 1.0
+ fortran.constraint_variables.tcycmn = 0.0
+ fortran.constraint_variables.t_current_ramp_up_min = 1.0
+ fortran.constraint_variables.vvhealw = 1.0
+ fortran.constraint_variables.walalw = 1.0
+ fortran.constraint_variables.f_alpha_energy_confinement_min = 5.0
+ fortran.constraint_variables.falpha_energy_confinement = 1.0
+ fortran.constraint_variables.fniterpump = 1.0
+ fortran.constraint_variables.zeffmax = 3.6
+ fortran.constraint_variables.fpoloidalpower = 1.0
+ fortran.constraint_variables.fpsep = 1.0
+ fortran.constraint_variables.fcqt = 1.0
+ fortran.constraint_variables.fecrh_ignition = 1.0
diff --git a/process/current_drive.py b/process/current_drive.py
index e6e80f6f04..015612653d 100644
--- a/process/current_drive.py
+++ b/process/current_drive.py
@@ -10,10 +10,8 @@
heat_transport_variables,
physics_variables,
)
-from process.fortran import (
- error_handling as eh,
-)
from process.plasma_profiles import PlasmaProfile
+from process.warning_handler import WarningManager
class NeutralBeam:
@@ -793,7 +791,6 @@ def legend(self, zlocal, arg):
Abramowitz and Stegun, equation 8.12.1
"""
if abs(arg) > (1.0e0 + 1.0e-10):
- eh.fdiags[0] = arg
raise ProcessValueError("Invalid argument", arg=arg)
arg2 = min(arg, (1.0e0 - 1.0e-10))
@@ -1073,7 +1070,9 @@ def lhrad(self):
rat0 = rat1
else:
- eh.report_error(16)
+ WarningManager.create_warning(
+ "LH penetration radius not found after lapno iterations, using 0.8*rminor"
+ )
rat0 = 0.8e0
return rat0
diff --git a/process/final.py b/process/final.py
index 34304a85e1..4db43f35fa 100644
--- a/process/final.py
+++ b/process/final.py
@@ -2,13 +2,13 @@
from tabulate import tabulate
+import process.constraints as constraints
from process import output as op
from process import (
process_output as po,
)
from process.fortran import (
constants,
- constraints,
numerics,
)
from process.objectives import objective_function
@@ -67,20 +67,15 @@ def output_once_through():
f2py_compatible_to_string(i)
for i in numerics.lablcc[numerics.icc[: numerics.neqns + numerics.nineqns] - 1]
]
- units = [f2py_compatible_to_string(i) for i in units]
- physical_constraint = [
- f"{c} {u}" for c, u in zip(value.tolist(), units, strict=False)
- ]
- physical_residual = [
- f"{c} {u}" for c, u in zip(residual.tolist(), units, strict=False)
- ]
+ physical_constraint = [f"{c} {u}" for c, u in zip(value, units, strict=False)]
+ physical_residual = [f"{c} {u}" for c, u in zip(residual, units, strict=False)]
table_data = {
"Constraint Name": labels,
- "Constraint Type": symbols.tolist(),
+ "Constraint Type": symbols,
"Physical constraint": physical_constraint,
"Constraint residual": physical_residual,
- "Normalised residual": residual_error.tolist(),
+ "Normalised residual": residual_error,
}
po.write(constants.nout, tabulate(table_data, headers="keys"))
diff --git a/process/fw.py b/process/fw.py
index 81112008db..c5b5562375 100644
--- a/process/fw.py
+++ b/process/fw.py
@@ -8,13 +8,10 @@
blanket_library,
build_variables,
constants,
- error_handling,
fwbs_variables,
)
-from process.fortran import (
- error_handling as eh,
-)
from process.utilities.f2py_string_patch import f2py_compatible_to_string
+from process.warning_handler import WarningManager
class Fw:
@@ -141,10 +138,12 @@ def fw_temp(
# Print debug info if temperature too low/high or NaN/Inf
if np.isnan(temp_k):
- eh.report_error(223)
+ WarningManager.create_warning("NaN first wall temperature")
elif (temp_k <= 100) or (temp_k > 1500):
- eh.fdiags[0] = temp_k
- eh.report_error(224)
+ WarningManager.create_warning(
+ "First wall temperature (temp_k) out of range : [100-1500] K",
+ temp_k=temp_k,
+ )
# Thermal conductivity of first wall material (W/m.K)
tkfw = self.fw_thermal_conductivity(temp_k)
@@ -397,17 +396,19 @@ def heat_transfer(self, masflx, rhof, radius, cf, viscf, kf):
# Check that Reynolds number is in valid range for the Gnielinski correlation
if (reynolds <= 3000.0) or (reynolds > 5.0e6):
- error_handling.fdiags[0] = reynolds
- error_handling.report_error(225)
+ WarningManager.create_warning(
+ "Reynolds number out of range : [3e3-5000e3]", reynolds=reynolds
+ )
# Check that Prandtl number is in valid range for the Gnielinski correlation
if (pr < 0.5) or (pr > 2000.0):
- error_handling.fdiags[0] = pr
- error_handling.report_error(226)
+ WarningManager.create_warning(
+ "Prandtl number out of range : [0.5-2000]", pr=pr
+ )
# Check that the Darcy friction factor is in valid range for the Gnielinski correlation
if f <= 0.0:
- error_handling.report_error(227)
+ WarningManager.create_warning("Negative Darcy friction factor (f)")
return heat_transfer
diff --git a/process/hcpb.py b/process/hcpb.py
index ad44d24bed..a8261b54d2 100644
--- a/process/hcpb.py
+++ b/process/hcpb.py
@@ -18,9 +18,7 @@
primary_pumping_variables,
tfcoil_variables,
)
-from process.fortran import (
- error_handling as eh,
-)
+from process.warning_handler import WarningManager
class CCFE_HCPB:
@@ -578,11 +576,13 @@ def nuclear_heating_blanket(self):
)
if fwbs_variables.p_blkt_nuclear_heat_total_mw < 1:
- eh.fdiags[0] = fwbs_variables.p_blkt_nuclear_heat_total_mw
- eh.fdiags[1] = ccfe_hcpb_module.exp_blanket
- eh.fdiags[2] = physics_variables.fusion_power
- eh.fdiags[3] = mass
- eh.report_error(274)
+ WarningManager.create_warning(
+ "Blanket heating is <1 MW or NaN. Is something wrong?",
+ p_blkt_nuclear_heat_total_mw=fwbs_variables.p_blkt_nuclear_heat_total_mw,
+ exp_blanket=ccfe_hcpb_module.exp_blanket,
+ fusion_power=physics_variables.fusion_power,
+ mass=mass,
+ )
def nuclear_heating_shield(self):
"""Nuclear heating in the shield for CCFE HCPB model
diff --git a/process/init.py b/process/init.py
index bef72bbde6..3978557148 100644
--- a/process/init.py
+++ b/process/init.py
@@ -12,6 +12,7 @@
from process.blanket_library import init_blanket_library, init_primary_pumping_variables
from process.build import init_build_variables
from process.buildings import init_buildings_variables
+from process.constraints import init_constraint_variables
from process.costs import init_cost_variables
from process.cs_fatigue import init_cs_fatigue_variables
from process.current_drive import init_current_drive_variables
@@ -43,6 +44,7 @@
from process.tf_coil import init_tfcoil_variables
from process.utilities.f2py_string_patch import f2py_compatible_to_string
from process.vacuum import init_vacuum_variables
+from process.warning_handler import WarningManager
from process.water_use import init_watuse_variables
@@ -54,9 +56,6 @@ def init_process():
the default values for the global variables, reads in data from
the input file, and checks the run parameters for consistency.
"""
- # Initialise error handling
- fortran.error_handling.initialise_error_list()
-
# Initialise the program variables
iteration_variables.initialise_iteration_variables()
@@ -240,11 +239,11 @@ def init_all_module_vars():
run. This matters ever since Process is used as a shared library, rather
than a 'run-once' executable.
"""
+ WarningManager.reinitialise()
fortran.numerics.init_numerics()
init_buildings_variables()
init_cost_variables()
init_divertor_variables()
- fortran.error_handling.init_error_handling()
init_fwbs_variables()
fortran.global_variables.init_global_variables()
init_ccfe_hcpb_module()
@@ -268,7 +267,7 @@ def init_all_module_vars():
init_vacuum_variables()
init_pf_power_variables()
init_build_variables()
- fortran.constraint_variables.init_constraint_variables()
+ init_constraint_variables()
init_pulse_variables()
init_rebco_variables()
init_reinke_variables()
@@ -277,8 +276,6 @@ def init_all_module_vars():
init_blanket_library()
init_dcll_module()
- fortran.init_module.init_fortran_modules()
-
def check_process(inputs): # noqa: ARG001
"""Routine to reset specific variables if certain options are
diff --git a/process/io/process_config.py b/process/io/process_config.py
index 85ad9c87f5..98682b077c 100644
--- a/process/io/process_config.py
+++ b/process/io/process_config.py
@@ -125,10 +125,7 @@ def error_status2readme(self, directory="."):
if os.path.isfile("MFILE.DAT"):
m_file = MFile(filename=directory + "/MFILE.DAT")
- error_status = (
- f"Error status: {m_file.data['error_status'].get_scan(-1)} "
- f"Error ID: {m_file.data['error_id'].get_scan(-1)}\n"
- )
+ error_status = f"Error status: {m_file.data['error_status'].get_scan(-1)} "
if self.comment != "":
with open(directory + "/README.txt", "a") as readme:
@@ -1151,7 +1148,7 @@ def write_error_summary(self, sample_index):
header += f" n_{label.replace('_(range_normalised)', ''):8s}"
# error status, id and ifail
- header += " error_status error_id ifail\n"
+ header += " error_status ifail\n"
with open(self.wdir + "/UQ_error_summary.txt", "w") as err_summary:
err_summary.write(header)
@@ -1165,10 +1162,7 @@ def write_error_summary(self, sample_index):
output += f" {m_file.data[f'nitvar{i:03}'].get_scan(-1):10f}"
# error status and id
- output += (
- f" {int(m_file.data['error_status'].get_scan(-1)):13d} "
- f"{int(m_file.data['error_id'].get_scan(-1)):8d}"
- )
+ output += f" {int(m_file.data['error_status'].get_scan(-1)):13d} "
# ifail
if m_file.data["error_status"].get_scan(-1) < 3:
output += f" {int(m_file.data['ifail'].get_scan(-1)):5d}\n"
diff --git a/process/io/process_funcs.py b/process/io/process_funcs.py
index 37ddb54dd7..6b0ebd74ac 100644
--- a/process/io/process_funcs.py
+++ b/process/io/process_funcs.py
@@ -12,7 +12,6 @@
import logging
from os.path import join as pjoin
-from pathlib import Path
from sys import stderr
from time import sleep
@@ -222,7 +221,6 @@ def check_logfile(logfile="process.log"):
Checks the log file of the PROCESS output.
Stops, if an error occured that needs to be
fixed before rerunning.
- XXX should be deprecated!! and replaced by check_input_error!
"""
with open(logfile) as outlogfile:
@@ -236,34 +234,6 @@ def check_logfile(logfile="process.log"):
exit()
-def check_input_error(wdir="."):
- """
- Checks, if an input error has occurred.
- Stops as a consequence.
- Will also fail if the MFILE.DAT isn't found.
- """
- try:
- mfile_path = Path(wdir) / "MFILE.DAT"
-
- if mfile_path.exists():
- mfile_path_str = str(mfile_path)
- mfile = MFile(filename=mfile_path_str)
- else:
- raise FileNotFoundError("MFile doesn't exist")
-
- error_id = mfile.data["error_id"].get_scan(-1)
-
- if error_id == 130:
- print(
- "Error in input file. Please check OUT.DAT for more information.",
- file=stderr,
- )
- exit()
- except Exception:
- logger.exception("Check input error exception")
- raise
-
-
########################################
diff --git a/process/main.py b/process/main.py
index 8b3c8c5312..d62f3441e6 100644
--- a/process/main.py
+++ b/process/main.py
@@ -79,7 +79,6 @@
# For VaryRun
from process.io.process_config import RunProcessConfig
from process.io.process_funcs import (
- check_input_error,
get_neqns_itervars,
get_variable_range,
no_unfeasible_mfile,
@@ -102,6 +101,7 @@
from process.tf_coil import TFCoil
from process.utilities.f2py_string_patch import string_to_f2py_compatible
from process.vacuum import Vacuum
+from process.warning_handler import WarningManager
from process.water_use import WaterUse
os.environ["PYTHON_PROCESS_ROOT"] = os.path.join(os.path.dirname(__file__))
@@ -312,12 +312,6 @@ def run(self):
neqns, itervars = get_neqns_itervars()
lbs, ubs = get_variable_range(itervars, config.factor)
- # If config file contains WDIR, use that. Otherwise, use the directory
- # containing the config file (used when running regression tests in
- # temp dirs)
- # TODO Not sure this is required any more
- wdir = config.wdir if config.wdir else Path(self.config_file).parent
-
# Check IN.DAT exists
if not input_path.exists():
raise FileNotFoundError
@@ -336,8 +330,6 @@ def run(self):
# Run process on an IN.DAT file
config.run_process(input_path, self.solver)
- check_input_error(wdir=wdir)
-
if not process_stopped():
no_unfeasible = no_unfeasible_mfile()
if no_unfeasible <= config.no_allowed_unfeasible:
@@ -500,7 +492,7 @@ def run_scan(self, solver):
def show_errors(self):
"""Report all informational/error messages encountered."""
- fortran.error_handling.show_errors()
+ WarningManager.show_errors(fortran.constants.nout)
def finish(self):
"""Run the finish subroutine to close files open in the Fortran.
diff --git a/process/output.py b/process/output.py
index 98c95c4060..d936926f38 100644
--- a/process/output.py
+++ b/process/output.py
@@ -11,11 +11,6 @@ def write(models, _outfile):
:param outfile: Fortran output unit identifier
:type outfile: int
"""
- # Turn on error reporting
- # (warnings etc. encountered in previous iterations may have cleared themselves
- # during the solution process)
- ft.error_handling.errors_on = True
-
# Call stellarator output routine instead if relevant
if ft.stellarator_variables.istell != 0:
models.stellarator.run(output=True)
diff --git a/process/pfcoil.py b/process/pfcoil.py
index 7b13760b45..e6cc925d39 100644
--- a/process/pfcoil.py
+++ b/process/pfcoil.py
@@ -15,7 +15,6 @@
from process.fortran import constants, numerics
from process.fortran import constraint_variables as ctv
from process.fortran import cs_fatigue_variables as csfv
-from process.fortran import error_handling as eh
from process.fortran import fwbs_variables as fwbsv
from process.fortran import pfcoil_module as pf
from process.fortran import pfcoil_variables as pfv
@@ -24,6 +23,7 @@
from process.fortran import tfcoil_variables as tfv
from process.fortran import times_variables as tv
from process.utilities.f2py_string_patch import f2py_compatible_to_string
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
@@ -514,7 +514,9 @@ def pfcoil(self):
else:
dics = 0.0e0
pfv.f_j_cs_start_end_flat_top = 1.0e0
- eh.report_error(71)
+ WarningManager.create_warning(
+ "OH coil not present; check volt-second calculations..."
+ )
# Split groups of coils into one set containing ncl coils
ncl = 0
@@ -1032,7 +1034,9 @@ def tf_pf_collision_detector(self):
pf_tf_collision += 1
if pf_tf_collision >= 1:
- eh.report_error(277)
+ WarningManager.create_warning(
+ "One or more collision between TF and PF coils. Check PF placement."
+ )
def solv(self, n_pf_groups_max, n_pf_coil_groups, nrws, gmat, bvec):
"""Solve a matrix using singular value decomposition.
@@ -1844,10 +1848,12 @@ def induct(self, output):
)
if noh > nohmax:
- eh.idiags[0] = noh
- eh.idiags[1] = nohmax
- eh.fdiags[0] = bv.dr_cs
- eh.report_error(73)
+ WarningManager.create_warning(
+ "Max no. of segments noh for OH coil > nohmax; increase dr_cs lower bound",
+ noh=noh,
+ nohmax=nohmax,
+ dr_cs=bv.dr_cs,
+ )
noh = min(noh, nohmax)
@@ -1892,9 +1898,6 @@ def induct(self, output):
if bv.dr_cs >= delzoh:
deltar = math.sqrt((bv.dr_cs**2 - delzoh**2) / 12.0e0)
else:
- # eh.fdiags[0] = bv.dr_cs
- # eh.fdiags[1] = delzoh
- # eh.report_error(74)
# Set deltar to something small and +ve instead; allows solver
# to continue and hopefully be constrained away from this point
deltar = 1.0e-6
@@ -2419,13 +2422,19 @@ def outpf(self):
if pfv.temp_cs_margin < 1.01e0 * tfv.tmargmin_cs:
pf.cslimit = True
if not pf.cslimit:
- eh.report_error(135)
+ WarningManager.create_warning(
+ "CS not using max current density: further optimisation may be possible"
+ )
# Check whether CS coil currents are feasible from engineering POV
if ctv.fjohc > 0.7:
- eh.report_error(286)
+ WarningManager.create_warning(
+ "fjohc shouldn't be above 0.7 for engineering reliability"
+ )
if ctv.fjohc0 > 0.7:
- eh.report_error(287)
+ WarningManager.create_warning(
+ "fjohc0 shouldn't be above 0.7 for engineering reliability"
+ )
# REBCO fractures in strains above ~+/- 0.7%
if (
@@ -2433,14 +2442,18 @@ def outpf(self):
or pfv.i_cs_superconductor == 8
or pfv.i_cs_superconductor == 9
) and abs(tfv.str_cs_con_res) > 0.7e-2:
- eh.report_error(262)
+ WarningManager.create_warning(
+ "Non physical strain used in CS. Use superconductor strain < +/- 0.7%"
+ )
if (
pfv.i_pf_superconductor == 6
or pfv.i_pf_superconductor == 8
or pfv.i_pf_superconductor == 9
) and abs(tfv.str_pf_con_res) > 0.7e-2:
- eh.report_error(263)
+ WarningManager.create_warning(
+ "Non physical strain used in PF. Use superconductor strain < +/- 0.7%"
+ )
else:
op.ocmmnt(self.outfile, "Resistive central solenoid")
diff --git a/process/physics.py b/process/physics.py
index 22fdb94d0e..c34bd17b43 100644
--- a/process/physics.py
+++ b/process/physics.py
@@ -21,7 +21,6 @@
constraint_variables,
current_drive_variables,
divertor_variables,
- error_handling,
fwbs_variables,
impurity_radiation_module,
numerics,
@@ -33,6 +32,7 @@
times_variables,
)
from process.utilities.f2py_string_patch import f2py_compatible_to_string
+from process.warning_handler import WarningManager
@nb.jit(nopython=True, cache=True)
@@ -2653,7 +2653,9 @@ def physics(self):
)
if reinke_variables.fzmin >= 1.0e0:
- error_handling.report_error(217)
+ WarningManager.create_warning(
+ "fzmin is greater than or equal to 1.0, this is at least notable"
+ )
po.write(
self.outfile,
@@ -2664,8 +2666,8 @@ def physics(self):
),
)
- @staticmethod
def calculate_density_limit(
+ self,
bt: float,
i_density_limit: int,
p_plasma_separatrix_mw: float,
@@ -2758,9 +2760,10 @@ def calculate_density_limit(
denom = (zeff - 1.0) * (1.0 - 4.0 / (3.0 * qcyl))
if denom <= 0.0:
if i_density_limit == 4:
- error_handling.fdiags[0] = denom
- error_handling.fdiags[1] = qcyl
- error_handling.report_error(80)
+ WarningManager.create_warning(
+ "qcyl < 4/3; dlimit(4) set to zero; model 5 will be enforced instead.",
+ qcyl=qcyl,
+ )
i_density_limit = 5
dlimit[3] = 0.0
@@ -3177,8 +3180,8 @@ def phyaux(
f_alpha_energy_confinement,
)
- @staticmethod
def plasma_ohmic_heating(
+ self,
f_c_plasma_inductive: float,
kappa95: float,
plasma_current: float,
@@ -3239,9 +3242,11 @@ def plasma_ohmic_heating(
# Check to see if plasma resistance is negative
# (possible if aspect ratio is too high)
if res_plasma <= 0.0:
- error_handling.fdiags[0] = res_plasma
- error_handling.fdiags[1] = physics_variables.aspect
- error_handling.report_error(83)
+ WarningManager.create_warning(
+ "Negative plasma resistance res_plasma",
+ res_plasma=res_plasma,
+ aspect=physics_variables.aspect,
+ )
# Ohmic heating power per unit volume
# Corrected from: pden_plasma_ohmic_mw = (f_c_plasma_inductive*plasma_current)**2 * ...
@@ -4344,7 +4349,9 @@ def outplas(self):
if physics_variables.ipedestal >= 1:
if physics_variables.ne0 < physics_variables.neped:
- error_handling.report_error(213)
+ WarningManager.create_warning(
+ "Central density is less than pedestal density"
+ )
po.ocmmnt(self.outfile, "Pedestal profiles are used.")
po.ovarrf(
@@ -5052,8 +5059,10 @@ def outplas(self):
)
if physics_variables.p_plasma_separatrix_mw <= 0.001e0:
- error_handling.fdiags[0] = physics_variables.p_plasma_separatrix_mw
- error_handling.report_error(87)
+ WarningManager.create_warning(
+ "Possible problem with high radiation power, forcing p_plasma_separatrix_mw to odd values",
+ p_plasma_separatrix_mw=physics_variables.p_plasma_separatrix_mw,
+ )
po.oblnkl(self.outfile)
po.ocmmnt(
self.outfile, " BEWARE: possible problem with high radiation power"
@@ -5275,13 +5284,15 @@ def outplas(self):
self.outfile,
"(physics_variables.bt outside Snipes 2000 fitted range)",
)
- error_handling.report_error(201)
+ WarningManager.create_warning("bt outside Snipes 2000 fitted range")
if (physics_variables.rminor < 0.15e0) or (
physics_variables.rminor > 1.15e0
):
po.ocmmnt(self.outfile, "(rminor outside Snipes 2000 fitted range)")
- error_handling.report_error(202)
+ WarningManager.create_warning(
+ "rminor outside Snipes 2000 fitted range"
+ )
if (physics_variables.rmajor < 0.55e0) or (
physics_variables.rmajor > 3.37e0
@@ -5290,7 +5301,9 @@ def outplas(self):
self.outfile,
"(physics_variables.rmajor outside Snipes 2000 fitted range)",
)
- error_handling.report_error(203)
+ WarningManager.create_warning(
+ "rmajor outside Snipes 2000 fitted range"
+ )
if (physics_variables.dnla < 0.09e20) or (
physics_variables.dnla > 3.16e20
@@ -5299,7 +5312,9 @@ def outplas(self):
self.outfile,
"(physics_variables.dnla outside Snipes 2000 fitted range)",
)
- error_handling.report_error(204)
+ WarningManager.create_warning(
+ "dnla outside Snipes 2000 fitted range"
+ )
if (physics_variables.kappa < 1.0e0) or (
physics_variables.kappa > 2.04e0
@@ -5308,13 +5323,17 @@ def outplas(self):
self.outfile,
"(physics_variables.kappa outside Snipes 2000 fitted range)",
)
- error_handling.report_error(205)
+ WarningManager.create_warning(
+ "kappa outside Snipes 2000 fitted range"
+ )
if (physics_variables.triang < 0.07e0) or (
physics_variables.triang > 0.74e0
):
po.ocmmnt(self.outfile, "(triang outside Snipes 2000 fitted range)")
- error_handling.report_error(206)
+ WarningManager.create_warning(
+ "triang outside Snipes 2000 fitted range"
+ )
po.oblnkl(self.outfile)
@@ -5324,7 +5343,9 @@ def outplas(self):
"(L-H threshold for closed divertor only. Limited data used in Snipes fit)",
)
po.oblnkl(self.outfile)
- error_handling.report_error(207)
+ WarningManager.create_warning(
+ "Closed divertor only. Limited data used in Snipes fit"
+ )
if (numerics.ioptimz > 0) and (numerics.active_constraints[14]):
po.ovarre(
@@ -5790,11 +5811,13 @@ def outplas(self):
)
# Error to catch if bootstap fraction limit has been enforced
if physics_module.err242 == 1:
- error_handling.report_error(242)
+ WarningManager.create_warning("Bootstrap fraction upper limit enforced")
# Error to catch if self-driven current fraction limit has been enforced
if physics_module.err243 == 1:
- error_handling.report_error(243)
+ WarningManager.create_warning(
+ "Predicted plasma driven current is more than upper limit on non-inductive fraction"
+ )
if current_drive_variables.f_c_plasma_bootstrap_max < 0.0e0:
po.ocmmnt(
@@ -5860,7 +5883,10 @@ def outplas(self):
)
# Error to show if diamagnetic current is above 1% but not used
if current_drive_variables.f_c_plasma_diamagnetic_scene > 0.01e0:
- error_handling.report_error(244)
+ WarningManager.create_warning(
+ "Diamagnetic fraction is more than 1%, but not calculated. "
+ "Consider using i_diamagnetic_current=2 and i_pfirsch_schluter_current=1."
+ )
elif physics_variables.i_diamagnetic_current == 1:
po.ocmmnt(
diff --git a/process/power.py b/process/power.py
index 05a020a3f1..75d7d18591 100644
--- a/process/power.py
+++ b/process/power.py
@@ -12,7 +12,6 @@
constraint_variables,
cost_variables,
current_drive_variables,
- error_handling,
fwbs_variables,
heat_transport_variables,
numerics,
@@ -25,6 +24,7 @@
times_variables,
)
from process.variables import AnnotatedVariable
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
@@ -2474,9 +2474,10 @@ def plant_thermal_efficiency(self, etath):
if (heat_transport_variables.tturb < 657.0e0) or (
heat_transport_variables.tturb > 915.0e0
):
- error_handling.idiags[0] = 2
- error_handling.fdiags[0] = heat_transport_variables.tturb
- error_handling.report_error(166)
+ WarningManager.create_warning(
+ "Turbine temperature tturb out of range of validity",
+ tturb=heat_transport_variables.tturb,
+ )
etath = (
0.1802e0 * np.log(heat_transport_variables.tturb)
@@ -2493,9 +2494,10 @@ def plant_thermal_efficiency(self, etath):
if (heat_transport_variables.tturb < 657.0e0) or (
heat_transport_variables.tturb > 915.0e0
):
- error_handling.idiags[0] = 2
- error_handling.fdiags[0] = heat_transport_variables.tturb
- error_handling.report_error(166)
+ WarningManager.create_warning(
+ "Turbine temperature tturb out of range of validity",
+ tturb=heat_transport_variables.tturb,
+ )
etath = (
0.1802e0 * np.log(heat_transport_variables.tturb)
@@ -2522,9 +2524,10 @@ def plant_thermal_efficiency(self, etath):
if (heat_transport_variables.tturb < 408.0e0) or (
heat_transport_variables.tturb > 1023.0e0
):
- error_handling.idiags[0] = 3
- error_handling.fdiags[0] = heat_transport_variables.tturb
- error_handling.report_error(166)
+ WarningManager.create_warning(
+ "Turbine temperature tturb out of range of validity",
+ tturb=heat_transport_variables.tturb,
+ )
etath = 0.4347e0 * np.log(heat_transport_variables.tturb) - 2.5043e0
@@ -2551,9 +2554,10 @@ def plant_thermal_efficiency_2(self, etath_liq):
if (heat_transport_variables.tturb < 408.0e0) or (
heat_transport_variables.tturb > 1023.0e0
):
- error_handling.idiags[0] = 3
- error_handling.fdiags[0] = heat_transport_variables.tturb
- error_handling.report_error(166)
+ WarningManager.create_warning(
+ "Turbine temperature tturb out of range of validity",
+ tturb=heat_transport_variables.tturb,
+ )
return 0.4347e0 * np.log(heat_transport_variables.tturb) - 2.5043e0
diff --git a/process/profiles.py b/process/profiles.py
index 758ee1930f..a2db4c8347 100644
--- a/process/profiles.py
+++ b/process/profiles.py
@@ -4,7 +4,8 @@
import numpy as np
import scipy as sp
-from process.fortran import error_handling, physics_variables
+from process.fortran import physics_variables
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
# Logging handler for console output
@@ -211,7 +212,9 @@ def ncore(
if ncore < 0.0:
# Allows solver to continue and
# warns the user to raise the lower bound on dene if the run did not converge
- error_handling.report_error(282)
+ WarningManager.create_warning(
+ "ncore is going negative when solving. Please raise the value of dene and or its lower limit"
+ )
ncore = 1.0e-6
return ncore
diff --git a/process/pulse.py b/process/pulse.py
index 669dac2069..e9ddb9fa2e 100644
--- a/process/pulse.py
+++ b/process/pulse.py
@@ -2,7 +2,6 @@
from process.fortran import (
constants,
constraint_variables,
- error_handling,
numerics,
pf_power_variables,
pfcoil_variables,
@@ -10,6 +9,7 @@
pulse_variables,
times_variables,
)
+from process.warning_handler import WarningManager
class Pulse:
@@ -185,12 +185,13 @@ def burn(self, output: bool):
tb = vsmax / v_plasma_loop_burn - times_variables.t_fusion_ramp
if tb < 0.0e0:
- error_handling.fdiags[0] = tb
- error_handling.fdiags[1] = vsmax
- error_handling.fdiags[2] = v_plasma_loop_burn
- error_handling.fdiags[3] = times_variables.t_fusion_ramp
- error_handling.report_error(93)
-
+ WarningManager.create_warning(
+ "Negative burn time available; reduce t_fusion_ramp or raise PF coil V-s capability",
+ tb=tb,
+ vsmax=vsmax,
+ v_plasma_loop_burn=v_plasma_loop_burn,
+ t_fusion_ramp=times_variables.t_fusion_ramp.item(),
+ )
times_variables.t_burn = max(0.0e0, tb)
# Output section
diff --git a/process/resistive_tf_coil.py b/process/resistive_tf_coil.py
index 5426dd60d9..c8f5598433 100644
--- a/process/resistive_tf_coil.py
+++ b/process/resistive_tf_coil.py
@@ -6,13 +6,13 @@
from process.fortran import (
build_variables,
constants,
- error_handling,
pfcoil_variables,
physics_variables,
sctfcoil_module,
tfcoil_variables,
)
from process.tf_coil import TFCoil
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
@@ -233,7 +233,9 @@ def run(self, output: bool):
)
except ValueError as e:
if e.args[1] == 245 and e.args[2] == 0:
- error_handling.report_error(245)
+ WarningManager.create_warning(
+ "Invalid stress model (r_tf_inboard = 0), stress constraint switched off"
+ )
tfcoil_variables.sig_tf_case = 0.0e0
tfcoil_variables.sig_tf_wp = 0.0e0
@@ -367,13 +369,16 @@ def res_tf_internal_geom(self):
# Reporting negative WP areas issues
if sctfcoil_module.awpc < 0.0e0:
- error_handling.fdiags[0] = sctfcoil_module.awpc
- error_handling.fdiags[0] = tfcoil_variables.dr_tf_wp
- error_handling.report_error(99)
+ WarningManager.create_warning(
+ "Winding pack cross-section problem",
+ awpc=sctfcoil_module.awpc,
+ dr_tf_wp=tfcoil_variables.dr_tf_wp,
+ )
elif sctfcoil_module.awptf < 0.0e0:
- error_handling.fdiags[0] = sctfcoil_module.awptf
- error_handling.report_error(101)
+ WarningManager.create_warning(
+ "Negative cable space dimension", awptf=sctfcoil_module.awptf
+ )
def tf_res_heating(self) -> None:
"""
@@ -612,38 +617,6 @@ def cpost(
r_tfin_inleg = r_tf_inboard_in + cas_in_th + gr_ins_th
# -#
- # Error traps
- # ------------
- # if rtop <= 0.0e0:
- # error_handling.fdiags[0] = rtop
- # error_handling.report_error(115)
-
- # if ztop <= 0.0e0:
- # error_handling.fdiags[0] = ztop
- # error_handling.report_error(116)
-
- # if rmid <= 0.0e0:
- # error_handling.fdiags[0] = rmid
- # error_handling.report_error(117)
-
- # if build_variables.hmax <= 0.0e0:
- # error_handling.fdiags[0] = build_variables.hmax
- # error_handling.report_error(118)
-
- # if (fcool < 0.0e0) or (fcool > 1.0e0):
- # error_handling.fdiags[0] = fcool
- # error_handling.report_error(119)
-
- # if rtop < rmid:
- # error_handling.fdiags[0] = rtop
- # error_handling.fdiags[1] = rmid
- # error_handling.report_error(120)
-
- # if build_variables.hmax < ztop:
- # error_handling.fdiags[0] = build_variables.hmax
- # error_handling.fdiags[1] = ztop
- # error_handling.report_error(121)
-
# ------------
# Mid-plane area calculations
@@ -688,7 +661,6 @@ def cpost(
vol_case_cp = 0.0e0
vol_gr_ins_cp = 0.0e0
vol_ins_cp = 0.0e0
- # error_handling.report_error(122)
return (
a_cp_cool,
vol_cond_cp,
@@ -748,14 +720,6 @@ def cpost(
r = rc - np.sqrt((rc - rmid) ** 2 - z * z)
- # if r <= 0.0e0:
- # error_handling.fdiags[0] = r
- # error_handling.fdiags[1] = rc
- # error_handling.fdiags[2] = rmid
- # error_handling.fdiags[3] = z
-
- # error_handling.report_error(123)
-
# Insulation cross-sectional area at z
yy_ins[ii] = (
np.pi * ((r_tfin_inleg + ins_th) ** 2 - r_tfin_inleg**2)
diff --git a/process/scan.py b/process/scan.py
index 7f85c67681..d1634b393a 100644
--- a/process/scan.py
+++ b/process/scan.py
@@ -3,6 +3,7 @@
import numpy as np
from tabulate import tabulate
+import process.constraints as constraints
import process.process_output as process_output
from process.caller import write_output_files
from process.exceptions import ProcessValueError
@@ -10,12 +11,10 @@
build_variables,
constants,
constraint_variables,
- constraints,
cost_variables,
cs_fatigue_variables,
current_drive_variables,
divertor_variables,
- error_handling,
fwbs_variables,
global_variables,
heat_transport_variables,
@@ -32,6 +31,7 @@
f2py_compatible_to_string,
string_to_f2py_compatible,
)
+from process.warning_handler import WarningManager
@dataclass
@@ -164,13 +164,10 @@ def run_scan(self):
number of output variable values are written to the MFILE.DAT file at
each scan point, for plotting or other post-processing purposes.
"""
- # Turn off error reporting (until next output)
- error_handling.errors_on = False
-
if scan_module.isweep == 0:
ifail = self.doopt()
write_output_files(models=self.models, ifail=ifail)
- error_handling.show_errors()
+ WarningManager.show_errors(constants.nout)
return
if scan_module.isweep > scan_module.ipnscns:
@@ -203,8 +200,6 @@ def post_optimise(self, ifail: int):
"""
numerics.sqsumsq = (numerics.rcm[: numerics.neqns] ** 2).sum() ** 0.5
- error_handling.errors_on = True
-
process_output.oheadr(constants.nout, "Numerics")
process_output.ocmmnt(
constants.nout, "PROCESS has performed a VMCON (optimisation) run."
@@ -216,8 +211,9 @@ def post_optimise(self, ifail: int):
)
process_output.oblnkl(constants.iotty)
- error_handling.idiags[0] = ifail
- error_handling.report_error(132)
+ WarningManager.create_warning(
+ "Optimisation solver VMCON returns with ifail != 1", ifail=ifail
+ )
self.verror(ifail)
process_output.oblnkl(constants.nout)
@@ -259,8 +255,9 @@ def post_optimise(self, ifail: int):
)
process_output.oblnkl(constants.iotty)
- error_handling.fdiags[0] = numerics.sqsumsq
- error_handling.report_error(134)
+ WarningManager.create_warning(
+ "High final VMCON constraint residues", sqsumsq=numerics.sqsumsq
+ )
process_output.ovarin(
constants.nout, "Number of iteration variables", "(nvar)", numerics.nvar
@@ -429,8 +426,8 @@ def post_optimise(self, ifail: int):
equality_constraint_table.append([
name,
sym[i],
- f"{con2[i]} {f2py_compatible_to_string(lab[i])}",
- f"{err[i]} {f2py_compatible_to_string(lab[i])}",
+ f"{con2[i]} {lab[i]}",
+ f"{err[i]} {lab[i]}",
con1[i],
])
process_output.ovarre(
@@ -468,8 +465,8 @@ def post_optimise(self, ifail: int):
inequality_constraint_table.append([
name,
sym[i],
- f"{con2[i]} {f2py_compatible_to_string(lab[i])}",
- f"{err[i]} {f2py_compatible_to_string(lab[i])}",
+ f"{con2[i]} {lab[i]}",
+ f"{err[i]} {lab[i]}",
])
process_output.ovarre(
constants.mfile,
@@ -656,8 +653,8 @@ def scan_1d(self):
scan_1d_ifail_dict[iscan] = ifail
write_output_files(models=self.models, ifail=ifail)
- error_handling.show_errors()
- error_handling.init_error_handling()
+ WarningManager.show_errors(constants.nout)
+ WarningManager.reinitialise()
# outvar now contains results
self.scan_1d_write_plot()
@@ -709,8 +706,8 @@ def scan_2d(self):
write_output_files(models=self.models, ifail=ifail)
- error_handling.show_errors()
- error_handling.init_error_handling()
+ WarningManager.show_errors(constants.nout)
+ WarningManager.reinitialise()
scan_2d_ifail_list[iscan_1][iscan_2] = ifail
iscan = iscan + 1
diff --git a/process/stellarator.py b/process/stellarator.py
index eb86689989..843c5b7b00 100644
--- a/process/stellarator.py
+++ b/process/stellarator.py
@@ -19,7 +19,6 @@
cost_variables,
current_drive_variables,
divertor_variables,
- error_handling,
fwbs_variables,
global_variables,
heat_transport_variables,
@@ -43,6 +42,7 @@
from process.physics import rether
from process.stellarator_config import load_stellarator_config
from process.utilities.f2py_string_patch import f2py_compatible_to_string
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
# Logging handler for console output
@@ -3378,11 +3378,13 @@ def intersect(self, x1, y1, x2, y2, xin):
xmax = min(np.max(x1), np.amax(x2))
if xmin >= xmax:
- error_handling.fdiags[0] = np.amin(x1)
- error_handling.fdiags[1] = np.amin(x2)
- error_handling.fdiags[2] = np.amax(x1)
- error_handling.fdiags[3] = np.amax(x2)
- error_handling.report_error(111)
+ WarningManager.create_warning(
+ "X ranges not overlapping",
+ x1_min=np.amin(x1),
+ x1_max=np.amax(x1),
+ x2_min=np.amin(x2),
+ x2_max=np.amax(x2),
+ )
# Ensure input guess for x is within this range
@@ -3430,20 +3432,24 @@ def intersect(self, x1, y1, x2, y2, xin):
x = x - 2.0e0 * dx * y / (yright - yleft)
if x < xmin:
- error_handling.fdiags[0] = x
- error_handling.fdiags[1] = xmin
- error_handling.report_error(112)
+ WarningManager.create_warning(
+ "X has dropped below Xmin; X has been set equal to Xmin",
+ x=x,
+ xmin=xmin,
+ )
x = xmin
break
if x > xmax:
- error_handling.fdiags[0] = x
- error_handling.fdiags[1] = xmax
- error_handling.report_error(113)
+ WarningManager.create_warning(
+ "X has risen above Xmax; X has been set equal to Xmax",
+ x=x,
+ xmax=xmax,
+ )
x = xmax
break
else:
- error_handling.report_error(114)
+ WarningManager.create_warning("Convergence too slow; X may be wrong")
return x
diff --git a/process/superconducting_tf_coil.py b/process/superconducting_tf_coil.py
index 4071bf282a..298d1f2818 100644
--- a/process/superconducting_tf_coil.py
+++ b/process/superconducting_tf_coil.py
@@ -11,7 +11,6 @@
build_variables,
constants,
divertor_variables,
- error_handling,
global_variables,
numerics,
pfcoil_variables,
@@ -22,6 +21,7 @@
)
from process.tf_coil import TFCoil
from process.utilities.f2py_string_patch import f2py_compatible_to_string
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
@@ -260,7 +260,9 @@ def run(self, output: bool):
)
except ValueError as e:
if e.args[1] == 245 and e.args[2] == 0:
- error_handling.report_error(245)
+ WarningManager.create_warning(
+ "Invalid stress model (r_tf_inboard = 0), stress constraint switched off"
+ )
tfcoil_variables.sig_tf_case = 0.0e0
tfcoil_variables.sig_tf_wp = 0.0e0
peaktfflag = 0
@@ -902,8 +904,11 @@ def supercon(
tc0m = 16.06e0
# If strain limit achieved, throw a warning and use the lower strain
if abs(strain) > 0.5e-2:
- error_handling.fdiags[0] = strain
- error_handling.report_error(261)
+ WarningManager.create_warning(
+ "TF strain was outside the region of applicability. Used lower strain",
+ strain=strain,
+ )
+
strain = np.sign(strain) * 0.5e-2
# j_crit_sc returned by superconductors.itersc is the critical current density in the
@@ -955,8 +960,10 @@ def supercon(
tc0m = tcritsc
# If strain limit achieved, throw a warning and use the lower strain
if abs(strain) > 0.5e-2:
- error_handling.fdiags[0] = strain
- error_handling.report_error(261)
+ WarningManager.create_warning(
+ "TF strain was outside the region of applicability. Used lower strain",
+ strain=strain,
+ )
strain = np.sign(strain) * 0.5e-2
j_crit_sc, _, _ = superconductors.itersc(thelium, bmax, strain, bc20m, tc0m)
@@ -974,8 +981,10 @@ def supercon(
tc0m = 16.06e0
# If strain limit achieved, throw a warning and use the lower strain
if abs(strain) > 0.5e-2:
- error_handling.fdiags[0] = strain
- error_handling.report_error(261)
+ WarningManager.create_warning(
+ "TF strain was outside the region of applicability. Used lower strain",
+ strain=strain,
+ )
strain = np.sign(strain) * 0.5e-2
# j_crit_sc returned by superconductors.itersc is the critical current density in the
@@ -1015,8 +1024,10 @@ def supercon(
tc0m = 185
# If strain limit achieved, throw a warning and use the lower strain
if abs(strain) > 0.7e-2:
- error_handling.fdiags[0] = strain
- error_handling.report_error(261)
+ WarningManager.create_warning(
+ "TF strain was outside the region of applicability. Used lower strain",
+ strain=strain,
+ )
strain = np.sign(strain) * 0.7e-2
j_crit_sc, _, _ = superconductors.gl_rebco(
@@ -1039,8 +1050,10 @@ def supercon(
tc0m = 92
# If strain limit achieved, throw a warning and use the lower strain
if abs(strain) > 0.7e-2:
- error_handling.fdiags[0] = strain
- error_handling.report_error(261)
+ WarningManager.create_warning(
+ "TF strain was outside the region of applicability. Used lower strain",
+ strain=strain,
+ )
strain = np.sign(strain) * 0.7e-2
# 'high current density' as per parameterisation described in Wolf,
@@ -1087,7 +1100,9 @@ def supercon(
# REBCO measurements from 2 T to 14 T, extrapolating outside this
if (isumat == 8) and (tfcoil_variables.bmaxtfrp >= 14):
- error_handling.report_error(266)
+ WarningManager.create_warning(
+ "Field on superconductor > 14 T (outside of interpolation range)"
+ )
# Temperature margin (already calculated in superconductors.bi2212 for isumat=2)
@@ -1742,15 +1757,17 @@ def sc_tf_internal_geom(self, i_tf_wp_geom, i_tf_case_geom, i_tf_turns_integer):
or sctfcoil_module.a_tf_ins <= 0.0e0
or sctfcoil_module.f_tf_ins <= 0.0e0
):
- error_handling.fdiags[0] = tfcoil_variables.acond
- error_handling.fdiags[1] = tfcoil_variables.avwp
- error_handling.fdiags[2] = tfcoil_variables.a_tf_coil_wp_turn_insulation
- error_handling.fdiags[3] = tfcoil_variables.aswp
- error_handling.fdiags[4] = sctfcoil_module.a_tf_steel
- error_handling.fdiags[5] = sctfcoil_module.f_tf_steel
- error_handling.fdiags[6] = sctfcoil_module.a_tf_ins
- error_handling.fdiags[7] = sctfcoil_module.f_tf_ins
- error_handling.report_error(276)
+ WarningManager.create_warning(
+ "One of the areas or fractions is negative in the internal SC TF coil geometry",
+ acond=tfcoil_variables.acond,
+ avwp=tfcoil_variables.avwp,
+ a_tf_coil_wp_turn_insulation=tfcoil_variables.a_tf_coil_wp_turn_insulation,
+ aswp=tfcoil_variables.aswp,
+ a_tf_steel=sctfcoil_module.a_tf_steel,
+ f_tf_steel=sctfcoil_module.f_tf_steel,
+ a_tf_ins=sctfcoil_module.a_tf_ins,
+ f_tf_ins=sctfcoil_module.f_tf_ins,
+ )
def tf_wp_geom(self, i_tf_wp_geom):
"""
@@ -1911,9 +1928,11 @@ def tf_wp_geom(self, i_tf_wp_geom):
# Negative WP area error reporting
if sctfcoil_module.awptf <= 0.0e0 or sctfcoil_module.awpc <= 0.0e0:
- error_handling.fdiags[0] = sctfcoil_module.awptf
- error_handling.fdiags[1] = sctfcoil_module.awpc
- error_handling.report_error(99)
+ WarningManager.create_warning(
+ "Winding pack cross-section problem",
+ awptf=sctfcoil_module.awptf,
+ awpc=sctfcoil_module.awpc,
+ )
def tf_case_geom(self, i_tf_wp_geom, i_tf_case_geom):
"""
@@ -1952,9 +1971,11 @@ def tf_case_geom(self, i_tf_wp_geom, i_tf_case_geom):
# Report error if the casing area is negative
if tfcoil_variables.acasetf <= 0.0e0 or tfcoil_variables.acasetfo <= 0.0e0:
- error_handling.fdiags[0] = tfcoil_variables.acasetf
- error_handling.fdiags[1] = tfcoil_variables.acasetfo
- error_handling.report_error(99)
+ WarningManager.create_warning(
+ "Winding pack cross-section problem",
+ acasetf=tfcoil_variables.acasetf,
+ acasetfo=tfcoil_variables.acasetfo,
+ )
# Average lateral casing thickness
# --------------
@@ -1998,10 +2019,12 @@ def tf_integer_turn_geom(self, n_layer, n_pancake, thwcndut, thicndut):
) / n_layer
if sctfcoil_module.t_turn_radial <= (2.0e0 * thicndut + 2.0e0 * thwcndut):
- error_handling.fdiags[0] = sctfcoil_module.t_turn_radial
- error_handling.fdiags[1] = thicndut
- error_handling.fdiags[2] = thwcndut
- error_handling.report_error(100)
+ WarningManager.create_warning(
+ "Negative cable space dimension; reduce conduit thicknesses or raise cpttf",
+ t_turn_radial=sctfcoil_module.t_turn_radial,
+ thicndut=thicndut,
+ thwcndut=thwcndut,
+ )
# Toroidal turn dimension [m]
sctfcoil_module.t_turn_toroidal = (
@@ -2010,10 +2033,12 @@ def tf_integer_turn_geom(self, n_layer, n_pancake, thwcndut, thicndut):
) / n_pancake
if sctfcoil_module.t_turn_toroidal <= (2.0e0 * thicndut + 2.0e0 * thwcndut):
- error_handling.fdiags[0] = sctfcoil_module.t_turn_toroidal
- error_handling.fdiags[1] = thicndut
- error_handling.fdiags[2] = thwcndut
- error_handling.report_error(100)
+ WarningManager.create_warning(
+ "Negative cable space dimension; reduce conduit thicknesses or raise cpttf",
+ t_turn_radial=sctfcoil_module.t_turn_toroidal,
+ thicndut=thicndut,
+ thwcndut=thwcndut,
+ )
tfcoil_variables.t_turn_tf = np.sqrt(
sctfcoil_module.t_turn_radial * sctfcoil_module.t_turn_toroidal
@@ -2057,15 +2082,19 @@ def tf_integer_turn_geom(self, n_layer, n_pancake, thwcndut, thicndut):
if (sctfcoil_module.t_cable_radial < 0.0e0) or (
sctfcoil_module.t_cable_toroidal < 0.0e0
):
- error_handling.fdiags[0] = acstf
- error_handling.fdiags[1] = sctfcoil_module.t_cable_radial
- error_handling.fdiags[2] = sctfcoil_module.t_cable_toroidal
- error_handling.report_error(101)
+ WarningManager.create_warning(
+ "Negative cable space dimension",
+ acstf=acstf,
+ t_cable_radial=sctfcoil_module.t_cable_radial,
+ t_cable_toroidal=sctfcoil_module.t_cable_toroidal,
+ )
else:
- error_handling.fdiags[0] = acstf
- error_handling.fdiags[1] = sctfcoil_module.t_cable_radial
- error_handling.fdiags[1] = sctfcoil_module.t_cable_toroidal
- error_handling.report_error(102)
+ WarningManager.create_warning(
+ "Cable space area problem; artificially set rounded corner radius to 0",
+ acstf=acstf,
+ t_cable_radial=sctfcoil_module.t_cable_radial,
+ t_cable_toroidal=sctfcoil_module.t_cable_toroidal,
+ )
sctfcoil_module.rbcndut = 0.0e0
acstf = (
sctfcoil_module.t_cable_radial * sctfcoil_module.t_cable_toroidal
@@ -3124,9 +3153,13 @@ def outtf(self, peaktfflag):
po.osubhd(self.outfile, "Ripple information:")
if tfcoil_variables.i_tf_shape == 1:
if peaktfflag == 1:
- error_handling.report_error(144)
+ WarningManager.create_warning(
+ "(TF coil peak field calculation) Winding pack width out of fitted range"
+ )
elif peaktfflag == 2:
- error_handling.report_error(145)
+ WarningManager.create_warning(
+ "(TF coil peak field calculation) Winding pack radial thickness out of fitted range"
+ )
po.ovarre(
self.outfile,
@@ -3579,13 +3612,17 @@ def tf_averaged_turn_geom(self, j_tf_wp, thwcndut, thicndut, i_tf_sc_mat):
if acstf <= 0.0e0:
if tfcoil_variables.t_conductor < 0.0e0:
- error_handling.fdiags[0] = acstf
- error_handling.fdiags[1] = sctfcoil_module.t_cable
- error_handling.report_error(101)
+ WarningManager.create_warning(
+ "Negative cable space dimension",
+ acstf=acstf,
+ t_cable=sctfcoil_module.t_cable,
+ )
else:
- error_handling.fdiags[0] = acstf
- error_handling.fdiags[1] = sctfcoil_module.t_cable
- error_handling.report_error(102)
+ WarningManager.create_warning(
+ "Cable space area problem; artificially set rounded corner radius to 0",
+ acstf=acstf,
+ t_cable=sctfcoil_module.t_cable,
+ )
rbcndut = 0.0e0
acstf = sctfcoil_module.t_cable**2
diff --git a/process/superconductors.py b/process/superconductors.py
index 06fa825623..8999451738 100644
--- a/process/superconductors.py
+++ b/process/superconductors.py
@@ -4,8 +4,8 @@
from scipy import optimize
from process.exceptions import ProcessValueError
-from process.fortran import error_handling as eh
from process.fortran import rebco_variables
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
@@ -545,18 +545,20 @@ def bottura_scaling(
tc0eps = tc0max * strfun ** (1 / 3)
if temperature / tc0eps >= 1.0:
- eh.fdiags[0] = temperature
- eh.fdiags[1] = tc0eps
- eh.report_error(159)
+ WarningManager.create_warning(
+ "Reduced temperature t artificially lowered",
+ temperature=temperature,
+ tc0eps=tc0eps,
+ )
# Reduced temperature at zero field, corrected for strain
# t > 1 is permitted, indicating the temperature is above the critical value at zero field.
t = temperature / tc0eps
if bmax / bc20eps >= 1.0:
- eh.fdiags[0] = bmax
- eh.fdiags[1] = bc20eps
- eh.report_error(160)
+ WarningManager.create_warning(
+ "Reduced field bzero artificially lowered", bmax=bmax, bc20eps=bc20eps
+ )
# Reduced field at zero temperature, taking account of strain
bzero = bmax / bc20eps
diff --git a/process/tf_coil.py b/process/tf_coil.py
index 7c2bb767dc..d4f0386bd7 100644
--- a/process/tf_coil.py
+++ b/process/tf_coil.py
@@ -11,7 +11,6 @@
from process.fortran import (
build_variables,
constants,
- error_handling,
fwbs_variables,
global_variables,
pfcoil_variables,
@@ -20,13 +19,13 @@
tfcoil_variables,
)
from process.fortran import build_variables as bv
-from process.fortran import error_handling as eh
from process.fortran import fwbs_variables as fwbsv
from process.fortran import tfcoil_variables as tfv
from process.utilities.f2py_string_patch import (
f2py_compatible_to_string,
string_to_f2py_compatible,
)
+from process.warning_handler import WarningManager
RMU0 = constants.rmu0
@@ -247,7 +246,9 @@ def run(self, output):
)
except ValueError as e:
if e.args[1] == 245 and e.args[2] == 0:
- error_handling.report_error(245)
+ WarningManager.create_warning(
+ "Invalid stress model (r_tf_inboard = 0), stress constraint switched off"
+ )
tfcoil_variables.sig_tf_case = 0.0e0
tfcoil_variables.sig_tf_wp = 0.0e0
@@ -719,12 +720,18 @@ def cntrpst(self):
# If the average conductor temperature difference is negative, set it to 0
if dtcncpav < 0.0e0:
- eh.report_error(249)
+ WarningManager.create_warning(
+ "Negative conductor average temperature difference, set to 0",
+ dtcncpav=dtcncpav,
+ )
dtcncpav = 0.0e0
# If the average conductor temperature difference is negative, set it to 0
if dtconcpmx < 0.0e0:
- eh.report_error(250)
+ WarningManager.create_warning(
+ "Negative conductor peak temperature difference, set to 0",
+ dtconcpmx=dtconcpmx,
+ )
dtconcpmx = 0.0e0
# Average conductor temperature
@@ -1008,8 +1015,10 @@ def he_density(temp: float) -> float:
# Fit range validation
if temp < 4.0e0 or temp > 50.0e0:
- eh.fdiags[0] = temp
- eh.report_error(257)
+ WarningManager.create_warning(
+ "Helium temperature out of helium property fiting range [4-50] K",
+ temp=temp,
+ )
# Oder 3 polynomial fit
if temp < 29.5e0:
@@ -1047,8 +1056,10 @@ def he_cp(temp: float) -> float:
# Fit range validation
if temp < 4.0e0 or temp > 50.0e0:
- eh.fdiags[0] = temp
- eh.report_error(257)
+ WarningManager.create_warning(
+ "Helium temperature out of helium property fiting range [4-50] K",
+ temp=temp,
+ )
# Order 3 polynomial fit in [4-30] K on the dimenion [K/(g.K)]
if temp < 29.5e0:
@@ -1088,8 +1099,10 @@ def he_visco(temp: float) -> float:
"""
if temp < 4.0e0 or temp > 50.0e0:
- eh.fdiags[0] = temp
- eh.report_error(257)
+ WarningManager.create_warning(
+ "Helium temperature out of helium property fiting range [4-50] K",
+ temp=temp,
+ )
# Order 4 polynomial exponential fit in [4-25] K
if temp < 22.5e0:
@@ -1127,8 +1140,10 @@ def he_th_cond(temp: float) -> float:
# Fit range validation
if temp < 4.0e0 or temp > 50.0e0:
- eh.fdiags[0] = temp
- eh.report_error(257)
+ WarningManager.create_warning(
+ "Helium temperature out of helium property fiting range [4-50] K",
+ temp=temp,
+ )
# Order 4 polynomial fit
if temp < 24.0e0:
@@ -1177,8 +1192,10 @@ def al_th_cond(temp: float) -> float:
# Fiting range verification
if temp < 15.0e0 or temp > 150.0e0:
- eh.fdiags[0] = temp
- eh.report_error(258)
+ WarningManager.create_warning(
+ "Aluminium temperature out of the th conductivity fit range [15-60] K",
+ temp=temp,
+ )
# fit 15 < T < 60 K (order 3 poly)
if temp < 60.0e0:
diff --git a/process/uncertainties/evaluate_uncertainties.py b/process/uncertainties/evaluate_uncertainties.py
index 378ac952b6..ec953e28db 100644
--- a/process/uncertainties/evaluate_uncertainties.py
+++ b/process/uncertainties/evaluate_uncertainties.py
@@ -35,7 +35,6 @@
from process.io.in_dat import InDat
from process.io.process_config import UncertaintiesConfig
from process.io.process_funcs import (
- check_input_error,
get_neqns_itervars,
get_variable_range,
no_unfeasible_mfile,
@@ -110,9 +109,6 @@ def run_monte_carlo(args):
input_path = Path(config.wdir) / "IN.DAT"
config.run_process(input_path)
- # Check for produced MFILE.DAT in working dir
- check_input_error(wdir=config.wdir)
-
if not process_stopped():
no_unfeasible = no_unfeasible_mfile()
diff --git a/process/utilities/errorlist.json b/process/utilities/errorlist.json
deleted file mode 100644
index 63e10f361e..0000000000
--- a/process/utilities/errorlist.json
+++ /dev/null
@@ -1,1459 +0,0 @@
-{
- "version": [
- "PROCESS version 2.1.1"
- ],
- "comment1": [
- "List of error messages used via error_handling.f90"
- ],
- "comment2": [
- "Increment n_errortypes if an error is added to this list"
- ],
- "n_errortypes": 289,
- "errors": [
- {
- "no": 1,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 7 if i_plasma_ignited=1"
- },
- {
- "no": 2,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 3,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 4,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 28 if i_plasma_ignited=1"
- },
- {
- "no": 5,
- "level": 3,
- "message": "CONSTRAINTS: temp_fw_peak = 0 ==> i_pulsed_plant=0; do not use constraint 39 if i_pulsed_plant=0"
- },
- {
- "no": 6,
- "level": 3,
- "message": "CONSTRAINTS: tcycmn = 0 ==> i_pulsed_plant=0; do not use constraint 42 if i_pulsed_plant=0"
- },
- {
- "no": 7,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 43 if itart=0"
- },
- {
- "no": 8,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 44 if itart=0"
- },
- {
- "no": 9,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 45 if itart=0"
- },
- {
- "no": 10,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 46 if itart=0"
- },
- {
- "no": 11,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 12,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 50 if ife=0"
- },
- {
- "no": 13,
- "level": 3,
- "message": "CONSTRAINTS: No such constraint equation number..."
- },
- {
- "no": 14,
- "level": 3,
- "message": "CONSTRAINTS: NaN/infty error for constraint equation..."
- },
- {
- "no": 15,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 16,
- "level": 2,
- "message": "LHRAD: LH penetration radius not found after lapno iterations, using 0.8*rminor"
- },
- {
- "no": 17,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 18,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 19,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 20,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 21,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 22,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 23,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 24,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 25,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 26,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 27,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 28,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 29,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 30,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 31,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 32,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 33,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 34,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 35,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 36,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 37,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 38,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 39,
- "level": 1,
- "message": "REMOVED"
- },
- {
- "no": 40,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 41,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 42,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 43,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 44,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 45,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 46,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 47,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 48,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 49,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 50,
- "level": 3,
- "message": "REMOVED"
- },
- {
- "no": 51,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 52,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 53,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 54,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 55,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 56,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 57,
- "level": 3,
- "message": "CONVXC: Illegal iteration variable number"
- },
- {
- "no": 58,
- "level": 3,
- "message": "CONVXC: Iteration variable is zero; change its initial value or lower bound"
- },
- {
- "no": 59,
- "level": 3,
- "message": "CONVXC: NaN error for iteration variable"
- },
- {
- "no": 60,
- "level": 3,
- "message": "CONVXC: scale(i) = 0 for iteration variable"
- },
- {
- "no": 61,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 62,
- "level": 1,
- "message": "RADIALB: Ripple result may be inaccurate, as the fit has been extrapolated"
- },
- {
- "no": 63,
- "level": 2,
- "message": "PORTSZ: Max beam tangency radius set =0 temporarily; change beamwd"
- },
- {
- "no": 64,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 65,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 66,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 67,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 68,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 69,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 70,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 71,
- "level": 2,
- "message": "PFCOIL: OH coil not present; check volt-second calculations..."
- },
- {
- "no": 72,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 73,
- "level": 2,
- "message": "INDUCT: Max no. of segments noh for OH coil > nohmax; increase dr_cs lower bound"
- },
- {
- "no": 74,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 75,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 76,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 77,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 78,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 79,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 80,
- "level": 2,
- "message": "calculate_density_limit: qcyl < 4/3; dlimit(4) set to zero; model 5 will be enforced instead"
- },
- {
- "no": 81,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 82,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 83,
- "level": 2,
- "message": "POHM: Negative plasma resistance res_plasma"
- },
- {
- "no": 84,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 85,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 86,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 87,
- "level": 2,
- "message": "OUTPLAS: Possible problem with high radiation power, forcing p_plasma_separatrix_mw to odd values"
- },
- {
- "no": 88,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 89,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 90,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 91,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 92,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 93,
- "level": 2,
- "message": "BURN: Negative burn time available; reduce t_fusion_ramp or raise PF coil V-s capability"
- },
- {
- "no": 94,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 95,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 96,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 97,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 98,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 99,
- "level": 2,
- "message": "SCTFCOIL: Winding pack cross-section problem..."
- },
- {
- "no": 100,
- "level": 2,
- "message": "SCTFCOIL: Negative cable space dimension; reduce conduit thicknesses or raise cpttf"
- },
- {
- "no": 101,
- "level": 2,
- "message": "SCTFCOIL: Negative cable space dimension"
- },
- {
- "no": 102,
- "level": 2,
- "message": "SCTFCOIL: Cable space area problem; artificially set rounded corner radius to 0"
- },
- {
- "no": 103,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 104,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 105,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 106,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 107,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 108,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 109,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 110,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 111,
- "level": 2,
- "message": "INTERSECT: X ranges not overlapping"
- },
- {
- "no": 112,
- "level": 2,
- "message": "INTERSECT: X has dropped below Xmin; X has been set equal to Xmin"
- },
- {
- "no": 113,
- "level": 2,
- "message": "INTERSECT: X has risen above Xmax; X has been set equal to Xmax"
- },
- {
- "no": 114,
- "level": 2,
- "message": "INTERSECT: Convergence too slow; X may be wrong..."
- },
- {
- "no": 115,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 116,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 117,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 118,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 119,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 120,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 121,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 122,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 123,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 124,
- "level": 2,
- "message": "VACUUM: Newton's method not converging; check fusion power, te"
- },
- {
- "no": 125,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 126,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 127,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 128,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 129,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 130,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 131,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 132,
- "level": 2,
- "message": "DOOPT: Optimisation solver VMCON returns with ifail /= 1"
- },
- {
- "no": 133,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 134,
- "level": 2,
- "message": "DOOPT: High final VMCON constraint residues"
- },
- {
- "no": 135,
- "level": 1,
- "message": "OUTPF: CS not using max current density: further optimisation may be possible"
- },
- {
- "no": 136,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 137,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 138,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 139,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 140,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 141,
- "level": 1,
- "message": "RADIALB: (TF coil ripple calculation) Dimensionless coil width X out of fitted range"
- },
- {
- "no": 142,
- "level": 1,
- "message": "RADIALB: (TF coil ripple calculation) No of TF coils not between 16 and 20 inclusive"
- },
- {
- "no": 143,
- "level": 1,
- "message": "RADIALB: (TF coil ripple calculation) (R+a)/rtot out of fitted range"
- },
- {
- "no": 144,
- "level": 1,
- "message": "OUTTF: (TF coil peak field calculation) Winding pack width out of fitted range"
- },
- {
- "no": 145,
- "level": 1,
- "message": "OUTTF: (TF coil peak field calculation) Winding pack radial thickness out of fitted range"
- },
- {
- "no": 146,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 147,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 148,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 149,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 150,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 151,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 152,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 153,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 154,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 155,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 156,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 157,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 158,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 159,
- "level": 2,
- "message": "ITERSC: Reduced temperature t artificially lowered"
- },
- {
- "no": 160,
- "level": 2,
- "message": "ITERSC: Reduced field bzero artificially lowered"
- },
- {
- "no": 161,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 162,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 163,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 164,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 165,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 166,
- "level": 2,
- "message": "PLANT_THERMAL_EFFICIENCY: Turbine temperature tturb out of range of validity"
- },
- {
- "no": 167,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 168,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 169,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 170,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 171,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 172,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 173,
- "level": 3,
- "message": "CONSTRAINTS: Do not use constraint 55 if i_blanket_type != 2"
- },
- {
- "no": 174,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 175,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 176,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 177,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 178,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 179,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 180,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 181,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 182,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 183,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 184,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 185,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 186,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 187,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 188,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 189,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 190,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 191,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 192,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 193,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 194,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 195,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 196,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 197,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 198,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 199,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 200,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 201,
- "level": 1,
- "message": "LH THRESHOLD: bt outside Snipes 2000 fitted range"
- },
- {
- "no": 202,
- "level": 1,
- "message": "LH THRESHOLD: rminor outside Snipes 2000 fitted range"
- },
- {
- "no": 203,
- "level": 1,
- "message": "LH THRESHOLD: rmajor outside Snipes 2000 fitted range"
- },
- {
- "no": 204,
- "level": 1,
- "message": "LH THRESHOLD: dnla outside Snipes 2000 fitted range"
- },
- {
- "no": 205,
- "level": 1,
- "message": "LH THRESHOLD: kappa outside Snipes 2000 fitted range"
- },
- {
- "no": 206,
- "level": 1,
- "message": "LH THRESHOLD: triang outside Snipes 2000 fitted range"
- },
- {
- "no": 207,
- "level": 1,
- "message": "LH THRESHOLD: Closed divertor only. Limited data used in Snipes fit"
- },
- {
- "no": 208,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 209,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 210,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 211,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 212,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 213,
- "level": 2,
- "message": "NCORE: Central density is less than pedestal density"
- },
- {
- "no": 214,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 215,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 216,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 217,
- "level": 2,
- "message": "REINKE IMPURITY MODEL: fzmin is greater than or equal to 1.0, this is at least notable"
- },
- {
- "no": 218,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 219,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 220,
- "level": 1,
- "message": "REMOVED"
- },
- {
- "no": 221,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 222,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 223,
- "level": 2,
- "message": "[fw.f90][fw_temp] NaN first wall temperature "
- },
- {
- "no": 224,
- "level": 2,
- "message": "[fw.f90][fw_temp] First wall temperature (temp_k) out of range : [100-1500] K"
- },
- {
- "no": 225,
- "level": 2,
- "message": "[fw.f90][heat_transfer] Reynolds number out of range : [3e3-5000e3]"
- },
- {
- "no": 226,
- "level": 2,
- "message": "[fw.f90][heat_transfer] Prandtl number out of range : [0.5-2000]"
- },
- {
- "no": 227,
- "level": 2,
- "message": "[fw.f90][heat_transfer] Negative Darcy friction factor (f)"
- },
- {
- "no": 228,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 229,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 230,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 231,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 232,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 233,
- "level": 1,
- "message": "REMOVED"
- },
- {
- "no": 234,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 235,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 236,
- "level": 3,
- "message": "CHECK: The 1/R toroidal B field dependency constraint is being depreciated."
- },
- {
- "no": 237,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 238,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 239,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 240,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 241,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 242,
- "level": 2,
- "message": "PHYSICS: Bootstrap fraction upper limit enforced"
- },
- {
- "no": 243,
- "level": 2,
- "message": "PHYSICS: Predicted plasma driven current is more than upper limit on non-inductive fraction."
- },
- {
- "no": 244,
- "level": 2,
- "message": "PHYSICS: Diamagnetic fraction is more than 1%, but not calculated. Consider using i_diamagnetic_current=2 and i_pfirsch_schluter_current=1."
- },
- {
- "no": 245,
- "level": 2,
- "message": "[sctfcoil] Invalid stress model (r_tf_inboard = 0), stress constraint switched off"
- },
- {
- "no": 246,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 247,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 248,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 249,
- "level": 2,
- "message": "[tfcoil][cntrpst] Negative conductor average temperature difference, set to 0"
- },
- {
- "no": 250,
- "level": 2,
- "message": "[tfcoil][cntrpst] Negative conductor peak temperature difference, set to 0"
- },
- {
- "no": 251,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 252,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 253,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 254,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 255,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 256,
- "level": 2,
- "message": "CHECK: Top CP radius larger that its value determined with plasma shape"
- },
- {
- "no": 257,
- "level": 2,
- "message": "[tfcoil] : Helium temperature out of helium property fiting range [4-50] K"
- },
- {
- "no": 258,
- "level": 2,
- "message": "[tfcoil] : Aluminium temperature out of the th conductivity fit range [15-60] K"
- },
- {
- "no": 259,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 260,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 261,
- "level": 2,
- "message": "[tfcoil]: TF strain was outside the region of applicability. Used lower strain"
- },
- {
- "no": 262,
- "level": 2,
- "message": "REBCO : Non physical strain used in CS. Use superconductor strain < +/- 0.7% "
- },
- {
- "no": 263,
- "level": 2,
- "message": "REBCO : Non physical strain used in PF. Use superconductor strain < +/- 0.7% "
- },
- {
- "no": 264,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 265,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 266,
- "level": 2,
- "message": "REBCO : Field on superconductor > 14 T (outside of interpolation range)"
- },
- {
- "no": 267,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 268,
- "level": 2,
- "message": "[radialb]: TF CP top radius (r_cp_top) replaced by 1.01*r_tf_inboard_out -> potential top rbuild issue"
- },
- {
- "no": 269,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 270,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 271,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 272,
- "level": 2,
- "message": "REMOVED"
- },
- {
- "no": 273,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 274,
- "level": 2,
- "message": "[HCPB]: Blanket heating is <1 MW or NaN. Is something wrong?"
- },
- {
- "no": 275,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 276,
- "level": 2,
- "message": "CHECK: One of the areas or fractions is negative in the internal SC TF coil geometry"
- },
- {
- "no": 277,
- "level": 2,
- "message": "CHECK: One or more collision between TF and PF coils. Check PF placement."
- },
- {
- "no": 278,
- "level": 2,
- "message": "CHECK: Your blanket modules are too small for the Liquid Metal pipes."
- },
- {
- "no": 279,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 280,
- "level": 2,
- "message": "CHECK: Outside temperature limit for one or more liquid metal breeder properties."
- },
- {
- "no": 281,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 282,
- "level": 2,
- "message": "CHECK: ncore is going negative when solving. Please raise the value of dene and or its lower limit."
- },
- {
- "no": 283,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 284,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 285,
- "level": 2,
- "message": "[tfcoil]: fiooic shouldn't be above 0.7 for engineering reliability."
- },
- {
- "no": 286,
- "level": 2,
- "message": "[pfcoil][cntrpost]: fjohc shouldn't be above 0.7 for engineering reliability."
- },
- {
- "no": 287,
- "level": 2,
- "message": "[pfcoil][cntrpost]: fjohc0 shouldn't be above 0.7 for engineering reliability."
- },
- {
- "no": 288,
- "level": 0,
- "message": "REMOVED"
- },
- {
- "no": 289,
- "level": 3,
- "message": "CHECK: Constraint equation 22 has been deprecated."
- }
- ]
-}
diff --git a/process/vacuum.py b/process/vacuum.py
index 1ec7163b05..d1873c02aa 100644
--- a/process/vacuum.py
+++ b/process/vacuum.py
@@ -6,7 +6,6 @@
from process import process_output as po
from process.fortran import build_variables as buv
from process.fortran import constants
-from process.fortran import error_handling as eh
from process.fortran import physics_variables as pv
from process.fortran import tfcoil_variables as tfv
from process.fortran import times_variables as tv
@@ -15,6 +14,7 @@
f2py_compatible_to_string,
string_to_f2py_compatible,
)
+from process.warning_handler import WarningManager
logger = logging.getLogger(__name__)
@@ -448,9 +448,11 @@ def vacuum(
break
else:
- eh.fdiags[0] = pv.fusion_power
- eh.fdiags[1] = pv.te
- eh.report_error(124)
+ WarningManager.create_warning(
+ "Newton's method not converging; check fusion power, te",
+ fusion_power=pv.fusion_power,
+ te=pv.te,
+ )
theta = math.pi / ntf
diff --git a/process/warning_handler.py b/process/warning_handler.py
new file mode 100644
index 0000000000..b0c6b33f28
--- /dev/null
+++ b/process/warning_handler.py
@@ -0,0 +1,65 @@
+import inspect
+from dataclasses import dataclass
+from typing import ClassVar
+
+import process.process_output as process_output
+from process.fortran import constants
+
+
+class ProcessUserWarning(UserWarning):
+ pass
+
+
+@dataclass
+class ProcessWarning:
+ msg: str
+ location: str
+
+
+class WarningManager:
+ _warnings: ClassVar[list[ProcessWarning]] = []
+
+ def __init__(self):
+ raise NotImplementedError(f"{self.__class__.__name__} cannot be instantiated.")
+
+ @classmethod
+ def create_warning(cls, msg: str, **diagnostics):
+ caller = inspect.stack()[1]
+ warning = ProcessWarning(
+ f"{msg}"
+ + ("\n\t" if diagnostics else "")
+ + "\n\t".join([f"\t{d}: {v!r}" for d, v in [*diagnostics.items()]]),
+ f"{caller.filename} @ line {caller.lineno}",
+ )
+
+ # will not add the warning again if it has the same message and line number
+ if warning not in cls._warnings:
+ cls._warnings.append(warning)
+
+ @classmethod
+ def warnings(cls):
+ return cls._warnings
+
+ @classmethod
+ def reinitialise(cls):
+ cls._warnings = []
+
+ @classmethod
+ def show_errors(cls, file_unit: int):
+ warning_string = (
+ "******************************************** Errors and Warnings *********************************************"
+ f"\n{cls.warning_string()}"
+ )
+ print(warning_string)
+ process_output.write(file_unit, warning_string)
+ process_output.ovarre(
+ constants.mfile,
+ "Error status",
+ "error_status",
+ 2 if len(WarningManager.warnings()) > 0 else 0,
+ "OP",
+ )
+
+ @classmethod
+ def warning_string(cls):
+ return "\n\n".join([f"({w.location}) {w.msg}" for w in cls._warnings])
diff --git a/source/fortran/constraint_equations.f90 b/source/fortran/constraint_equations.f90
deleted file mode 100755
index 1a0bb7e502..0000000000
--- a/source/fortran/constraint_equations.f90
+++ /dev/null
@@ -1,3442 +0,0 @@
-! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-module constraints
- !! author: J Morris
- !!
- !! Module defining the constraint equations and the routine that evaluates the
- !! constraint equations.
-
- ! Import modules
-#ifndef dp
- use, intrinsic :: iso_fortran_env, only: dp=>real64
-#endif
- use error_handling, only: report_error, idiags, fdiags
-
- implicit none
-
- public :: constraint_eqns
-
-! type constraint_args_type
-! real(dp) :: cc
-! !! Residual error in constraint equation
-! real(dp) :: con
-! !! constraint value for constraint equation in physical units
-! real(dp) :: err
-! !! residual error in constraint equation in physical units
-! character(len=1) :: symbol
-! !! `=<`, `>`, `<` symbol for constraint equation denoting its type
-! character(len=10) :: units
-! !! constraint units in constraint equation
-! end type
-
-contains
-
- subroutine constraint_eqns(m,ieqn,cc,con,err,symbol,units)
- !! Routine that formulates the constraint equations
- !!
- !! **author: P J Knight** (UKAEA)
- !!
- !! **author: J Morris** (UKAEA)
- !!
- !! if `ieqn` is zero or negative, evaluate all the constraint equations, otherwise
- !! evaluate only the `ieqn`th equation. The code attempts to make \( cc(i) = 0 \) for all
- !! \( i \in \{1,\cdots,m\} \) equations. All relevant consistency equations should be active in
- !! order to make a self-consistent machine.
- !!
- !! **References**
- !!
- !! 1.
- use numerics, only: icc
-
- implicit none
-
- integer, intent(in) :: m
- !! Number of constraint equations to solve
-
- integer, intent(in) :: ieqn
- !! Switch for constraint equations to evaluate;
-
- real(dp), dimension(m), intent(out) :: cc
- !! Residual error in equation i
-
- real(dp), optional, dimension(m), intent(out) :: con
- !! constraint value for equation i in physical units
-
- real(dp), optional, dimension(m), intent(out) :: err
- !! residual error in equation i in physical units
-
- character(len=1), optional, dimension(m), intent(out) :: symbol
- !! `=<`, `>`, `<` symbol for equation i denoting its type
-
- character*10, optional, dimension(m), intent(out) :: units
- !! constraint units in equation i
-
- ! Local variables
- integer :: i, i1, i2
-
- real(dp) :: tmp_cc = 0
- !! Residual error in constraint equation
- real(dp) :: tmp_con = 0
- !! constraint value for constraint equation in physical units
- real(dp) :: tmp_err = 0
- !! residual error in constraint equation in physical units
- character(len=1) :: tmp_symbol = ' '
- !! `=<`, `>`, `<` symbol for constraint equation denoting its type
- character(len=10) :: tmp_units = ' '
- !! constraint units in constraint equation
-
-
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
- ! If ieqn is positive, only evaluate the 'ieqn'th constraint residue,
- ! otherwise evaluate all m constraint residues
- if (ieqn > 0) then
- i1 = ieqn ; i2 = ieqn
- else
- i1 = 1 ; i2 = m
- end if
-
- ! Consistency (equality) constraints should converge to zero.
- do i = i1,i2
-
- ! The constraint value in physical units is
- ! a) for consistency equations, the quantity to be equated, or
- ! b) for limit equations, the limiting value.
- ! The symbol is = for a consistency equation, < for an upper limit
- ! or > for a lower limit.
- select case (icc(i))
-
- ! Relationship between beta, temperature (keV) and density
- case (1); call constraint_eqn_001(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Global plasma power balance equation
- case (2); call constraint_eqn_002(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Global power balance equation for ions
- case (3); call constraint_eqn_003(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Global power balance equation for electrons
- case (4); call constraint_eqn_004(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for density upper limit
- case (5); call constraint_eqn_005(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for epsilon beta-poloidal upper limit
- case (6); call constraint_eqn_006(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for hot beam ion density
- case (7); call constraint_eqn_007(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for neutron wall load upper limit
- case (8); call constraint_eqn_008(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for fusion power upper limit
- case (9); call constraint_eqn_009(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Obsolete
- case (10); call constraint_eqn_010(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for radial build
- case (11); call constraint_eqn_011(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for volt-second capability lower limit
- case (12); call constraint_eqn_012(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for burn time lower limit
- case (13); call constraint_eqn_013(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation to fix number of NBI decay lengths to plasma centre
- case (14); call constraint_eqn_014(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for L-H power threshold limit
- case (15); call constraint_eqn_015(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for net electric power lower limit
- case (16); call constraint_eqn_016(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for radiation power upper limit
- case (17); call constraint_eqn_017(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for divertor heat load upper limit
- case (18); call constraint_eqn_018(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for MVA upper limit
- case (19); call constraint_eqn_019(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for neutral beam tangency radius upper limit
- case (20); call constraint_eqn_020(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for minor radius lower limit
- case (21); call constraint_eqn_021(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Obsolete
- case (22); call constraint_eqn_022(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for conducting shell radius / rminor upper limit
- case (23); call constraint_eqn_023(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for beta upper limit
- case (24); call constraint_eqn_024(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for peak toroidal field upper limit
- case (25); call constraint_eqn_025(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for Central Solenoid current density upper limit at EOF
- case (26); call constraint_eqn_026(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for Central Solenoid current density upper limit at BOP
- case (27); call constraint_eqn_027(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for fusion gain (big Q) lower limit
- case (28); call constraint_eqn_028(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for inboard major radius
- case (29); call constraint_eqn_029(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for injection power upper limit
- case (30); call constraint_eqn_030(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil case stress upper limit (SCTF)
- case (31); call constraint_eqn_031(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil conduit stress upper limit (SCTF)
- case (32); call constraint_eqn_032(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil operating/critical J upper limit (SCTF)
- case (33); call constraint_eqn_033(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil dump voltage upper limit (SCTF)
- case (34); call constraint_eqn_034(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil J_wp/J_prot upper limit (SCTF)
- case (35); call constraint_eqn_035(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil s/c temperature margin lower limit (SCTF)
- case (36); call constraint_eqn_036(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for current drive gamma upper limit
- case (37); call constraint_eqn_037(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for first wall temperature upper limit
- case (39); call constraint_eqn_039(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for auxiliary power lower limit
- case (40); call constraint_eqn_040(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for plasma current ramp-up time lower limit
- case (41); call constraint_eqn_041(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for cycle time lower limit
- case (42); call constraint_eqn_042(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for average centrepost temperature
- case (43); call constraint_eqn_043(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for centrepost temperature upper limit (TART)
- case (44); call constraint_eqn_044(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for edge safety factor lower limit (TART)
- case (45); call constraint_eqn_045(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for Ip/Irod upper limit (TART)
- case (46); call constraint_eqn_046(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for TF coil toroidal thickness upper limit
- case (47); call constraint_eqn_047(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for poloidal beta upper limit
- case (48); call constraint_eqn_048(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Issue #508 Remove RFP option Equation to ensure reversal parameter F is negative
- case (49); call constraint_eqn_049(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! IFE option: Equation for repetition rate upper limit
- case (50); call constraint_eqn_050(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation to enforce startup flux = available startup flux
- case (51); call constraint_eqn_051(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for tritium breeding ratio lower limit
- case (52); call constraint_eqn_052(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for fast neutron fluence on TF coil upper limit
- case (53); call constraint_eqn_053(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for peak TF coil nuclear heating upper limit
- case (54); call constraint_eqn_054(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for helium concentration in vacuum vessel upper limit
- case (55); call constraint_eqn_055(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for power through separatrix / major radius upper limit
- case (56); call constraint_eqn_056(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Obsolete
- case (57); call constraint_eqn_057(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Obsolete
- case (58); call constraint_eqn_058(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for neutral beam shine-through fraction upper limit
- case (59); call constraint_eqn_059(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for Central Solenoid s/c temperature margin lower limit
- case (60); call constraint_eqn_060(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for availability limit
- case (61); call constraint_eqn_061(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Lower limit on f_alpha_energy_confinement the ratio of alpha particle to energy confinement times
- case (62); call constraint_eqn_062(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Upper limit on niterpump (vacuum_model = simple)
- case (63); call constraint_eqn_063(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Upper limit on Zeff
- case (64); call constraint_eqn_064(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Limit TF dump time to calculated quench time
- case (65); call constraint_eqn_065(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Limit on rate of change of energy in poloidal field
- case (66); call constraint_eqn_066(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Simple upper limit on radiation wall load
- case (67); call constraint_eqn_067(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! New Psep scaling (PsepB/qAR)
- case (68); call constraint_eqn_068(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Ensure separatrix power is less than value from Kallenbach divertor
- case (69); call constraint_eqn_069(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Separatrix temperature consistency
- case (70); call constraint_eqn_070(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Separatrix density consistency
- case (71); call constraint_eqn_071(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Central Solenoid Tresca yield criterion
- case (72); call constraint_eqn_072(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! ensure separatrix power is greater than the L-H power + auxiliary power
- case (73); call constraint_eqn_073(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! ensure TF coil quench temperature < tmax_croco
- case (74); call constraint_eqn_074(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! ensure that TF coil current / copper area < Maximum value ONLY used for croco HTS coil
- case (75); call constraint_eqn_075(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Eich critical separatrix density model
- case (76); call constraint_eqn_076(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for maximum TF current per turn upper limit
- case (77); call constraint_eqn_077(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for Reinke criterion, divertor impurity fraction lower limit
- case (78); call constraint_eqn_078(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Equation for maximum CS field
- case (79); call constraint_eqn_079(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Lower limit p_plasma_separatrix_mw
- case (80); call constraint_eqn_080(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint equation making sure that ne(0) > ne(ped)
- case (81); call constraint_eqn_081(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint equation making sure that stellarator coils dont touch in toroidal direction
- case (82); call constraint_eqn_082(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint ensuring radial build consistency for stellarators
- case (83); call constraint_eqn_083(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for lower limit of beta
- case (84); call constraint_eqn_084(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for CP lifetime
- case (85); call constraint_eqn_085(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for turn dimension
- case (86); call constraint_eqn_086(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for cryogenic power
- case (87); call constraint_eqn_087(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for TF coil strain
- case (88); call constraint_eqn_088(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! ensure that OH coil current / copper area < Maximum value ONLY used for croco HTS coil
- case (89); call constraint_eqn_089(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for minimum CS stress load cycles
- case (90); call constraint_eqn_090(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for indication of ECRH ignitability
- case (91); call constraint_eqn_091(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- ! Constraint for D/T ratio
- case (92); call constraint_eqn_092(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- case default
-
- idiags(1) = icc(i)
- call report_error(13)
- tmp_cc = 0
- tmp_con = 0
- tmp_err = 0
- tmp_symbol = ' '
- tmp_units = ' '
-
- end select
-
- cc(i) = tmp_cc
- if (present(con)) then
- con(i) = tmp_con
- if (present(err)) err(i) = tmp_err
- if (present(symbol)) symbol(i) = tmp_symbol
- if (present(units)) units(i) = tmp_units
- end if
-
- if (isnan(cc(i)).or.(abs(cc(i))>9.99D99)) then
-
- ! Add debugging lines as appropriate...
- select case (icc(i))
-
- ! Relationship between beta, temperature (keV) and density (consistency equation)
- case (1); call constraint_err_001()
- ! Equation for net electric power lower limit
- case (16); call constraint_err_016()
- ! Equation for injection power upper limit
- case (30); call constraint_err_030()
- ! Limit on rate of change of energy in poloidal field
- case (66); call constraint_err_066()
-
- end select
-
- idiags(1) = icc(i) ; fdiags(1) = cc(i)
- call report_error(14)
-
- end if
-
- end do
- ! Issue 505 Reverse the sign so it works as an inequality constraint (cc(i) > 0)
- ! This will have no effect if it is used as an equality constraint because it will be squared.
- cc = -cc
-
- end subroutine constraint_eqns
-
- !--- Error-handling routines
-
- subroutine constraint_err_001()
- !! Error in: Relationship between beta, temperature (keV) and density (consistency equation)
- !! author: P B Lloyd, CCFE, Culham Science Centre
- use physics_variables, only: beta_fast_alpha, beta_beam, dene, ten, nd_ions_total, tin, btot, beta
- write(*,*) 'beta_fast_alpha = ', beta_fast_alpha
- write(*,*) 'beta_beam = ', beta_beam
- write(*,*) 'dene = ', dene
- write(*,*) 'ten = ', ten
- write(*,*) 'nd_ions_total = ', nd_ions_total
- write(*,*) 'tin = ', tin
- write(*,*) 'btot = ',btot
- write(*,*) 'beta = ', beta
- end subroutine
-
- subroutine constraint_err_016()
- !! Error in: Equation for net electric power lower limit
- !! author: P B Lloyd, CCFE, Culham Science Centre
- use constraint_variables, only: fpnetel, pnetelin
- use heat_transport_variables, only: pnetelmw
- implicit none
- write(*,*) 'fpnetel = ', fpnetel
- write(*,*) 'pnetelmw = ', pnetelmw
- write(*,*) 'pnetelin = ', pnetelin
- end subroutine
-
- subroutine constraint_err_030()
- !! Error in: Equation for injection power upper limit
- !! author: P B Lloyd, CCFE, Culham Science Centre
- use current_drive_variables, only: p_hcd_injected_total_mw, p_hcd_injected_max
- use constraint_variables, only: fpinj
- implicit none
- write(*,*) 'fpinj = ', fpinj
- write(*,*) 'p_hcd_injected_max = ', p_hcd_injected_max
- write(*,*) 'p_hcd_injected_total_mw = ', p_hcd_injected_total_mw
- end subroutine
-
- subroutine constraint_err_066()
- !! Error in: Limit on rate of change of energy in poloidal field
- !! author: P B Lloyd, CCFE, Culham Science Centre
- use constraint_variables, only: fpoloidalpower
- use pf_power_variables, only: maxpoloidalpower, peakpoloidalpower
- implicit none
- write(*,*) 'fpoloidalpower = ', fpoloidalpower
- write(*,*) 'maxpoloidalpower = ', maxpoloidalpower
- write(*,*) 'peakpoloidalpower = ', peakpoloidalpower
- end subroutine constraint_err_066
-
- !---
-
- subroutine constraint_eqn_001(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- !! author: J Morris
- !! category: equality constraint
- !!
- !! Relationship between beta, temperature (keV) and density
- !!
- !! \begin{equation}
- !! c_i = 1 - \frac{1}{\beta}\left( \beta_{ft} + \beta_{NBI} + 2 \times 10^3 \mu_0 e
- !! \left( \frac{n_e T_e + n_i T_i}{B_{tot}^2} \right) \right)
- !! \end{equation}
- !!
- !! - \( \beta \) -- total plasma beta
- !! - \( \beta_{ft} \) -- fast alpha beta component
- !! - \( \beta_{NBI} \) -- neutral beam beta component
- !! - \( n_e \) -- electron density [m\(^3\)]
- !! - \( n_i \) -- total ion density [m\(^3\)]
- !! - \( T_e \) -- density weighted average electron temperature [keV]
- !! - \( T_i \) -- density weighted average ion temperature [keV]
- !! - \( B_{tot} \) -- total toroidal + poloidal field [T]
-
- use physics_variables, only: beta_fast_alpha, beta_beam, dene, ten, nd_ions_total, tin, btot, beta
- use constants, only: electron_charge,rmu0
-
- implicit none
-
- ! type(constraint_args_type), intent(out) :: args
- real(dp), intent(out) :: tmp_cc
- real(dp), intent(out) :: tmp_con
- real(dp), intent(out) :: tmp_err
- character(len=1), intent(out) :: tmp_symbol
- character(len=10), intent(out) :: tmp_units
-
- !! constraint derived type
-
- tmp_cc = 1.0D0 - (beta_fast_alpha + beta_beam + &
- 2.0D3*rmu0*electron_charge * (dene*ten + nd_ions_total*tin)/btot**2 )/beta
- tmp_con = beta * (1.0D0 - tmp_cc)
- tmp_err = beta * tmp_cc
- tmp_symbol = '='
- tmp_units = ''
-
- end subroutine constraint_eqn_001
-
- subroutine constraint_eqn_002(tmp_cc, tmp_con, tmp_err, tmp_symbol, tmp_units)
- !! author: J. Morris
- !! category: equality constraint
- !!
- !! Global plasma power balance equation
- !!
- !! \begin{equation} c_i =
- !! \end{equation}
- !!
- !! i_rad_loss : input integer : switch for radiation loss term usage in power balance (see User Guide):device.dat):All possible informational/error messages are initialised
- !! via a call to
- !! initialise_error_list. Thereafter, any routine
- !! that needs to flag a message should call
- !! report_error with the relevant error identifier as
- !! the argument. Up to eight integer and eight floating-point diagnostic
- !! values may be saved by the user in arrays idiags and
- !! fdiags, respectively, for debugging purposes.
- !!
The list of messages reported during the course of a run
- !! may be displayed by calling routine
- !! show_errors.
- !!
The error_status variable returns the highest severity
- !! level that has been encountered; if a severe error is flagged
- !! (level 3) the program is terminated immediately.
- !! None
- !
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-#ifndef dp
- use, intrinsic :: iso_fortran_env, only: dp=>real64
-#endif
- implicit none
-
- private
- public :: errors_on, error_status, idiags, fdiags
- public :: initialise_error_list, report_error, show_errors, &
- init_error_handling
-
- ! Switch to turn error handling on
- ! Error reporting is turned off, until either a severe error is found, or
- ! during an output step. Warnings during intermediate iteration steps
- ! may be premature and might clear themselves at the final converged point.
-
- logical :: errors_on
-
- ! Levels of severity
-
- integer, parameter :: ERROR_OKAY = 0
- integer, parameter :: ERROR_INFO = 1
- integer, parameter :: ERROR_WARN = 2
- integer, parameter :: ERROR_SEVERE = 3
-
- integer :: error_id
- !! error_id : identifier for final message encountered
-
- ! Overall status
-
- integer :: error_status
- !! error_status : overall status flag for a run; on exit:
The error messages are read in from a JSON-format file. - !! None - ! - ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - - use fson_library, only: fson_parse, fson_value, fson_get, fson_destroy - implicit none - - ! Arguments - - ! Local variables - - integer :: n_errortypes - character(len=180) :: filename - type(fson_value), pointer :: errorfile - - ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - - ! Parse the json file - character(len=200) :: process_dir - CALL get_environment_variable("PYTHON_PROCESS_ROOT", process_dir) - if (process_dir == "") then - filename = INSTALLDIR//'/process//utilities/errorlist.json' - else - filename = trim(process_dir)//'/utilities/errorlist.json' - end if - errorfile => fson_parse(trim(filename)) - - ! Allocate memory for error_type array contents - - call fson_get(errorfile, "n_errortypes", n_errortypes) - - ! Guard against re-allocation - if (allocated(error_type)) deallocate(error_type) - allocate(error_type(n_errortypes)) - - error_type(:)%level = ERROR_OKAY - error_type(:)%message = ' ' - - ! Extract information arrays from the file - - call fson_get(errorfile, "errors", "level", error_type%level) - call fson_get(errorfile, "errors", "message", error_type%message) - - ! Clean up - - call fson_destroy(errorfile) - - end subroutine initialise_error_list - - ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - - subroutine report_error(error_id) - - !! Adds the latest error message to the list already specified - !! author: P J Knight, CCFE, Culham Science Centre - !! error_id : input integer : identifier (error_type element number) - !! for the relevant error - !! This routine should be called if a informational, warning - !! or error message needs to be flagged. - !! It uses a linked list (see references) to provide - !! an audit trail for any such messages during the program - !! execution. - !!
Up to eight integer and eight floating-point diagnostic
- !! values may be saved by the user in arrays idiags and
- !! fdiags, respectively, for debugging; these arrays must
- !! be assigned with up to eight values each prior to calling this routine.
- !!
The error_status variable returns the highest severity
- !! level that has been encountered; if a severe error is flagged
- !! (level 3) the program is terminated immediately.
- !! Introduction to Fortran 90/95, Stephen J, Chapman, pp.467-472,
- !! McGraw-Hill, ISBN 0-07-115896-0
- !
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
- implicit none
-
- ! Arguments
-
- integer, intent(in) :: error_id
-
- ! Local variables
-
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
- ! Turn on error handling if a severe error has been encountered
-
- if (error_type(error_id)%level == ERROR_SEVERE) errors_on = .true.
-
- ! Error handling is only turned on during an output step, not during
- ! intermediate iteration steps
-
- if (.not.errors_on) then
- idiags = INT_DEFAULT
- fdiags = FLT_DEFAULT
- return
- end if
-
- if (.not.associated(error_head)) then
- allocate(error_head)
- error_tail => error_head
- else
-
-! SJP Issue #867
-! Remove consecutive identical error messages
-
- if (error_tail%id == error_id) return
-
- allocate(error_tail%ptr)
- error_tail => error_tail%ptr
- end if
-
- error_tail%id = error_id
- error_tail%data%level = error_type(error_id)%level
- error_tail%data%message = error_type(error_id)%message
- error_tail%data%idiags = idiags ; idiags = INT_DEFAULT
- error_tail%data%fdiags = fdiags ; fdiags = FLT_DEFAULT
-
- nullify (error_tail%ptr)
-
- ! Update the overall error status (highest severity level encountered)
- ! and stop the program if a severe error has occurred
-
- error_status = max(error_status, error_type(error_id)%level)
-
- if (error_status == ERROR_SEVERE) then
- call show_errors
- write(*,*) 'PROCESS stopping.'
- stop 1
- end if
-
- end subroutine report_error
-
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
- subroutine show_errors
-
- !! Reports all informational/error messages encountered
- !! author: P J Knight, CCFE, Culham Science Centre
- !! None
- !! This routine provides a summary audit trail of all the errors
- !! encountered during the program's execution.
- !! Introduction to Fortran 90/95, Stephen J, Chapman, pp.467-472,
- !! McGraw-Hill, ISBN 0-07-115896-0
- !
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
- use constants, only: iotty, nout
- use process_output_fortran, only: oblnkl, oheadr, ocmmnt, ovarin, ostars
- implicit none
-
- ! Arguments
-
- ! Local variables
-
- type (error_list_item), pointer :: ptr
- integer :: i
- character(len=50) :: status_message
-
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
- call oheadr(iotty,'Errors and Warnings')
- call oheadr(nout,'Errors and Warnings')
- call ocmmnt(nout,'(See top of file for solver errors and warnings.)')
-
- select case (error_status)
- case (0)
- status_message = 'No messages'
- case (1)
- status_message = 'Information messages only'
- case (2)
- status_message = 'Warning messages'
- case (3)
- status_message = 'Errors'
- case default
- status_message = 'Incorrect value of error_status'
- end select
-
- call ocmmnt(nout,'PROCESS status flag: '//status_message)
- write(*,*) 'PROCESS status flag: '//status_message
- call oblnkl(iotty)
- call ovarin(nout,'PROCESS error status flag','(error_status)',error_status)
-
- ptr => error_head
-
- if (.not.associated(ptr)) then
- call ovarin(nout,'Final error/warning identifier','(error_id)',error_id)
- return
- end if
-
- write(*,*) 'ID LEVEL MESSAGE'
-
- output: do
- if (.not.associated(ptr)) exit output
-
- error_id = ptr%id
- write(nout,'(i3,t7,i3,t13,a80)') ptr%id,ptr%data%level,ptr%data%message
- write(*, '(i3,t7,i3,t13,a80)') ptr%id,ptr%data%level,ptr%data%message
-
- if (any(ptr%data%idiags /= INT_DEFAULT)) then
- write(*,*) 'Integer diagnostic values for this error:'
- do i = 1,8
- if (ptr%data%idiags(i) /= INT_DEFAULT) then
- write(*,'(i4,a,i14)') i,') ',ptr%data%idiags(i)
- write(nout,'(i4,a,i14)') i,') ',ptr%data%idiags(i)
- endif
- end do
- end if
- if (any(ptr%data%fdiags /= FLT_DEFAULT)) then
- write(*,*) 'Floating point diagnostic values for this error:'
- do i = 1,8
- if (ptr%data%fdiags(i) /= FLT_DEFAULT) then
- write(*,'(i4,a,1pe14.5)') i,') ',ptr%data%fdiags(i)
- write(nout,'(i4,a,1pe14.5)') i,') ',ptr%data%fdiags(i)
- endif
- end do
- end if
- write(*,*) ' '
-
- ptr => ptr%ptr
- end do output
- call ostars(iotty, 110)
- call oblnkl(iotty)
-
- call ovarin(nout,'Final error identifier','(error_id)',error_id)
-
- end subroutine show_errors
-
-end module error_handling
diff --git a/source/fortran/fson_library.f90 b/source/fortran/fson_library.f90
deleted file mode 100644
index f1010b224a..0000000000
--- a/source/fortran/fson_library.f90
+++ /dev/null
@@ -1,1794 +0,0 @@
-! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-!+PJK Need to convert to autodoc style + usual PROCESS standard layout
-
-!-----------------------------------------------------------------------------------------
-!
-! FSON library
-!
-! Extracted from https://github.com/josephalevin/fson on 9th July 2014
-!
-! with selected updates taken from the version forked as
-!
-! https://github.com/jmozmoz/fson/commit/b210a9011bb804957546e2a4b6eade578e7035ef
-!
-! plus some improvements to help with array handling and double precision, by
-! P J Knight, 17th July 2014
-!
-! Comprises the following original FSON files:
-! string.f95, value_m.f95, fson_path_m.f95, fson.f95
-!
-!-----------------------------------------------------------------------------------------
-!
-! Copyright (c) 2012 Joseph A. Levin
-!
-! Permission is hereby granted, free of charge, to any person obtaining a copy of this
-! software and associated documentation files (the "Software"), to deal in the Software
-! without restriction, including without limitation the rights to use, copy, modify, merge,
-! publish, distribute, sublicense, and/or sell copies of the Software, and to permit
-! persons to whom the Software is furnished to do so, subject to the following conditions:
-!
-! The above copyright notice and this permission notice shall be included in all copies or
-! substantial portions of the Software.
-!
-! THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED,
-! INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
-! PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
-! LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT
-! OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
-! DEALINGS IN THE SOFTWARE.
-!
-!-----------------------------------------------------------------------------------------
-
-module fson_string_m
-
- private
-
- public :: fson_string, fson_string_create, fson_string_destroy, fson_string_length, &
- fson_string_append, fson_string_clear
- public :: fson_string_equals, fson_string_copy
-
- integer, parameter :: BLOCK_SIZE = 32
-
- type fson_string
- character (len = BLOCK_SIZE) :: chars
- integer :: index = 0
- type(fson_string), pointer :: next => null()
- end type fson_string
-
- interface fson_string_append
- module procedure append_chars, append_string
- end interface
-
- interface fson_string_copy
- module procedure copy_chars
- end interface
-
- interface fson_string_equals
- module procedure equals_string
- end interface
-
- interface fson_string_length
- module procedure string_length
- end interface
-
-contains
-
- !
- ! FSON STRING CREATE
- !
- function fson_string_create(chars) result(new)
- character(len=*), optional :: chars
- type(fson_string), pointer :: new
-
- allocate(new)
-
- ! append chars if available
- if (present(chars)) then
- call append_chars(new, chars)
- end if
-
- end function fson_string_create
-
- !
- ! FSON STRING DESTROY
- !
- recursive subroutine fson_string_destroy(this)
- type(fson_string), pointer :: this
-
- if (associated(this % next)) then
- call fson_string_destroy(this % next)
- end if
-
- nullify (this % next)
- nullify (this)
- end subroutine fson_string_destroy
-
- !
- ! ALLOCATE BLOCK
- !
- subroutine allocate_block(this)
- type(fson_string), pointer :: this
- type(fson_string), pointer :: new
-
- if (.not.associated(this % next)) then
- allocate(new)
- this % next => new
- end if
-
- end subroutine allocate_block
-
- !
- ! APPEND_STRING
- !
- subroutine append_string(str1, str2)
- type(fson_string), pointer :: str1, str2
- integer length, i
-
- length = string_length(str2)
-
- do i = 1, length
- call append_char(str1, get_char_at(str2, i))
- end do
-
- end subroutine append_string
-
- !
- ! APPEND_CHARS
- !
- subroutine append_chars(str, c)
- type(fson_string), pointer :: str
- character (len = *), intent(in) :: c
- integer length, i
-
- length = len(c)
-
- do i = 1, length
- call append_char(str, c(i:i))
- end do
-
- end subroutine append_chars
-
- !
- ! APPEND_CHAR
- !
- recursive subroutine append_char(str, c)
- type(fson_string), pointer :: str
- character, intent(in) :: c
-
- if (str % index >= BLOCK_SIZE) then
- !set down the chain
- call allocate_block(str)
- call append_char(str % next, c)
-
- else
- ! set local
- str % index = str % index + 1
- str % chars(str % index:str % index) = c
- end if
-
- end subroutine append_char
-
- !
- ! COPY CHARS
- !
- subroutine copy_chars(this, to)
- type(fson_string), pointer :: this
- character(len = *), intent(inout) :: to
- integer :: length
-
- length = min(string_length(this), len(to))
-
- do i = 1, length
- to(i:i) = get_char_at(this, i)
- end do
-
- ! pad with nothing
- do i = length + 1, len(to)
- to(i:i) = ""
- end do
-
- end subroutine copy_chars
-
- !
- ! CLEAR
- !
- recursive subroutine string_clear(this)
- type(fson_string), pointer :: this
-
- if (associated(this % next)) then
- call string_clear(this % next)
- deallocate(this % next)
- nullify (this % next)
- end if
-
- this % index = 0
-
- end subroutine string_clear
-
- !
- ! SIZE
- !
- recursive integer function string_length(str) result(count)
- type(fson_string), pointer :: str
-
- count = str % index
-
- if (str % index == BLOCK_SIZE .AND. associated(str % next)) then
- count = count + string_length(str % next)
- end if
-
- end function string_length
-
- !
- ! GET CHAR AT
- !
- recursive character function get_char_at(this, i) result(c)
- type(fson_string), pointer :: this
- integer, intent(in) :: i
-
- if (i <= this % index) then
- c = this % chars(i:i)
- else
- c = get_char_at(this % next, i - this % index)
- end if
-
- end function get_char_at
-
- !
- ! EQUALS STRING
- !
- logical function equals_string(this, other) result(equals)
- type(fson_string), pointer :: this, other
- integer :: i
- equals = .false.
-
- if (fson_string_length(this) /= fson_string_length(other)) then
- equals = .false.
- return
- else if (fson_string_length(this) == 0) then
- equals = .true.
- return
- end if
-
- do i=1, fson_string_length(this)
- if (get_char_at(this, i) /= get_char_at(other, i)) then
- equals = .false.
- return
- end if
- end do
-
- equals = .true.
-
- end function equals_string
-
-end module fson_string_m
-
-! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-module fson_value_m
- use fson_string_m, only: fson_string
- implicit none
-
- private
-
- public :: fson_value, fson_value_create, fson_value_destroy, fson_value_add, fson_value_get, &
- fson_value_count, fson_value_print
-
- ! constants for the value types
- integer, public, parameter :: TYPE_UNKNOWN = -1
- integer, public, parameter :: TYPE_NULL = 0
- integer, public, parameter :: TYPE_OBJECT = 1
- integer, public, parameter :: TYPE_ARRAY = 2
- integer, public, parameter :: TYPE_STRING = 3
- integer, public, parameter :: TYPE_INTEGER = 4
- integer, public, parameter :: TYPE_REAL = 5
- integer, public, parameter :: TYPE_LOGICAL = 6
-
- !
- ! FSON VALUE
- !
- type fson_value
- type(fson_string), pointer :: name => null()
- integer :: value_type = TYPE_UNKNOWN
- logical :: value_logical
- integer :: value_integer
- real :: value_real
- !+PJK
- real(kind(1.0D0)) :: value_double
- !-PJK
- type(fson_string), pointer :: value_string => null()
- type(fson_value), pointer :: next => null()
- type(fson_value), pointer :: parent => null()
- type(fson_value), pointer :: children => null()
- end type fson_value
-
- !
- ! FSON VALUE GET
- !
- ! Use either a 1 based index or member name to get the value.
- interface fson_value_get
- module procedure get_by_index
- module procedure get_by_name_chars
- module procedure get_by_name_string
- end interface
-
-contains
-
- !
- ! FSON VALUE CREATE
- !
- function fson_value_create() result(new)
- type(fson_value), pointer :: new
-
- allocate(new)
-
- end function fson_value_create
-
- !
- ! FSON VALUE DESTROY
- !
- recursive subroutine fson_value_destroy(this)
- use fson_string_m, only: fson_string_destroy
- implicit none
-
- type(fson_value), pointer :: this
-
- if (associated(this % children)) then
- call fson_value_destroy(this % children)
- nullify(this % children)
- end if
-
- if (associated(this % next)) then
- call fson_value_destroy(this % next)
- nullify (this % next)
- end if
-
- if (associated(this % name)) then
- call fson_string_destroy(this % name)
- nullify (this % name)
- end if
-
- if (associated(this % value_string)) then
- call fson_string_destroy(this % value_string)
- nullify (this % value_string)
- end if
-
- nullify(this)
-
- end subroutine fson_value_destroy
-
- !
- ! FSON VALUE ADD
- !
- ! Adds the member to the linked list
- subroutine fson_value_add(this, member)
- type(fson_value), pointer :: this, member, p
-
- ! associate the parent
- member % parent => this
-
- ! add to linked list
- if (associated(this % children)) then
- ! get to the tail of the linked list
- p => this % children
- do while (associated(p % next))
- p => p % next
- end do
-
- p % next => member
- else
- this % children => member
- end if
-
- end subroutine fson_value_add
-
- !
- ! FSON_VALUE_COUNT
- !
- integer function fson_value_count(this) result(count)
- type(fson_value), pointer :: this, p
-
- count = 0
-
- p => this % children
-
- do while (associated(p))
- count = count + 1
- p => p % next
- end do
-
- end function fson_value_count
-
- !
- ! GET BY INDEX
- !
- function get_by_index(this, index) result(p)
- type(fson_value), pointer :: this, p
- integer, intent(in) :: index
- integer :: i
-
- p => this % children
-
- do i = 1, index - 1
- p => p % next
- end do
-
- end function get_by_index
-
- !
- ! GET BY NAME CHARS
- !
- function get_by_name_chars(this, name) result(p)
- use fson_string_m, only: fson_string, fson_string_create
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), intent(in) :: name
-
- type(fson_string), pointer :: string
-
- ! convert the char array into a string
- string => fson_string_create(name)
-
- p => get_by_name_string(this, string)
-
- end function get_by_name_chars
-
- !
- ! GET BY NAME STRING
- !
- function get_by_name_string(this, name) result(p)
- use fson_string_m, only: fson_string, fson_string_equals
- implicit none
-
- type(fson_value), pointer :: this, p
- type(fson_string), pointer :: name
- integer :: i
-
- if (this % value_type /= TYPE_OBJECT) then
- nullify(p)
- return
- end if
-
- do i=1, fson_value_count(this)
- p => fson_value_get(this, i)
- if (fson_string_equals(p%name, name)) then
- return
- end if
- end do
-
- ! didn't find anything
- nullify(p)
-
- end function get_by_name_string
-
- !
- ! FSON VALUE PRINT
- !
- recursive subroutine fson_value_print(this, indent)
- use fson_string_m, only: fson_string_copy
- implicit none
-
- type(fson_value), pointer :: this, element
- integer, optional, intent(in) :: indent
- character(len=1024) :: tmp_chars
- integer :: tab, i, count, spaces
-
- !+PJK
- if (.not.associated(this)) return
- !-PJK
-
- if (present(indent)) then
- tab = indent
- else
- tab = 0
- end if
-
- spaces = tab * 2
-
- select case (this % value_type)
- case(TYPE_OBJECT)
- print *, repeat(" ", spaces), "{"
- count = fson_value_count(this)
- do i = 1, count
- ! get the element
- element => fson_value_get(this, i)
- ! get the name
- call fson_string_copy(element % name, tmp_chars)
- ! print the name
- print *, repeat(" ", spaces), '"', trim(tmp_chars), '":'
- ! recursive print of the element
- call fson_value_print(element, tab + 1)
- ! print the separator if required
- if (i < count) then
- print *, repeat(" ", spaces), ","
- end if
- end do
-
- print *, repeat(" ", spaces), "}"
- case (TYPE_ARRAY)
- print *, repeat(" ", spaces), "["
- count = fson_value_count(this)
- do i = 1, count
- ! get the element
- element => fson_value_get(this, i)
- ! recursive print of the element
- call fson_value_print(element, tab + 1)
- ! print the separator if required
- if (i < count) then
- print *, ","
- end if
- end do
- print *, repeat(" ", spaces), "]"
- case (TYPE_NULL)
- print *, repeat(" ", spaces), "null"
- case (TYPE_STRING)
- call fson_string_copy(this % value_string, tmp_chars)
- print *, repeat(" ", spaces), '"', trim(tmp_chars), '"'
- case (TYPE_LOGICAL)
- if (this % value_logical) then
- print *, repeat(" ", spaces), "true"
- else
- print *, repeat(" ", spaces), "false"
- end if
- case (TYPE_INTEGER)
- print *, repeat(" ", spaces), this % value_integer
- case (TYPE_REAL)
- print *, repeat(" ", spaces), this % value_real ! N.B. doubles will be shown as single precision
- end select
- end subroutine fson_value_print
-
-end module fson_value_m
-
-! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-! Modifications by P J Knight to handle arrays more conveniently than by providing
-! an array_callback argument to get_array, and to handle double precision values better
-
-module fson_path_m
-
- private
-
- public :: fson_path_get, array_callback
-
- interface fson_path_get
- module procedure get_by_path
- module procedure get_integer
- module procedure get_real
- module procedure get_double
- module procedure get_logical
- module procedure get_chars
- module procedure get_array
- !+PJK
- module procedure get_int_array
- module procedure get_real_array
- module procedure get_double_array
- module procedure get_string_array
- module procedure get_int_array_in_struct
- module procedure get_real_array_in_struct
- module procedure get_double_array_in_struct
- module procedure get_string_array_in_struct
- !-PJK
- end interface
-
-contains
-
- !
- ! GET BY PATH
- !
- ! $ = root
- ! @ = this
- ! . = child object member
- ! [] = child array element
- !
- recursive subroutine get_by_path(this, path, p)
- use fson_value_m, only: fson_value, fson_value_get
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), intent(inout) :: path
- integer :: i, length, child_i
- character :: c
- logical :: array
-
- ! default to assuming relative to this
- p => this
-
- child_i = 1
-
- array = .false.
-
- length = len_trim(path)
-
- do i=1, length
- c = path(i:i)
- select case (c)
- case ("$")
- ! root
- do while (associated (p % parent))
- p => p % parent
- end do
- child_i = i + 1
- case ("@")
- ! this
- p => this
- child_i = i + 1
- case (".", "[")
- ! get child member from p
- if (child_i < i) then
- p => fson_value_get(p, path(child_i:i-1))
- else
- child_i = i + 1
- cycle
- end if
-
- if (.not.associated(p)) then
- return
- end if
-
- child_i = i + 1
-
- ! check if this is an array
- ! if so set the array flag
- if (c == "[") then
- ! start looking for the array element index
- array = .true.
- end if
- case ("]")
- if (.not.array) then
- print *, "ERROR: Unexpected ], not missing preceding ["
- call exit(1)
- end if
- array = .false.
- child_i = parse_integer(path(child_i:i-1))
- p => fson_value_get(p, child_i)
-
- child_i = i + 1
- end select
- end do
-
- ! grab the last child if present in the path
- if (child_i <= length) then
- p => fson_value_get(p, path(child_i:i-1))
- if (.not.associated(p)) then
- return
- else
- end if
- end if
-
- end subroutine get_by_path
-
- !
- ! PARSE INTEGER
- !
- integer function parse_integer(chars) result(integral)
- character(len=*) :: chars
- character :: c
- integer :: tmp, i
-
- integral = 0
- do i=1, len_trim(chars)
- c = chars(i:i)
- select case(c)
- case ("0":"9")
- ! digit
- read (c, '(i1)') tmp
- ! shift
- if (i > 1) then
- integral = integral * 10
- end if
- ! add
- integral = integral + tmp
-
- case default
- return
- end select
- end do
-
- end function parse_integer
-
- !
- ! GET INTEGER
- !
- subroutine get_integer(this, path, value)
- use fson_value_m, only: type_integer, type_real, type_logical, fson_value
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), optional :: path
- integer :: value
-
- nullify(p)
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_INTEGER) then
- value = p % value_integer
- else if (p % value_type == TYPE_REAL) then
- value = p % value_real
- else if (p % value_type == TYPE_LOGICAL) then
- if (p % value_logical) then
- value = 1
- else
- value = 0
- end if
- else
- print *, "Unable to resolve value to integer: ", path
- call exit(1)
- end if
-
- end subroutine get_integer
-
- !
- ! GET REAL
- !
- subroutine get_real(this, path, value)
- use fson_value_m, only: type_integer, type_real, type_logical, fson_value
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), optional :: path
- real :: value
-
- nullify(p)
-
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_INTEGER) then
- value = p % value_integer
- else if (p % value_type == TYPE_REAL) then
- value = p % value_real
- else if (p % value_type == TYPE_LOGICAL) then
- if (p % value_logical) then
- value = 1
- else
- value = 0
- end if
- else
- print *, "Unable to resolve value to real: ", path
- call exit(1)
- end if
-
- end subroutine get_real
-
- !
- ! GET DOUBLE
- !
- subroutine get_double(this, path, value)
- use fson_value_m, only: type_integer, type_real, type_logical, fson_value
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), optional :: path
- real(kind(1.0D0)) :: value
-
- nullify(p)
-
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_INTEGER) then
- value = p % value_integer
- else if (p % value_type == TYPE_REAL) then
- value = p % value_double ! PJK from value_real
- else if (p % value_type == TYPE_LOGICAL) then
- if (p % value_logical) then
- value = 1
- else
- value = 0
- end if
- else
- print *, "Unable to resolve value to double: ", path
- call exit(1)
- end if
-
- end subroutine get_double
-
- !
- ! GET LOGICAL
- !
- subroutine get_logical(this, path, value)
- use fson_value_m, only: type_integer, type_logical, fson_value
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), optional :: path
- logical :: value
-
- nullify(p)
-
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_INTEGER) then
- value = (p % value_integer > 0)
- else if (p % value_type == TYPE_LOGICAL) then
- value = p % value_logical
- else
- print *, "Unable to resolve value to logical: ", path
- call exit(1)
- end if
-
- end subroutine get_logical
-
- !
- ! GET CHARS
- !
- subroutine get_chars(this, path, value)
- use fson_value_m, only: type_string, fson_value
- use fson_string_m, only: fson_string_copy
- implicit none
-
- type(fson_value), pointer :: this, p
- character(len=*), optional :: path
- character(len=*) :: value
-
- nullify(p)
-
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_STRING) then
- call fson_string_copy(p % value_string, value)
- else
- print *, "Unable to resolve value to characters: ", path
- call exit(1)
- end if
-
- end subroutine get_chars
-
- !
- ! GET ARRAY (original version using array_callback)
- !
- subroutine get_array(this, path, array_callback)
- use fson_value_m, only: type_array, fson_value_get, fson_value_count, &
- fson_value
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*), optional :: path
- integer :: index, count
-
- ! ELEMENT CALLBACK (PJK: Added example comments)
- interface
- subroutine array_callback(element, index, count)
- use fson_value_m, only: fson_value
- implicit none
-
- ! In the actual routine add a second 'use' line as follows:
- !use shared_data ! contains declarations for the array(s) to be populated
-
- type(fson_value), pointer :: element
- integer :: index, count
-
- ! Example usage in the actual routine:
- ! call fson_get(element, "element_name", array_name(index))
-
- end subroutine array_callback
- end interface
-
- nullify(p)
-
- ! resolve the path to the value
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- call array_callback(element, index, count)
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_array
-
- !+PJK
- !
- ! GET INT ARRAY
- !
- subroutine get_int_array(this, path, array)
- use fson_value_m, only: type_array, fson_value_get, fson_value_count, &
- fson_value
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*), optional :: path
- integer :: index, count
- integer, dimension(:) :: array
-
- nullify(p)
-
- ! resolve the path to the value
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- array(index) = element%value_integer
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_int_array
-
- !
- ! GET REAL ARRAY
- !
- subroutine get_real_array(this, path, array)
- use fson_value_m, only: type_array, fson_value_get, fson_value_count, &
- fson_value
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*), optional :: path
- integer :: index, count
- real, dimension(:) :: array
-
- nullify(p)
-
- ! resolve the path to the value
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- array(index) = element%value_real
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_real_array
-
- !
- ! GET DOUBLE ARRAY
- !
- subroutine get_double_array(this, path, array)
- use fson_value_m, only: type_array, fson_value_get, fson_value_count, &
- fson_value
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*), optional :: path
- integer :: index, count
- real(kind(1.0D0)), dimension(:) :: array
-
- nullify(p)
-
- ! resolve the path to the value
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- array(index) = element%value_double
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_double_array
-
- !
- ! GET STRING ARRAY
- !
- subroutine get_string_array(this, path, array)
- use fson_value_m, only: fson_value_count, fson_value_get, fson_value, &
- TYPE_ARRAY
- use fson_string_m, only: fson_string_copy
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*), optional :: path
- integer :: index, count
- character(len=*), dimension(:) :: array
-
- nullify(p)
-
- ! resolve the path to the value
- if (present(path)) then
- call get_by_path(this=this, path=path, p=p)
- else
- p => this
- end if
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- call fson_string_copy(element%value_string, array(index))
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_string_array
-
- !
- ! GET INT ARRAY IN STRUCTURE
- !
- subroutine get_int_array_in_struct(this, path, subpath, array)
- use fson_value_m, only: fson_value_count, fson_value_get, fson_value, &
- TYPE_ARRAY
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*) :: path, subpath
- integer, dimension(:), intent(out) :: array
- integer :: index, count
-
- nullify(p)
-
- ! resolve the path to the value
- call get_by_path(this=this, path=path, p=p)
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- call get_integer(element, subpath, array(index))
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_int_array_in_struct
-
- !
- ! GET REAL ARRAY IN STRUCTURE
- !
- subroutine get_real_array_in_struct(this, path, subpath, array)
- use fson_value_m, only: fson_value_count, fson_value_get, fson_value, &
- TYPE_ARRAY
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*) :: path, subpath
- real, dimension(:), intent(out) :: array
- integer :: index, count
-
- nullify(p)
-
- ! resolve the path to the value
- call get_by_path(this=this, path=path, p=p)
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- call get_real(element, subpath, array(index))
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_real_array_in_struct
-
- !
- ! GET DOUBLE ARRAY IN STRUCTURE
- !
- subroutine get_double_array_in_struct(this, path, subpath, array)
- use fson_value_m, only: fson_value_count, fson_value_get, fson_value, &
- TYPE_ARRAY
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*) :: path, subpath
- real(kind(1.0D0)), dimension(:), intent(out) :: array
- integer :: index, count
-
- nullify(p)
-
- ! resolve the path to the value
- call get_by_path(this=this, path=path, p=p)
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- call get_double(element, subpath, array(index))
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_double_array_in_struct
-
- !
- ! GET STRING ARRAY IN STRUCTURE
- !
- subroutine get_string_array_in_struct(this, path, subpath, array)
- use fson_value_m, only: fson_value_count, fson_value_get, fson_value, &
- TYPE_ARRAY
- implicit none
-
- type(fson_value), pointer :: this, p, element
- character(len=*) :: path, subpath
- character(len=*), dimension(:), intent(out) :: array
- integer :: index, count
-
- nullify(p)
-
- ! resolve the path to the value
- call get_by_path(this=this, path=path, p=p)
-
- if (.not.associated(p)) then
- print *, "Unable to resolve path: ", path
- call exit(1)
- end if
-
- if (p % value_type == TYPE_ARRAY) then
- count = fson_value_count(p)
- do index = 1, count
- element => fson_value_get(p, index)
- call get_chars(element, subpath, array(index))
- end do
- else
- print *, "Resolved value is not an array. ", path
- call exit(1)
- end if
-
- end subroutine get_string_array_in_struct
- !-PJK
-
-end module fson_path_m
-
-! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-module fson_library
-
- !! JSON file reading module
- !! author: P J Knight, CCFE, Culham Science Centre
- !! N/A
- !! This module uses a local copy of the freely-available FSON library,
- !! written by Joseph A. Levin and distributed via github,
- !! to enable PROCESS to read in information from JSON-formatted files.
- !! None
- !
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
- use fson_value_m, only: fson_print => fson_value_print, &
- fson_destroy => fson_value_destroy, fson_value
- use fson_path_m, only: fson_get => fson_path_get
-
- implicit none
-
- private
-
- public :: fson_parse, fson_value, fson_get, fson_print, fson_destroy, &
- init_fson_library
-
- ! FILE IOSTAT CODES
- integer, parameter :: end_of_file = -1
- integer, parameter :: end_of_record = -2
-
- ! PARSING STATES
- integer, parameter :: STATE_LOOKING_FOR_VALUE = 1
- integer, parameter :: STATE_IN_OBJECT = 2
- integer, parameter :: STATE_IN_PAIR_NAME = 3
- integer, parameter :: STATE_IN_PAIR_VALUE = 4
-
- ! POP/PUSH CHARACTER
- integer :: pushed_index
- character (len = 10) :: pushed_char
-
-contains
-
- subroutine init_fson_library
- !! Initialise fson library module variables
- implicit none
-
- pushed_index = 0
- end subroutine init_fson_library
-
- ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
- !
- ! FSON PARSE
- !
- function fson_parse(file, unit) result(p)
- use fson_value_m, only: fson_value_create
- implicit none
-
- type(fson_value), pointer :: p
- integer, optional, intent(inout) :: unit
- character(len = *), intent(in) :: file
- logical :: unit_available
- integer :: u
- ! init the pointer to null
- nullify(p)
-
- ! select the file unit to use
- if (present(unit)) then
- u = unit
- else
- ! find the first available unit
- unit_available = .false.
- u = 20
-
- do while (.not.unit_available)
- inquire(unit=u, exist=unit_available)
- u = u + 1
- end do
- end if
-
- ! open the file
- open (unit=u, file=file, status="old", action="read", &
- form="formatted", position="rewind")
-
- ! create the value and associate the pointer
- p => fson_value_create()
-
- ! parse as a value
- call parse_value(unit=u, value=p)
-
- ! close the file
- if ( .not. present(unit)) then
- close (u)
- end if
-
- end function fson_parse
-
- !
- ! PARSE_VALUE
- !
- recursive subroutine parse_value(unit, value)
- use fson_value_m, only: TYPE_ARRAY, TYPE_LOGICAL, TYPE_NULL, TYPE_OBJECT, &
- TYPE_STRING
- implicit none
-
- integer, intent(inout) :: unit
- type(fson_value), pointer :: value
- logical :: eof
- character :: c
-
- ! for some unknown reason the next pointer is getting messed with the pop
- type(fson_value), pointer :: hack
-
- ! start the hack
- hack => value % next
-
- ! pop the next non whitespace character off the file
- c = pop_char(unit, eof=eof, skip_ws=.true.)
-
- ! finish the hack; set the next pointer to whatever it was before the pop
- value % next => hack
-
- if (eof) then
- return
- else
- select case (c)
- case ("{")
- ! start object
- value % value_type = TYPE_OBJECT
- call parse_object(unit, value)
- case ("[")
- ! start array
- value % value_type = TYPE_ARRAY
- call parse_array(unit, value)
- case ("]")
- ! end an empty array
- nullify(value)
- case ('"')
- ! string
- value % value_type = TYPE_STRING
- value % value_string => parse_string(unit)
- case ("t")
- ! true
- value % value_type = TYPE_LOGICAL
- call parse_for_chars(unit, "rue")
- value % value_logical = .true.
- case ("f")
- ! false
- value % value_type = TYPE_LOGICAL
- value % value_logical = .false.
- call parse_for_chars(unit, "alse")
- case ("n")
- value % value_type = TYPE_NULL
- call parse_for_chars(unit, "ull")
- case("-", "0": "9")
- call push_char(c)
- call parse_number(unit, value)
- case default
- print *, "ERROR: Unexpected character while parsing value. '", c, "' ASCII=", iachar(c)
- call exit (1)
- end select
- end if
-
- end subroutine parse_value
-
- !
- ! PARSE OBJECT
- !
- recursive subroutine parse_object(unit, parent)
- use fson_value_m, only: fson_value_create, fson_value_add
- implicit none
-
- integer, intent(inout) :: unit
- type(fson_value), pointer :: parent, pair
-
- logical :: eof
- character :: c
-
- ! pair name
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- print *, "ERROR: Unexpected end of file while parsing start of object."
- call exit (1)
- else if ("}" == c) then
- ! end of an empty object
- return
- else if ('"' == c) then
- pair => fson_value_create()
- pair % name => parse_string(unit)
- else
- print *, "ERROR: Expecting string: '", c, "'"
- call exit (1)
- end if
-
- ! pair value
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- print *, "ERROR: Unexpected end of file while parsing object member. 1"
- call exit (1)
- else if (":" == c) then
- ! parse the value
- call parse_value(unit, pair)
- call fson_value_add(parent, pair)
- else
- print *, "ERROR: Expecting : and then a value. ", c
- call exit (1)
- end if
-
- ! another possible pair
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- return
- else if ("," == c) then
- ! read the next member
- call parse_object(unit=unit, parent=parent)
- else if ("}" == c) then
- return
- else
- print *, "ERROR: Expecting end of object.", c
- call exit (1)
- end if
-
- end subroutine parse_object
-
- !
- ! PARSE ARRAY
- !
- recursive subroutine parse_array(unit, array)
- use fson_value_m, only: fson_value_create, fson_value_add
- implicit none
-
- integer, intent(inout) :: unit
- type(fson_value), pointer :: array, element
-
- logical :: eof
- character :: c
-
- ! try to parse an element value
- element => fson_value_create()
- call parse_value(unit, element)
-
- ! parse value will disassociate an empty array value
- if (associated(element)) then
- call fson_value_add(array, element)
- end if
-
- ! popped the next character
- c = pop_char(unit, eof=eof, skip_ws=.true.)
-
- if (eof) then
- return
- else if ("," == c) then
- ! parse the next element
- call parse_array(unit, array)
- else if ("]" == c) then
- ! end of array
- return
- end if
-
- end subroutine parse_array
-
- !
- ! PARSE STRING
- !
- function parse_string(unit) result(string)
- use fson_string_m, only: fson_string, fson_string_create, fson_string_append
- implicit none
-
- integer, intent(inout) :: unit
- type(fson_string), pointer :: string
-
- logical :: eof
- character :: c, last
-
- string => fson_string_create()
-
- do
- c = pop_char(unit, eof=eof, skip_ws=.false.)
- if (eof) then
- print *, "Expecting end of string"
- call exit(1)!
- else if ('"' == c .and. last /= '\') then !'
- exit
- else
- last = c
- call fson_string_append(string, c)
- end if
- end do
- end function parse_string
-
- !
- ! PARSE FOR CHARACTERS
- !
- subroutine parse_for_chars(unit, chars)
- integer, intent(in) :: unit
- character(len = *), intent(in) :: chars
- integer :: i, length
- logical :: eof
- character :: c
-
- length = len_trim(chars)
-
- do i = 1, length
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- print *, "ERROR: Unexpected end of file while parsing array."
- call exit (1)
- else if (c /= chars(i:i)) then
- print *, "ERROR: Unexpected character.'", c,"'", chars(i:i)
- call exit (1)
- end if
- end do
-
- end subroutine parse_for_chars
-
- !
- ! PARSE NUMBER
- !
- subroutine parse_number(unit, value)
- use fson_value_m, only: TYPE_INTEGER, TYPE_REAL
- implicit none
-
- integer, intent(inout) :: unit
- type(fson_value), pointer :: value
- logical :: eof, negative, decimal, scientific
- character :: c
- integer :: integral, exp, digit_count
- real(kind(1.0D0)) :: frac
-
- ! first character is either - or a digit
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- print *, "ERROR: Unexpected end of file while parsing number."
- call exit (1)
- else if ("-" == c) then
- negative = .true.
- else
- negative = .false.
- call push_char(c)
- end if
-
- ! parse the integral
- integral = parse_integer(unit)
-
- decimal = .false.
- scientific = .false.
-
- do
- ! first character is either - or a digit
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- print *, "ERROR: Unexpected end of file while parsing number."
- call exit (1)
- else
- select case (c)
- case (".")
- ! this is already fractional number
- if (decimal) then
- ! already found a decimal place
- print *, "ERROR: Unexpected second decimal place while parsing number."
- call exit(1)
- end if
- decimal = .true.
- frac = parse_integer(unit, digit_count)
- frac = frac / (10.0D0 ** digit_count)
- case ("e", "E")
- ! this is already an exponent number
- if (scientific) then
- ! already found a e place
- print *, "ERROR: Unexpected second exponent while parsing number."
- call exit(1)
- end if
- scientific = .true.
- ! this number has an exponent
- exp = parse_integer(unit)
-
- case default
- ! this is a integer
- if (decimal) then
-
- ! add the integral
- frac = frac + integral
-
- if (scientific) then
- ! apply exponent
- frac = frac * (10.0D0 ** exp)
- end if
-
- ! apply negative
- if (negative) then
- frac = frac * (-1)
- end if
-
- value % value_type = TYPE_REAL
- value % value_real = frac
- value % value_double = frac
- else
- if (scientific) then
- ! apply exponent
- integral = integral * (10.0D0 ** exp)
- end if
-
- ! apply negative
- if (negative) then
- integral = integral * (-1)
- end if
-
- value % value_type = TYPE_INTEGER
- value % value_integer = integral
- !+PJK
- ! Following two lines are included in case the decimal point of
- ! a floating point number has been (accidentally) left out
- value % value_real = integral
- value % value_double = integral
- !-PJK
- end if
- call push_char(c)
- exit
- end select
- end if
- end do
-
- end subroutine parse_number
-
- !
- ! PARSE INTEGER
- !
- integer(kind=8) function parse_integer(unit, digit_count) result(integral)
- integer, intent(in) :: unit
- integer, optional, intent(inout) :: digit_count
- logical :: eof, found_sign, found_digit
- character :: c
- integer :: tmp, icount, isign
- integer, parameter :: max_integer_length = 18
-
- icount = 0
- integral = 0
- isign = 1
- found_sign = .false.
- found_digit = .false.
- do
- c = pop_char(unit, eof=eof, skip_ws=.true.)
- if (eof) then
- print *, "ERROR: Unexpected end of file while parsing digit."
- call exit (1)
- else
- select case(c)
- case ("+")
- if (found_sign.or.found_digit) then
- print *, "ERROR: Miss formatted number."
- call exit(1)
- end if
- found_sign = .true.
- case ("-")
- if (found_sign.or.found_digit) then
- print *, "ERROR: Miss formatted number."
- call exit(1)
- end if
- found_sign = .true.
- isign = -1
- case ("0":"9")
- found_sign = .true.
- if (icount > max_integer_length) then
- print *, "ERROR: Too many digits for an integer."
- call exit(1)
- end if
- ! digit
- read (c, '(i1)') tmp
- ! shift
- if (icount > 0) then
- integral = integral * 10
- end if
- ! add
- integral = integral + tmp
-
- ! increase the icount
- icount = icount + 1
- case default
- if (present(digit_count)) then
- digit_count = icount
- end if
- call push_char(c)
- integral = isign * integral
- return
- end select
- end if
- end do
-
- end function parse_integer
-
- !
- ! POP CHAR
- !
- recursive character function pop_char(unit, eof, skip_ws) result(popped)
- integer, intent(in) :: unit
- logical, intent(out) :: eof
- logical, intent(in), optional :: skip_ws
-
- integer :: ios
- character :: c
- logical :: ignore
-
- !+PJK
- ios = 0
- !-PJK
- eof = .false.
- if (.not.present(skip_ws)) then
- ignore = .false.
- else
- ignore = skip_ws
- end if
-
- do
- if (pushed_index > 0) then
- ! there is a character pushed back on, most likely from the number parsing
- c = pushed_char(pushed_index:pushed_index)
- pushed_index = pushed_index - 1
- else
- read (unit=unit, fmt="(a)", advance="no", iostat=ios) c
- end if
- if (ios == end_of_record) then
- cycle
- else if (ios == end_of_file) then
- eof = .true.
- exit
- else if (iachar(c) <= 31) then ! PJK from 32
- ! non printing ascii characters
- cycle
- else if (ignore .and. c == " ") then
- cycle
- else
- popped = c
- exit
- end if
- end do
-
- end function pop_char
-
- !
- ! PUSH CHAR
- !
- subroutine push_char(c)
- character, intent(inout) :: c
- pushed_index = pushed_index + 1
- pushed_char(pushed_index:pushed_index) = c
-
- end subroutine push_char
-
-end module fson_library
diff --git a/source/fortran/init_module.f90 b/source/fortran/init_module.f90
index 59b2d69edb..c177a9fe7b 100644
--- a/source/fortran/init_module.f90
+++ b/source/fortran/init_module.f90
@@ -10,18 +10,6 @@ module init_module
contains
- subroutine init_fortran_modules
- !! Temporary routine to call initialisation routines for Fortran modules
- !! that are not wrapped by f2py and thus cannot be called from Python.
-
- use fson_library, only: init_fson_library
-
- implicit none
-
- call init_fson_library
-
- end subroutine init_fortran_modules
-
subroutine open_files
use global_variables, only: verbose, fileprefix, output_prefix
use constants, only: nout, mfile
diff --git a/tests/conftest.py b/tests/conftest.py
index 361df550f3..1a5fe21e91 100644
--- a/tests/conftest.py
+++ b/tests/conftest.py
@@ -10,7 +10,7 @@
from _pytest.fixtures import SubRequest
from system_check import system_compatible
-from process.fortran import error_handling as eh
+from process.warning_handler import WarningManager
def pytest_addoption(parser):
@@ -128,8 +128,7 @@ def initialise_error_module():
Initialise the error module initially otherwise segmentation faults can
occur when tested subroutines raise errors.
"""
- eh.init_error_handling()
- eh.initialise_error_list()
+ WarningManager.reinitialise()
@pytest.fixture
@@ -143,7 +142,7 @@ def reinitialise_error_module():
# TODO Perhaps this should be autoused by all tests? Specify use explicitly
# for now for known error-raisers
yield
- eh.init_error_handling()
+ WarningManager.reinitialise()
@pytest.fixture(autouse=True)
diff --git a/tests/integration/test_pfcoil_int.py b/tests/integration/test_pfcoil_int.py
index fe023c6e97..4ac16e0634 100644
--- a/tests/integration/test_pfcoil_int.py
+++ b/tests/integration/test_pfcoil_int.py
@@ -16,7 +16,6 @@
from process.cs_fatigue import CsFatigue
from process.fortran import build_variables as bv
from process.fortran import constants
-from process.fortran import error_handling as eh
from process.fortran import fwbs_variables as fwbsv
from process.fortran import pfcoil_module as pf
from process.fortran import pfcoil_variables as pfv
@@ -54,7 +53,6 @@ def test_pfcoil(monkeypatch, pfcoil):
monkeypatch.setattr(bv, "dr_tf_inboard", 1.4)
monkeypatch.setattr(bv, "r_tf_outboard_mid", 1.66e1)
monkeypatch.setattr(bv, "dr_bore", 2.15)
- monkeypatch.setattr(eh, "idiags", np.full(8, -999999))
monkeypatch.setattr(fwbsv, "denstl", 7.8e3)
monkeypatch.setattr(pfv, "rpf1", 0.0)
monkeypatch.setattr(pfv, "m_pf_coil_structure_total", 0.0)
@@ -190,7 +188,6 @@ def test_ohcalc(monkeypatch, reinitialise_error_module, pfcoil):
monkeypatch.setattr(bv, "hmax", 8.864)
monkeypatch.setattr(bv, "dr_cs", 6.510e-1)
monkeypatch.setattr(fwbsv, "denstl", 7.8e3)
- monkeypatch.setattr(eh, "idiags", np.full(8, 0))
monkeypatch.setattr(pfv, "n_cs_pf_coils", 5)
monkeypatch.setattr(pfv, "b_cs_peak_flat_top_end", 1.4e1)
monkeypatch.setattr(pfv, "i_cs_stress", 0)
@@ -262,7 +259,6 @@ def test_ohcalc(monkeypatch, reinitialise_error_module, pfcoil):
monkeypatch.setattr(tfv, "poisson_steel", 3.0e-1)
# Mocks for superconpf()
- monkeypatch.setattr(eh, "fdiags", np.full(8, -9.99999e5))
monkeypatch.setattr(tfv, "tmargmin_cs", 1.5)
monkeypatch.setattr(tfv, "temp_margin", 0.0)
monkeypatch.setattr(tfv, "b_crit_upper_nbti", 1.486e1)
@@ -2197,20 +2193,6 @@ def test_peakb(monkeypatch: pytest.MonkeyPatch, pfcoil: PFCoil):
monkeypatch.setattr(bv, "iohcl", 1)
monkeypatch.setattr(bv, "hmax", 9.0730900215620327)
monkeypatch.setattr(bv, "dr_cs", 0.55242000000000002)
- monkeypatch.setattr(
- eh,
- "idiags",
- np.array([
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- ]),
- )
monkeypatch.setattr(
pfv,
"r_pf_coil_inner",
@@ -2564,8 +2546,6 @@ def test_superconpf(monkeypatch: pytest.MonkeyPatch, pfcoil: PFCoil):
"""
# TODO This test would benefit from parameterisation for different SC
# materials (isumat)
- monkeypatch.setattr(eh, "fdiags", np.zeros(8))
- monkeypatch.setattr(eh, "idiags", np.zeros(8))
monkeypatch.setattr(tfv, "tmargmin_cs", 0.0)
monkeypatch.setattr(tfv, "temp_margin", 0.0)
monkeypatch.setattr(tfv, "b_crit_upper_nbti", 0.0)
@@ -2746,34 +2726,6 @@ def test_induct(pfcoil: PFCoil, monkeypatch: pytest.MonkeyPatch):
"""
monkeypatch.setattr(bv, "iohcl", 1)
monkeypatch.setattr(bv, "dr_cs", 0.55242000000000002)
- monkeypatch.setattr(
- eh,
- "fdiags",
- np.array([
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- ]),
- )
- monkeypatch.setattr(
- eh,
- "idiags",
- np.array([
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- -999999,
- ]),
- )
monkeypatch.setattr(pfv, "n_cs_pf_coils", 7)
monkeypatch.setattr(
pfv,
diff --git a/tests/integration/test_vmcon.py b/tests/integration/test_vmcon.py
index c5613e008c..abe636239f 100644
--- a/tests/integration/test_vmcon.py
+++ b/tests/integration/test_vmcon.py
@@ -13,16 +13,12 @@
import pytest
from process.evaluators import Evaluators
-from process.fortran import error_handling
from process.init import init_all_module_vars
from process.solver import get_solver
# Debug-level terminal output logging
logger = logging.getLogger(__name__)
-# Provide helpful errors in the event of a failed Vmcon test
-error_handling.initialise_error_list()
-
@pytest.fixture(autouse=True)
def reinit():
diff --git a/tests/unit/test_costs_1990.py b/tests/unit/test_costs_1990.py
index 0ff4454224..41d2bb2337 100644
--- a/tests/unit/test_costs_1990.py
+++ b/tests/unit/test_costs_1990.py
@@ -26,7 +26,6 @@
times_variables,
vacuum_variables,
)
-from process.fortran import error_handling as eh
from process.fortran import fwbs_variables as fv
from process.fortran import heat_transport_variables as htv
@@ -41,24 +40,6 @@ def costs():
return Costs()
-@pytest.fixture
-def initialise_error_module(monkeypatch):
- """pytest fixture to initialise error module
-
- Any routine which can raise an error should initialise
- the error module otherwise segmentation faults can occur.
-
- This fixture also resets the `fdiags` array to 0's.
-
- :param monkeypatch: Mock fixture
- :type monkeypatch: object
- """
- eh.init_error_handling()
- eh.initialise_error_list()
- # monkeypatch.setattr(eh, 'fdiags', np.zeros(8))
- # monkeypatch.setattr(eh, 'errors_on', False)
-
-
def acc2261_param(**kwargs):
"""Make parameters for a single acc2261() test.
@@ -5395,7 +5376,7 @@ class Acc2253Param(NamedTuple):
),
),
)
-def test_acc2253_urt(acc2253param, monkeypatch, costs, initialise_error_module):
+def test_acc2253_urt(acc2253param, monkeypatch, costs):
"""
Automatically generated Regression Unit Test for acc2253.