Add plotting and documentation for cyclotron radiation sink.

Includes new post-processed variable, leading to regeneration of all sim tests. PiperOrigin-RevId: 712888033
google-deepmind · Jan 8, 2025 · f7427f1 · f7427f1
1 parent 9aa2340
commit f7427f1
Show file tree

Hide file tree

Showing 66 changed files with 664 additions and 38 deletions.
diff --git a/docs/configuration.rst b/docs/configuration.rst
@@ -931,7 +931,7 @@ Time dependent Gaussian pellet source. No first-principle-based model is yet imp
   Total particle source in units of particles/s
 
 generic_particle_source
-^^^^^^^^^^^^^^^^^^^
+^^^^^^^^^^^^^^^^^^^^^^^
 
 Time dependent Gaussian particle source. No first-principle-based model is yet implemented in TORAX.
 
@@ -988,6 +988,26 @@ Bremsstrahlung model from Wesson, with an optional correction for relativistic e
 
 ``use_relativistic_correction`` (bool = False)
 
+cyclotron_radiation_heat_sink
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Cyclotron radiation model from Albajar NF 2001 with a deposition profile from Artaud NF 2018.
+
+``mode`` (str = 'model')
+
+``wall_reflection_coeff`` (float = 0.9)
+  Machine-dependent dimensionless parameter corresponding to the fraction of
+  cyclotron radiation reflected off the wall and reabsorbed by the plasma.
+
+``beta_min`` (float = 0.5)
+
+``beta_max`` (float = 8.0)
+
+``beta_grid_size`` (int = 32)
+  beta in this context is a variable in the temperature profile parameterization used
+  in the Albajar model. The parameter is fit with simple grid search performed over
+  the range ``[beta_min, beta_max]``, with ``beta_grid_size`` uniformly spaced steps.
+
 electron_cyclotron_source
 ^^^^^^^^^^^^^^^^^^^^^^^^^
 Electron-cyclotron heating and current drive, based on the local efficiency model in `Lin-Liu et al., 2003 <https://doi.org/10.1063/1.1610472>`_.

diff --git a/docs/physics_models.rst b/docs/physics_models.rst
@@ -409,17 +409,25 @@ Presently, TORAX provides three built-in formula-based particle sources for the
 
 Radiation
 ---------
-Currently, TORAX only has a dedicated model for Bremsstrahlung. Models for cyclotron radiation,
-recombination, and line radiation are still left for future work.
+Currently, TORAX has dedicated models for Bremsstrahlung and cyclotron radiation.
+Models for line radiation are left for future work.
 
 Bremsstrahlung
 ^^^^^^^^^^^^^^
 
 Uses the model from Wesson, John, and David J. Campbell. Tokamaks. Vol. 149.
 An optional correction for relativistic effects from Stott PPCF 2005 can be enabled with the flag "use_relativistic_correction".
 
-Ion Cyclotron Resonance
-^^^^^^^^^^^^^^^^^^^^^^^
+Cyclotron Radiation
+^^^^^^^^^^^^^^^^^^^
+
+Uses the total radiation power from `Albajar NF 2001 <https://doi.org/10.1088/0029-5515/41/6/301>`_
+with a deposition profile from `Artaud NF 2018 <https://doi.org/10.1088/1741-4326/aad5b1>`_.
+The Albajar model includes a parameterization of the temperature profile which in TORAX is fit by simple
+grid search for computational efficiency.
+
+Ion Cyclotron Resonance Heating
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
 Presently this source is implemented for a SPARC specific ICRH scenario.
 

diff --git a/torax/geometry/geometry.py b/torax/geometry/geometry.py
@@ -152,6 +152,7 @@ class Geometry:
   spr_cell: chex.Array
   spr_face: chex.Array
   delta_face: chex.Array
+  elongation: chex.Array
   elongation_face: chex.Array
   g0: chex.Array
   g0_face: chex.Array
@@ -273,6 +274,7 @@ class GeometryProvider:
   spr_cell: interpolated_param.InterpolatedVarSingleAxis
   spr_face: interpolated_param.InterpolatedVarSingleAxis
   delta_face: interpolated_param.InterpolatedVarSingleAxis
+  elongation: interpolated_param.InterpolatedVarSingleAxis
   elongation_face: interpolated_param.InterpolatedVarSingleAxis
   g0: interpolated_param.InterpolatedVarSingleAxis
   g0_face: interpolated_param.InterpolatedVarSingleAxis
@@ -416,6 +418,7 @@ class StandardGeometryProvider(GeometryProvider):
   jtot_face: interpolated_param.InterpolatedVarSingleAxis
   delta_upper_face: interpolated_param.InterpolatedVarSingleAxis
   delta_lower_face: interpolated_param.InterpolatedVarSingleAxis
