Merge branch 'main' into 2953-additional-functionality-to-plot_scanspy

CMakeLists.txt

@@ -123,6 +123,10 @@ LIST(APPEND PROCESS_SRCS
     reinke_variables.f90
     neoclassics_module.f90
     blanket_library.f90
+    stellarator_variables.f90
+    stellarator.f90
+    stellarator_configuration.f90
+    stellarator_fwbs.f90
 )
 
 PREPROCESS()

documentation/proc-pages/fusion-devices/stellarator.md

@@ -1,8 +1,5 @@
 # Stellarator Model
 
-!!! danger inline end "Stellarator removed" 
-    The Stellarator model was removed in issue #1853 because of concerns around licensing issues with a number of sections of code used in the model. Issue #1854 aims to reintroduce the model after these licensing issues are addressed.
-
 The code has the ability to perform calculations based on the physics and engineering of a stellarator, which, although being a toroidal device, is radically different in a number of ways from a tokamak.
 
 The model is largely based on W7-X and the HELIAS 5-B stellarator power plant design[^1] (Figure 1) and related modelling that has been performed by IPP Greifswald[^2] [^3] [^4] [^5]

process/caller.py

@@ -37,6 +37,13 @@ def call_models(self, xc):
         # Convert variables
         ft.define_iteration_variables.convxc(xc, nvars)
 
+        # Perform the various function calls
+        # Stellarator caller
+        if ft.stellarator_variables.istell != 0:
+            self.models.stellarator.run(output=False)
+            # TODO Is this return safe?
+            return
+
         # Inertial Fusion Energy calls
         if ft.ife_variables.ife != 0:
             self.models.ife.run(output=False)

process/evaluators.py

@@ -4,6 +4,7 @@
 from process.fortran import cost_variables as cv
 from process.fortran import numerics
 from process.fortran import physics_variables as pv
+from process.fortran import stellarator_variables as sv
 from process.fortran import times_variables as tv
 from process.fortran import function_evaluator
 import numpy as np
@@ -56,22 +57,22 @@ def fcnvmc1(self, n, m, xv, ifail):
         self.caller.call_models(xv)
 
         # Convergence loop to ensure burn time consistency
-
-        loop = 0
-        while (loop < 10) and (
-            abs((tv.tburn - tv.tburn0) / max(tv.tburn, 0.01)) > 0.001
-        ):
-            loop += 1
-            self.caller.call_models(xv)
-            if gv.verbose == 1:
-                print("Internal tburn consistency check: ", tv.tburn, tv.tburn0)
-
-        if loop >= 10:
-            print(
-                "Burn time values are not consistent in iteration: ",
-                numerics.nviter,
-            )
-            print("tburn, tburn0: ", tv.tburn, tv.tburn0)
+        if sv.istell == 0:
+            loop = 0
+            while (loop < 10) and (
+                abs((tv.tburn - tv.tburn0) / max(tv.tburn, 0.01)) > 0.001
+            ):
+                loop += 1
+                self.caller.call_models(xv)
+                if gv.verbose == 1:
+                    print("Internal tburn consistency check: ", tv.tburn, tv.tburn0)
+
+            if loop >= 10:
+                print(
+                    "Burn time values are not consistent in iteration: ",
+                    numerics.nviter,
+                )
+                print("tburn, tburn0: ", tv.tburn, tv.tburn0)
 
         # Evaluate figure of merit (objective function)
         objf = function_evaluator.funfom()

process/io/sankey_funcs.py

@@ -17,7 +17,6 @@
 def plot_full_sankey(
     mfilename="MFILE.DAT",
 ):  # Plots the power flow from PROCESS as a Sankey Diagram
-
     # ------------------------------- Pulling values from the MFILE -------------------------------
 
     m_file = MFile(mfilename)
@@ -98,7 +97,6 @@ def plot_full_sankey(
 
     # Loop 1 to get 'Plasma Heating' branch tip coords; loop 2 to match 'PLASMA' branch
     for _ in range(2):
-
         # The visual settings of the Sankey Plot
         plt.rcParams.update({"font.size": 9})
         fig = plt.figure()
@@ -451,7 +449,6 @@ def plot_full_sankey(
 
 
 def plot_sankey(mfilename="MFILE.DAT"):  # Plot simplified power flow Sankey Diagram
-
     # ------------------------------- Pulling values from the MFILE -------------------------------
 
     m_file = MFile(mfilename)
@@ -522,7 +519,7 @@ def plot_sankey(mfilename="MFILE.DAT"):  # Plot simplified power flow Sankey Dia
     pthermfw_blkt = m_file.data["pthermfw_blkt"].get_scan(
         -1
     )  # Heat for electricity (MW)
-    htpmw_fw_blkt = m_file.data["htpmw_blkt_tot"].get_scan(
+    htpmw_fw_blkt = m_file.data["htpmw_fw_blkt"].get_scan(
         -1
     )  # 1st wall & blanket pumping (MW)
     pthermmw_p = pthermfw_blkt - htpmw_fw_blkt  # Heat - pumping power (MW)
@@ -563,7 +560,6 @@ def plot_sankey(mfilename="MFILE.DAT"):  # Plot simplified power flow Sankey Dia
 
     # Loop 1 to get 'Plasma Heating' branch tip coords; loop 2 to match 'PLASMA' branch
     for _ in range(2):
-
         # ------------------------------------ Visual Settings ------------------------------------
 
         plt.rcParams.update({"font.size": 9})  # Setting font size to 9

process/main.py

@@ -49,6 +49,7 @@
 from process.pulse import Pulse
 from process.scan import Scan
 from process import final
+from process.stellarator import Stellarator
 from process.structure import Structure
 from process.build import Build
 from process.utilities.f2py_string_patch import string_to_f2py_compatible
@@ -580,11 +581,20 @@ def __init__(self):
         self.costs = Costs()
         self.ife = IFE(availability=self.availability, costs=self.costs)
         self.plasma_profile = PlasmaProfile()
-        self.costs_2015 = Costs2015()
-        self.physics = Physics(plasma_profile=self.plasma_profile)
         self.fw = Fw()
         self.blanket_library = BlanketLibrary(fw=self.fw)
         self.ccfe_hcpb = CCFE_HCPB(blanket_library=self.blanket_library)
+        self.stellarator = Stellarator(
+            availability=self.availability,
+            buildings=self.buildings,
+            vacuum=self.vacuum,
+            costs=self.costs,
+            power=self.power,
+            plasma_profile=self.plasma_profile,
+            hcpb=self.ccfe_hcpb,
+        )
+        self.costs_2015 = Costs2015()
+        self.physics = Physics(plasma_profile=self.plasma_profile)
         self.dcll = DCLL(blanket_library=self.blanket_library)
 
 

process/output.py

@@ -17,6 +17,11 @@ def write(models, outfile):
     # 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)
+        return
+
     #  Call IFE output routine instead if relevant
     if ft.ife_variables.ife != 0:
         models.ife.run(output=True)

process/solver.py

@@ -169,7 +169,11 @@ def solve(self) -> int:
             B = np.identity(numerics.nvar) * self.b
 
         def _solver_callback(i: int, _x, _result, convergence_param: float):
-            print(f"{i+1} | Convergence Parameter: {convergence_param:.3E}", end="\r")
+            print(
+                f"{i+1} | Convergence Parameter: {convergence_param:.3E}",
+                end="\r",
+                flush=True,
+            )
 
         try:
             x, _, _, res = solve(