From bab12e20109abfdd0761f30c3db520081832646e Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 6 Aug 2024 09:45:09 +0100 Subject: [PATCH 1/8] Added new bootstap model switch --- process/physics.py | 7 +++++++ source/fortran/input.f90 | 2 +- 2 files changed, 8 insertions(+), 1 deletion(-) diff --git a/process/physics.py b/process/physics.py index ac895a08cc..482b5561e5 100644 --- a/process/physics.py +++ b/process/physics.py @@ -1002,6 +1002,8 @@ def physics(self): current_drive_variables.bootipf = current_drive_variables.bscf_wilson elif physics_variables.ibss == 4: current_drive_variables.bootipf = current_drive_variables.bscf_sauter + elif physics_variables.ibss == 5: + current_drive_variables.bootipf = current_drive_variables.bscf_sakai else: error_handling.idiags[0] = physics_variables.ibss error_handling.report_error(75) @@ -4115,6 +4117,11 @@ def outplas(self): self.outfile, " (Sauter et al bootstrap current fraction model used)", ) + elif physics_variables.ibss == 5: + po.ocmmnt( + self.outfile, + " (Sakai et al bootstrap current fraction model used)", + ) if physics_variables.idia == 0: po.ocmmnt( diff --git a/source/fortran/input.f90 b/source/fortran/input.f90 index 9654d2e262..f5dd77b16b 100644 --- a/source/fortran/input.f90 +++ b/source/fortran/input.f90 @@ -619,7 +619,7 @@ subroutine parse_input_file(in_file,out_file,show_changes) call parse_real_variable('taumax', taumax, 0.1D0, 100.0D0, & 'Maximum allowed energy confinement time (s)') case ('ibss') - call parse_int_variable('ibss', ibss, 1, 4, & + call parse_int_variable('ibss', ibss, 1, 5, & 'Switch for bootstrap scaling') case ('iculbl') call parse_int_variable('iculbl', iculbl, 0, 3, & From bd71e4b693f8d87cd716b58e120f54475dcfca6d Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 6 Aug 2024 10:01:34 +0100 Subject: [PATCH 2/8] feat: Add Sakai bootstrap current fraction model This commit adds a new bootstrap current fraction model based on the Sakai formula. The `bootstrap_fraction_sakai` method is implemented in the `Physics` class in the `physics.py` file. This model calculates the bootstrap fraction using a set of parameters such as plasma beta, bootstrap current coefficient, safety factor, neoclassical viscosity parameters, and epsilon. The formula is derived from the work of Sakai et al. [^1^]. [^1^]: Sakai, R., Fujita, T., Okamoto, A. (2019). Derivation of bootstrap current fraction scaling formula for 0-D system code analysis. Fusion Engineering and Design, 149, 111322. https://doi.org/10.1016/j.fusengdes.2019.111322 --- process/physics.py | 58 +++++++++++++++++++++- source/fortran/current_drive_variables.f90 | 3 ++ 2 files changed, 60 insertions(+), 1 deletion(-) diff --git a/process/physics.py b/process/physics.py index 482b5561e5..1fd3fcf367 100644 --- a/process/physics.py +++ b/process/physics.py @@ -990,6 +990,10 @@ def physics(self): current_drive_variables.cboot * self.bootstrap_fraction_sauter() ) + current_drive_variables.bscf_sakai = ( + current_drive_variables.cboot * self.bootstrap_fraction_sakai() + ) + if current_drive_variables.bscfmax < 0.0e0: current_drive_variables.bootipf = abs(current_drive_variables.bscfmax) current_drive_variables.plasipf = current_drive_variables.bootipf @@ -4066,6 +4070,14 @@ def outplas(self): current_drive_variables.bscf_wilson, "OP ", ) + + po.ovarrf( + self.outfile, + "Bootstrap fraction (Sakai)", + "(bscf_sakai)", + current_drive_variables.bscf_sakai, + "OP ", + ) po.ovarrf( self.outfile, "Diamagnetic fraction (Hender)", @@ -4121,7 +4133,7 @@ def outplas(self): po.ocmmnt( self.outfile, " (Sakai et al bootstrap current fraction model used)", - ) + ) if physics_variables.idia == 0: po.ocmmnt( @@ -4669,6 +4681,50 @@ def bootstrap_fraction_sauter(): return np.sum(da * jboot, axis=0) / physics_variables.plascur + @staticmethod + def bootstrap_fraction_sakai( + betap: float, + cboot: float, + q95: float, + q0: float, + alphan: float, + alphat: float, + eps: float, + rli: float, + ) -> float: + """ + Calculate the bootstrap fraction using the Sakai formula. + + Parameters: + betap (float): Plasma beta parameter. + cboot (float): Bootstrap current coefficient. + plascur (float): Plasma current. + q95 (float): Safety factor at 95% of the plasma radius. + q0 (float): Safety factor at the magnetic axis. + alphan (float): Neoclassical toroidal viscosity parameter. + alphat (float): Neoclassical poloidal viscosity parameter. + eps (float): Epsilon parameter. + + Returns: + float: The calculated bootstrap fraction. + + Notes: + The profile assumed for the alphan anf alpat indexes is only a prabolic profile without a pedestal. + + References: + Ryosuke Sakai, Takaaki Fujita, Atsushi Okamoto, Derivation of bootstrap current fraction scaling formula for 0-D system code analysis, + Fusion Engineering and Design, Volume 149, 2019, 111322, ISSN 0920-3796, + https://doi.org/10.1016/j.fusengdes.2019.111322. + """ + return ( + 10 ** (0.951 * eps - 0.948) + * betap ** (1.226 * eps + 1.584) + * rli ** (-0.184 - 0.285) + * (q95 / q0) ** (-0.042 * eps - 0.02) + * alphan ** (0.13 * eps + 0.05) + * alphat ** (0.502 * eps - 0.273) + ) + def fhfac(self, is_): """Function to find H-factor for power balance author: P J Knight, CCFE, Culham Science Centre diff --git a/source/fortran/current_drive_variables.f90 b/source/fortran/current_drive_variables.f90 index 45b357e961..d8db4048a5 100644 --- a/source/fortran/current_drive_variables.f90 +++ b/source/fortran/current_drive_variables.f90 @@ -42,6 +42,9 @@ module current_drive_variables real(dp) :: bscf_wilson !! bootstrap current fraction, Wilson et al model + real(dp) :: bscf_sakai + !! Bootstrap current fraction, Sakai et al model + real(dp) :: cboot !! bootstrap current fraction multiplier (`ibss=1`) From 7fe362db902373f9ef41498c9b8b22cbc74ef4ec Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 6 Aug 2024 10:20:41 +0100 Subject: [PATCH 3/8] feat: Implement Sakai bootstrap current fraction model --- process/physics.py | 23 +++++++++++++++-------- 1 file changed, 15 insertions(+), 8 deletions(-) diff --git a/process/physics.py b/process/physics.py index 1fd3fcf367..473f69134f 100644 --- a/process/physics.py +++ b/process/physics.py @@ -991,7 +991,16 @@ def physics(self): ) current_drive_variables.bscf_sakai = ( - current_drive_variables.cboot * self.bootstrap_fraction_sakai() + current_drive_variables.cboot + * self.bootstrap_fraction_sakai( + betap=physics_variables.betap, + q95=physics_variables.q95, + q0=physics_variables.q0, + alphan=physics_variables.alphan, + alphat=physics_variables.alphat, + eps=physics_variables.eps, + rli=physics_variables.rli, + ) ) if current_drive_variables.bscfmax < 0.0e0: @@ -1008,6 +1017,7 @@ def physics(self): current_drive_variables.bootipf = current_drive_variables.bscf_sauter elif physics_variables.ibss == 5: current_drive_variables.bootipf = current_drive_variables.bscf_sakai + print(current_drive_variables.bscf_sakai) else: error_handling.idiags[0] = physics_variables.ibss error_handling.report_error(75) @@ -4070,7 +4080,7 @@ def outplas(self): current_drive_variables.bscf_wilson, "OP ", ) - + po.ovarrf( self.outfile, "Bootstrap fraction (Sakai)", @@ -4684,7 +4694,6 @@ def bootstrap_fraction_sauter(): @staticmethod def bootstrap_fraction_sakai( betap: float, - cboot: float, q95: float, q0: float, alphan: float, @@ -4697,12 +4706,10 @@ def bootstrap_fraction_sakai( Parameters: betap (float): Plasma beta parameter. - cboot (float): Bootstrap current coefficient. - plascur (float): Plasma current. q95 (float): Safety factor at 95% of the plasma radius. q0 (float): Safety factor at the magnetic axis. - alphan (float): Neoclassical toroidal viscosity parameter. - alphat (float): Neoclassical poloidal viscosity parameter. + alphan (float): Density profile index + alphat (float): Temperature profile index eps (float): Epsilon parameter. Returns: @@ -4712,7 +4719,7 @@ def bootstrap_fraction_sakai( The profile assumed for the alphan anf alpat indexes is only a prabolic profile without a pedestal. References: - Ryosuke Sakai, Takaaki Fujita, Atsushi Okamoto, Derivation of bootstrap current fraction scaling formula for 0-D system code analysis, + Ryosuke Sakai, Takaaki Fujita, Atsushi Okamoto, Derivation of bootstrap current fraction scaling formula for 0-D system code analysis, Fusion Engineering and Design, Volume 149, 2019, 111322, ISSN 0920-3796, https://doi.org/10.1016/j.fusengdes.2019.111322. """ From 671d886313b20ea30870288ed2bbda619d17fbd1 Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 6 Aug 2024 13:54:23 +0100 Subject: [PATCH 4/8] Update docstring --- process/physics.py | 12 +++++++++--- source/fortran/physics_variables.f90 | 1 + 2 files changed, 10 insertions(+), 3 deletions(-) diff --git a/process/physics.py b/process/physics.py index 473f69134f..edfb388198 100644 --- a/process/physics.py +++ b/process/physics.py @@ -4705,18 +4705,24 @@ def bootstrap_fraction_sakai( Calculate the bootstrap fraction using the Sakai formula. Parameters: - betap (float): Plasma beta parameter. + betap (float): Plasma poloidal beta. q95 (float): Safety factor at 95% of the plasma radius. q0 (float): Safety factor at the magnetic axis. alphan (float): Density profile index alphat (float): Temperature profile index - eps (float): Epsilon parameter. + eps (float): Inverse aspect ratio. Returns: float: The calculated bootstrap fraction. Notes: - The profile assumed for the alphan anf alpat indexes is only a prabolic profile without a pedestal. + The profile assumed for the alphan anf alpat indexes is only a prabolic profile without a pedestal (L-mode). + The Root Mean Squared Error for the fitting database of this formula was 0.025 + Concentrating on the positive shear plasmas using the ACCOME code equilibria with the fully non-inductively driven + conditions with neutral beam (NB) injection only are calculated. + The electron temperature and the ion temperature were assumed to be equal + This can be used for all apsect ratios. + The diamagnetic fraction is included in this formula. References: Ryosuke Sakai, Takaaki Fujita, Atsushi Okamoto, Derivation of bootstrap current fraction scaling formula for 0-D system code analysis, diff --git a/source/fortran/physics_variables.f90 b/source/fortran/physics_variables.f90 index 995df178f9..7999fc5b7d 100644 --- a/source/fortran/physics_variables.f90 +++ b/source/fortran/physics_variables.f90 @@ -249,6 +249,7 @@ module physics_variables !! - =2 for Nevins et al general scaling !! - =3 for Wilson et al numerical scaling !! - =4 for Sauter et al scaling + !! - =5 for Sakai et al scaling integer :: iculbl !! switch for beta limit scaling (`constraint equation 24`) From 801b8a78610a6d81821d308a45eeac0e93b27fe0 Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 6 Aug 2024 14:03:12 +0100 Subject: [PATCH 5/8] feat: Add Sakai bootstrap current fraction model docs --- documentation/proc-pages/physics-models/plasma_current.md | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/documentation/proc-pages/physics-models/plasma_current.md b/documentation/proc-pages/physics-models/plasma_current.md index a1395a0edd..d2204fbf06 100644 --- a/documentation/proc-pages/physics-models/plasma_current.md +++ b/documentation/proc-pages/physics-models/plasma_current.md @@ -68,6 +68,7 @@ existence of pedestals, whereas the Sauter et al. scaling | 2 | General scaling -- To use a more general scaling method, set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). | 3 | Numerically fitted scaling [^8] -- To use a numerically fitted scaling method, valid for all aspect ratios, set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). | 4 | Sauter, Angioni and Lin-Liu scaling [^9] [^10] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). +| 5 | Sakai, Fujita and Okamoto scaling [^11] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). !!! Note "Fixed Bootstrap Current" Direct input -- To input the bootstrap current fraction directly, set `bscfmax` @@ -114,4 +115,5 @@ Task meeting EU-JA #16, Fusion for Energy, Garching, 24--25 June 2014 [^7]: Menard et al. (2016), Nuclear Fusion, 56, 106023 [^8]: H.R. Wilson, Nuclear Fusion **32** (1992) 257 [^9]: O. Sauter, C. Angioni and Y.R. Lin-Liu, Physics of Plasmas **6** (1999) 2834 -[^10]: O. Sauter, C. Angioni and Y.R. Lin-Liu, Physics of Plasmas **9** (2002) 5140 \ No newline at end of file +[^10]: O. Sauter, C. Angioni and Y.R. Lin-Liu, Physics of Plasmas **9** (2002) 5140 +[^11]: Ryosuke Sakai, Takaaki Fujita, Atsushi Okamoto, Derivation of bootstrap current fraction scaling formula for 0-D system code analysis, Fusion Engineering and Design, Volume 149, 2019, 111322, ISSN 0920-3796, https://doi.org/10.1016/j.fusengdes.2019.111322. From e9a14d7482d9d3942e3df0b0c50f6f3f39051b7c Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 6 Aug 2024 14:53:54 +0100 Subject: [PATCH 6/8] Added unit test for large tokamak and ST case --- process/physics.py | 2 +- tests/unit/test_physics.py | 84 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 85 insertions(+), 1 deletion(-) diff --git a/process/physics.py b/process/physics.py index edfb388198..1e37148fb4 100644 --- a/process/physics.py +++ b/process/physics.py @@ -4732,7 +4732,7 @@ def bootstrap_fraction_sakai( return ( 10 ** (0.951 * eps - 0.948) * betap ** (1.226 * eps + 1.584) - * rli ** (-0.184 - 0.285) + * rli ** (-0.184 - 0.282) * (q95 / q0) ** (-0.042 * eps - 0.02) * alphan ** (0.13 * eps + 0.05) * alphat ** (0.502 * eps - 0.273) diff --git a/tests/unit/test_physics.py b/tests/unit/test_physics.py index 125569741a..32bbac8ce6 100644 --- a/tests/unit/test_physics.py +++ b/tests/unit/test_physics.py @@ -468,6 +468,90 @@ def test_bootstrap_fraction_sauter(bootstrapfractionsauterparam, monkeypatch, ph assert bfs == pytest.approx(bootstrapfractionsauterparam.expected_bfs) +class BootstrapFractionSakaiParam(NamedTuple): + betap: Any = None + + q95: Any = None + + q0: Any = None + + alphan: Any = None + + alphat: Any = None + + eps: Any = None + + rli: Any = None + + expected_bfs: Any = None + + +@pytest.mark.parametrize( + "bootstrapfractionsakaiparam", + ( + BootstrapFractionSakaiParam( + betap=1.3184383457774960, + q95=3.5151046634673557, + q0=1.0, + alphan=1.0, + alphat=1.45, + eps=1 / 3, + rli=1.2098126022585098, + expected_bfs=0.34204201506155418, + ), + BootstrapFractionSakaiParam( + betap=1.1701245502231756, + q95=5.1746754543339177, + q0=2.0, + alphan=0.9, + alphat=1.3999999999999999, + eps=1 / 1.8, + rli=0.3, + expected_bfs=0.90349498124262029, + ), + ), +) +def test_bootstrap_fraction_sakai(bootstrapfractionsakaiparam, monkeypatch, physics): + """ + Automatically generated Regression Unit Test for bootstrap_fraction_sakai. + + This test was generated using data from tests/regression/input_files/large_tokamak.IN.DAT. + and tests/regression/input_files/st_regression.IN.DAT. + + :param bootstrapfractionsauterparam: the data used to mock and assert in this test. + :type bootstrapfractionsauterparam: bootstrapfractionsauterparam + + :param monkeypatch: pytest fixture used to mock module/class variables + :type monkeypatch: _pytest.monkeypatch.monkeypatch + """ + + monkeypatch.setattr(physics_variables, "betap", bootstrapfractionsakaiparam.betap) + + monkeypatch.setattr(physics_variables, "q95", bootstrapfractionsakaiparam.q95) + + monkeypatch.setattr(physics_variables, "q0", bootstrapfractionsakaiparam.q0) + + monkeypatch.setattr(physics_variables, "alphan", bootstrapfractionsakaiparam.alphan) + + monkeypatch.setattr(physics_variables, "alphat", bootstrapfractionsakaiparam.alphat) + + monkeypatch.setattr(physics_variables, "eps", bootstrapfractionsakaiparam.eps) + + monkeypatch.setattr(physics_variables, "rli", bootstrapfractionsakaiparam.rli) + + bfs = physics.bootstrap_fraction_sakai( + betap=bootstrapfractionsakaiparam.betap, + q95=bootstrapfractionsakaiparam.q95, + q0=bootstrapfractionsakaiparam.q0, + alphan=bootstrapfractionsakaiparam.alphan, + alphat=bootstrapfractionsakaiparam.alphat, + eps=bootstrapfractionsakaiparam.eps, + rli=bootstrapfractionsakaiparam.rli, + ) + + assert bfs == pytest.approx(bootstrapfractionsakaiparam.expected_bfs) + + class CulcurParam(NamedTuple): normalised_total_beta: Any = None From f1414f47667d9b0258b9b31f597eab01524a8d9c Mon Sep 17 00:00:00 2001 From: mn3981 Date: Wed, 7 Aug 2024 09:51:47 +0100 Subject: [PATCH 7/8] feat: Update Sakai bootstrap current calculation This commit modifies the `Physics` class in the `physics.py` file to update the calculation of the Sakai bootstrap current. The Sakai formula states that the ACCOME dataset used has the toroidal diamagnetic current included in the bootstrap current. Therefore, the commit removes the PROCESS calculated diamagnetic current from the bootstrap current calculation. This change improves the accuracy of the bootstrap current fraction model. --- .../proc-pages/physics-models/plasma_current.md | 2 +- process/physics.py | 9 +++++++-- 2 files changed, 8 insertions(+), 3 deletions(-) diff --git a/documentation/proc-pages/physics-models/plasma_current.md b/documentation/proc-pages/physics-models/plasma_current.md index d2204fbf06..0a0a40a588 100644 --- a/documentation/proc-pages/physics-models/plasma_current.md +++ b/documentation/proc-pages/physics-models/plasma_current.md @@ -68,7 +68,7 @@ existence of pedestals, whereas the Sauter et al. scaling | 2 | General scaling -- To use a more general scaling method, set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). | 3 | Numerically fitted scaling [^8] -- To use a numerically fitted scaling method, valid for all aspect ratios, set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). | 4 | Sauter, Angioni and Lin-Liu scaling [^9] [^10] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). -| 5 | Sakai, Fujita and Okamoto scaling [^11] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). +| 5 | Sakai, Fujita and Okamoto scaling [^11] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). The model includes the toroidal diamagnetic current in the calculation due to the dataset, so the PROCESS diamagnetic current fraction is removed from the calculated value for greater self-consistency. !!! Note "Fixed Bootstrap Current" Direct input -- To input the bootstrap current fraction directly, set `bscfmax` diff --git a/process/physics.py b/process/physics.py index 1e37148fb4..372b71d3f3 100644 --- a/process/physics.py +++ b/process/physics.py @@ -1016,8 +1016,11 @@ def physics(self): elif physics_variables.ibss == 4: current_drive_variables.bootipf = current_drive_variables.bscf_sauter elif physics_variables.ibss == 5: - current_drive_variables.bootipf = current_drive_variables.bscf_sakai - print(current_drive_variables.bscf_sakai) + # Sakai states that the ACCOME dataset used has the toridal diamagnetic current included in the bootstrap current + # So we remove the PROCESS calculated diamagnetic current from the bootstrap current + current_drive_variables.bootipf = ( + current_drive_variables.bscf_sakai - current_drive_variables.diaipf + ) else: error_handling.idiags[0] = physics_variables.ibss error_handling.report_error(75) @@ -4729,6 +4732,8 @@ def bootstrap_fraction_sakai( Fusion Engineering and Design, Volume 149, 2019, 111322, ISSN 0920-3796, https://doi.org/10.1016/j.fusengdes.2019.111322. """ + # Sakai states that the ACCOME dataset used has the toridal diamagnetic current included in the bootstrap current + # So we remove the PROCESS calculated diamagnetic current at the function call point return ( 10 ** (0.951 * eps - 0.948) * betap ** (1.226 * eps + 1.584) From b51dbf4f9be3ce44aa1b1e12767ad3d2552c4d3e Mon Sep 17 00:00:00 2001 From: mn3981 Date: Thu, 8 Aug 2024 15:05:08 +0100 Subject: [PATCH 8/8] Added level 3 warning about calculating diamagnetic current when using the Sakai scaling --- documentation/proc-pages/physics-models/plasma_current.md | 2 +- process/physics.py | 8 +++----- process/utilities/errorlist.json | 7 ++++++- source/fortran/initial.f90 | 6 +++++- 4 files changed, 15 insertions(+), 8 deletions(-) diff --git a/documentation/proc-pages/physics-models/plasma_current.md b/documentation/proc-pages/physics-models/plasma_current.md index 0a0a40a588..1a2a3f3029 100644 --- a/documentation/proc-pages/physics-models/plasma_current.md +++ b/documentation/proc-pages/physics-models/plasma_current.md @@ -68,7 +68,7 @@ existence of pedestals, whereas the Sauter et al. scaling | 2 | General scaling -- To use a more general scaling method, set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). | 3 | Numerically fitted scaling [^8] -- To use a numerically fitted scaling method, valid for all aspect ratios, set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). | 4 | Sauter, Angioni and Lin-Liu scaling [^9] [^10] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). -| 5 | Sakai, Fujita and Okamoto scaling [^11] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). The model includes the toroidal diamagnetic current in the calculation due to the dataset, so the PROCESS diamagnetic current fraction is removed from the calculated value for greater self-consistency. +| 5 | Sakai, Fujita and Okamoto scaling [^11] -- Set `bscfmax` to the maximum required bootstrap current fraction ($\leq 1$). The model includes the toroidal diamagnetic current in the calculation due to the dataset, so `idia = 0` can only be used with it !!! Note "Fixed Bootstrap Current" Direct input -- To input the bootstrap current fraction directly, set `bscfmax` diff --git a/process/physics.py b/process/physics.py index 372b71d3f3..e616abff60 100644 --- a/process/physics.py +++ b/process/physics.py @@ -1017,10 +1017,8 @@ def physics(self): current_drive_variables.bootipf = current_drive_variables.bscf_sauter elif physics_variables.ibss == 5: # Sakai states that the ACCOME dataset used has the toridal diamagnetic current included in the bootstrap current - # So we remove the PROCESS calculated diamagnetic current from the bootstrap current - current_drive_variables.bootipf = ( - current_drive_variables.bscf_sakai - current_drive_variables.diaipf - ) + # So the diamagnetic current calculation should be turned off when using, (idia = 0). + current_drive_variables.bootipf = current_drive_variables.bscf_sakai else: error_handling.idiags[0] = physics_variables.ibss error_handling.report_error(75) @@ -4733,7 +4731,7 @@ def bootstrap_fraction_sakai( https://doi.org/10.1016/j.fusengdes.2019.111322. """ # Sakai states that the ACCOME dataset used has the toridal diamagnetic current included in the bootstrap current - # So we remove the PROCESS calculated diamagnetic current at the function call point + # So the diamganetic current should not be calculated with this. idia = 0 return ( 10 ** (0.951 * eps - 0.948) * betap ** (1.226 * eps + 1.584) diff --git a/process/utilities/errorlist.json b/process/utilities/errorlist.json index dafabd320e..69611229a3 100644 --- a/process/utilities/errorlist.json +++ b/process/utilities/errorlist.json @@ -8,7 +8,7 @@ "comment2": [ "Increment n_errortypes if an error is added to this list" ], - "n_errortypes": 283, + "n_errortypes": 284, "errors": [ { "no": 1, @@ -1424,6 +1424,11 @@ "no": 283, "level": 3, "message": "[tfcoil]: Can only have i_tf_turns_integer = 1 with i_tf_wp_geom = 0" + }, + { + "no": 284, + "level": 3, + "message": "CHECK: idia = 0 should be used with the Sakai plasma current scaling." } ] } diff --git a/source/fortran/initial.f90 b/source/fortran/initial.f90 index 9a1512d740..97f0b0b45d 100755 --- a/source/fortran/initial.f90 +++ b/source/fortran/initial.f90 @@ -258,7 +258,7 @@ subroutine check use physics_variables, only: aspect, fdeut, fgwped, fhe3, & fgwsep, ftrit, ibss, i_single_null, icurr, idivrt, ishape, & iradloss, isc, ipedestal, ilhthresh, itart, nesep, rhopedn, rhopedt, & - rnbeam, neped, te, tauee_in, tesep, teped, itartpf, ftar + rnbeam, neped, te, tauee_in, tesep, teped, itartpf, ftar, idia use pulse_variables, only: lpulse use reinke_variables, only: fzactual, impvardiv use tfcoil_variables, only: casthi, casthi_is_fraction, casths, i_tf_sup, & @@ -828,6 +828,10 @@ subroutine check call report_error(283) end if + if ( ibss == 5 .and. idia /= 0 ) then + call report_error(284) + end if + ! Setting t_cable_tf_is_input to true if t_cable_tf is an input if ( abs(t_cable_tf) < epsilon(t_cable_tf) ) then t_cable_tf_is_input = .false.