+  elongation: interpolated_param.InterpolatedVarSingleAxis
   elongation_face: interpolated_param.InterpolatedVarSingleAxis
 
   @functools.partial(jax_utils.jit, static_argnums=0)
@@ -636,6 +639,7 @@ def build_circular_geometry(
       Rout=Rout,
       Rout_face=Rout_face,
       # Set the circular geometry-specific params.
+      elongation=elongation,
       elongation_face=elongation_face,
       volume_hires=volume_hires,
       area_hires=area_hires,
@@ -1528,6 +1532,9 @@ def build_standard_geometry(
   delta_face = 0.5 * (delta_upper_face + delta_lower_face)
 
   # elongation
+  elongation = rhon_interpolation_func(
+      rho_norm, intermediate.elongation
+  )
   elongation_face = rhon_interpolation_func(
       rho_face_norm, intermediate.elongation
   )
@@ -1628,6 +1635,7 @@ def build_standard_geometry(
       jtot_face=jtot_face,
       delta_upper_face=delta_upper_face,
       delta_lower_face=delta_lower_face,
+      elongation=elongation,
       elongation_face=elongation_face,
       volume_hires=volume_hires,
       area_hires=area_hires,

diff --git a/torax/plotting/configs/default_plot_config.py b/torax/plotting/configs/default_plot_config.py
@@ -145,19 +145,32 @@
             suppress_zero_values=True,  # Do not plot all-zero data
         ),
         plotruns_lib.PlotProperties(
-            attrs=('q_brems', 'q_rad'),
-            labels=(r'$Q_\mathrm{brems}$', r'$Q_\mathrm{rad}$'),
+            attrs=('q_brems', 'q_rad', 'q_cycl'),
+            labels=(
+                r'$Q_\mathrm{brems}$',
+                r'$Q_\mathrm{rad}$',
+                r'$Q_\mathrm{cycl}$',
+            ),
             ylabel=r'Heat sink density $[MW~m^{-3}]$',
             suppress_zero_values=True,  # Do not plot all-zero data
         ),
         plotruns_lib.PlotProperties(
             plot_type=plotruns_lib.PlotType.TIME_SERIES,
-            attrs=('p_auxiliary', 'p_ohmic', 'p_alpha', 'p_sink'),
+            attrs=(
+                'p_auxiliary',
+                'p_ohmic',
+                'p_alpha',
+                'p_brems',
+                'p_rad',
+                'p_cycl',
+            ),
             labels=(
                 r'$P_\mathrm{aux}$',
                 r'$P_\mathrm{ohm}$',
                 r'$P_\mathrm{\alpha}$',
-                r'$P_\mathrm{sink}$',
+                r'$P_\mathrm{brems}$',
+                r'$P_\mathrm{rad}$',
+                r'$P_\mathrm{cycl}$',
             ),
             ylabel=r'Total heating/sink powers $[MW]$',
             legend_fontsize=6,  # Smaller fontsize for this plot

diff --git a/torax/plotting/configs/sources_plot_config.py b/torax/plotting/configs/sources_plot_config.py
@@ -58,12 +58,21 @@
         ),
         plotruns_lib.PlotProperties(
             plot_type=plotruns_lib.PlotType.TIME_SERIES,
-            attrs=('p_auxiliary', 'p_ohmic', 'p_alpha', 'p_sink'),
+            attrs=(
+                'p_auxiliary',
+                'p_ohmic',
+                'p_alpha',
+                'p_rad',
+                'p_brems',
+                'p_cycl',
+            ),
             labels=(
                 r'$P_\mathrm{aux}$',
                 r'$P_\mathrm{ohm}$',
                 r'$P_\mathrm{\alpha}$',
-                r'$P_\mathrm{sink}$',
+                r'$P_\mathrm{rad}$',
+                r'$P_\mathrm{brems}$',
+                r'$P_\mathrm{cycl}$',
             ),
             ylabel=r'Total heating/sink powers $[MW]$',
             legend_fontsize=8,  # Smaller fontsize for this plot
@@ -99,8 +108,12 @@
             suppress_zero_values=True,  # Do not plot all-zero data
         ),
         plotruns_lib.PlotProperties(
-            attrs=('q_brems', 'q_rad'),
-            labels=(r'$Q_\mathrm{brems}$', r'$Q_\mathrm{rad}$'),
+            attrs=('q_brems', 'q_rad', 'q_cycl'),
+            labels=(
+                r'$Q_\mathrm{brems}$',
+                r'$Q_\mathrm{rad}$',
+                r'$Q_\mathrm{cycl}$',
+            ),
             ylabel=r'Heat sink density $[MW~m^{-3}]$',
             suppress_zero_values=True,  # Do not plot all-zero data
         ),

diff --git a/torax/plotting/plotruns_lib.py b/torax/plotting/plotruns_lib.py
@@ -112,6 +112,7 @@ class PlotData:
   q_alpha_e: np.ndarray  # [MW/m^3]
   q_ohmic: np.ndarray  # [MW/m^3]
   q_brems: np.ndarray  # [MW/m^3]
+  q_cycl: np.ndarray  # [MW/m^3]
   q_ei: np.ndarray  # [MW/m^3]
   q_rad: np.ndarray  # [MW/m^3]
   Q_fusion: np.ndarray  # pylint: disable=invalid-name  # Dimensionless
@@ -126,6 +127,9 @@ class PlotData:
   p_ohmic: np.ndarray  # [MW]
   p_alpha: np.ndarray  # [MW]
   p_sink: np.ndarray  # [MW]
+  p_brems: np.ndarray  # [MW]
+  p_cycl: np.ndarray  # [MW]
+  p_rad: np.ndarray  # [MW]
   t: np.ndarray  # [s]
   rho_cell_coord: np.ndarray  # Normalized toroidal flux coordinate
   rho_face_coord: np.ndarray  # Normalized toroidal flux coordinate
@@ -180,12 +184,14 @@ def _transform_data(ds: xr.Dataset):
         'fusion_heat_source_el': 1e6,  # W/m^3 to MW/m^3
         'ohmic_heat_source': 1e6,  # W/m^3 to MW/m^3
         'bremsstrahlung_heat_sink': 1e6,  # W/m^3 to MW/m^3
+        'cyclotron_radiation_heat_sink': 1e6,  # W/m^3 to MW/m^3
         'impurity_radiation_heat_sink': 1e6,  # W/m^3 to MW/m^3
         'qei_source': 1e6,  # W/m^3 to MW/m^3
         'P_ohmic': 1e6,  # W to MW
         'P_external_tot': 1e6,  # W to MW
         'P_alpha_tot': 1e6,  # W to MW
         'P_brems': 1e6,  # W to MW
+        'P_cycl': 1e6,  # W to MW
         'P_ecrh': 1e6,  # W to MW
         'P_rad': 1e6,  # W to MW
         'I_ecrh': 1e6,  # A to MA
@@ -258,6 +264,9 @@ def _transform_data(ds: xr.Dataset):
       q_brems=get_optional_data(
           core_sources_dataset, 'bremsstrahlung_heat_sink', 'cell'
       ),
+      q_cycl=get_optional_data(
+          core_sources_dataset, 'cyclotron_radiation_heat_sink', 'cell'
+      ),
       q_rad=get_optional_data(
           core_sources_dataset, 'impurity_radiation_heat_sink', 'cell'
       ),
@@ -279,7 +288,11 @@ def _transform_data(ds: xr.Dataset):
       ).to_numpy(),
       p_alpha=post_processed_outputs_dataset['P_alpha_tot'].to_numpy(),
       p_sink=post_processed_outputs_dataset['P_brems'].to_numpy()
-      + post_processed_outputs_dataset['P_rad'].to_numpy(),
+      + post_processed_outputs_dataset['P_rad'].to_numpy()
+      + post_processed_outputs_dataset['P_cycl'].to_numpy(),
+      p_brems=post_processed_outputs_dataset['P_brems'].to_numpy(),
+      p_rad=post_processed_outputs_dataset['P_rad'].to_numpy(),
+      p_cycl=post_processed_outputs_dataset['P_cycl'].to_numpy(),
       te_volume_avg=post_processed_outputs_dataset['te_volume_avg'].to_numpy(),
       ti_volume_avg=post_processed_outputs_dataset['ti_volume_avg'].to_numpy(),
       ne_volume_avg=post_processed_outputs_dataset['ne_volume_avg'].to_numpy(),

diff --git a/torax/post_processing.py b/torax/post_processing.py
@@ -39,6 +39,7 @@
 EL_HEAT_SOURCE_TRANSFORMATIONS = {
     'ohmic_heat_source': 'P_ohmic',
     'bremsstrahlung_heat_sink': 'P_brems',
+    'cyclotron_radiation_heat_sink': 'P_cycl',
     'electron_cyclotron_source': 'P_ecrh',
     'impurity_radiation_heat_sink': 'P_rad',
 }

diff --git a/torax/sources/bootstrap_current_source.py b/torax/sources/bootstrap_current_source.py
@@ -97,14 +97,6 @@ class BootstrapCurrentSource(source.Source):
   def source_name(self) -> str:
     return self.SOURCE_NAME
 
-  @property
-  def supported_modes(self) -> tuple[runtime_params_lib.Mode, ...]:
-    return (
-        runtime_params_lib.Mode.ZERO,
-        runtime_params_lib.Mode.MODEL_BASED,
-        runtime_params_lib.Mode.PRESCRIBED,
-    )
-
   @property
   def output_shape_getter(self) -> source.SourceOutputShapeFunction:
     return _default_output_shapes