Add eich profile

process/core/io/plot/summary.py

@@ -79,6 +79,7 @@
     PlasmaGeometryModelType,
     PlasmaShapeModelType,
 )
+from process.models.physics.scrape_off_layer import ScrapeOffLayer
 from process.models.superconductors import SuperconductorModel
 from process.models.tfcoil.base import (
     TFCoilShapeModel,
@@ -9090,6 +9091,50 @@ def plot_sol_power_decay_length_comparison(axis: plt.Axes, mfile: MFile, scan: i
     axis.set_facecolor("#f0f0f0")
 
 
+def plot_midplane_near_sol_radial_profile(axis: plt.Axes, mfile: MFile, scan: int):
+    """Function to plot the radial profile of the near SOL at the midplane."""
+    rmajor = mfile.get("rmajor", scan=scan)
+    rminor = mfile.get("rminor", scan=scan)
+    len_plasma_sol_power_decay = mfile.get(
+        "len_plasma_sol_eich13_power_decay", scan=scan
+    )
+    r = np.linspace(
+        (rmajor + rminor), (rmajor + rminor) + (3 * len_plasma_sol_power_decay), 100
+    )
+
+    radial_profile = (
+        ScrapeOffLayer().calculate_outboard_midplane_near_sol_radial_profile(
+            rmajor=rmajor,
+            rminor=rminor,
+            len_plasma_sol_power_decay=mfile.get(
+                "len_plasma_sol_eich13_power_decay", scan=scan
+            ),
+            pflux_plasma_outboard_sol_parallel_mw=mfile.get(
+                "pflux_plasma_outboard_sol_eich13_parallel_mw", scan=scan
+            ),
+            r=r,
+        )
+    )
+
+    axis.axvline(
+        x=rmajor + rminor + len_plasma_sol_power_decay,
+        color="k",
+        linestyle="--",
+        label=r"$\lambda_q$",
+    )
+
+    axis.set_xlim([
+        rmajor + rminor,
+        (rmajor + rminor) + (3 * len_plasma_sol_power_decay),
+    ])
+    axis.plot(r, radial_profile)
+    axis.grid()
+    axis.legend()
+    axis.set_title(r"Upstream Near SOL $q_{\parallel}$ Radial Profile")
+    axis.set_xlabel("Radial Position [m]")
+    axis.set_ylabel(r"$q_{\parallel}$ [MW/m$^2$]")
+
+
 def plot_h_threshold_comparison(axis: plt.Axes, mfile: MFile, scan: int, u_seed=None):
     """Function to plot a scatter box plot of L-H threshold power comparisons.
 
@@ -16457,6 +16502,10 @@ def _add_page(name: str | None = None):
         _add_page("plasma_compare_3").add_subplot(221), m_file, scan
     )
 
+    plot_midplane_near_sol_radial_profile(
+        _add_page("midplane_near_sol_radial_profile").add_subplot(111), m_file, scan
+    )
+
     plot_debye_length_profile(
         _add_page("microscopic_quantities").add_subplot(232), m_file, scan
     )

process/models/physics/scrape_off_layer.py

@@ -89,8 +89,7 @@ def output(self) -> None:
         )
         po.ovarre(
             self.outfile,
-            "Plasma outboard midplane Eich 2013 SOL parallel power flux "
-            "(qₗₗ,ᵤ) [MW/m²]",
+            "Plasma outboard midplane Eich 2013 SOL parallel power flux (qₗₗ,ᵤ) [MW/m²]",
             "(pflux_plasma_outboard_sol_eich13_parallel_mw)",
             self.data.physics.pflux_plasma_outboard_sol_eich13_parallel_mw,
         )
@@ -251,3 +250,61 @@ def calculate_upstream_sol_outboard_parallel_area(
             * len_plasma_sol_power_decay
             * (b_plasma_surface_poloidal_average / b_plasma_outboard_total)
         )
+
+    @staticmethod
+    def calculate_outboard_midplane_near_sol_radial_profile(
+        rmajor: float,
+        rminor: float,
+        len_plasma_sol_power_decay: float,
+        pflux_plasma_outboard_sol_parallel_mw: float,
+        r: float | np.ndarray,
+    ) -> float | np.ndarray:
+        """Calculate the outboard midplane near SOL radial profile (qₗₗ(r)) [MW/m²].
+
+        Parameters
+        ----------
+        rmajor : float
+            Major radius of the plasma (R₀) [m]
+        rminor : float
+            Minor radius of the plasma (a) [m]
+        len_plasma_sol_power_decay : float
+            Power decay length (λ_q) [m]
+        pflux_plasma_outboard_sol_parallel_mw : float
+            Parallel power flux at the outboard midplane (qₗₗ,ᵤ) [MW/m²]
+        r : float|np.ndarray
+            Radial position(s) at which to calculate the SOL profile [m]
+
+        Returns
+        -------
+        float|np.ndarray
+            Outboard midplane SOL radial profile (qₗₗ(r)) [MW/m²]
+
+        Notes
+        -----
+        - The exponential model is highly valid in the "near-SOL" (typically the first
+        few millimeters to a centimeter outside the separatrix). In this region, parallel
+        heat transport is dominated by classical electron heat conduction
+        (Spitzer-Härm conductivity), which is vastly faster than perpendicular diffusion.
+        This competition between fast parallel conduction and slow perpendicular
+        diffusion naturally produces an exponential radial profile.
+
+        - The midplane exponential assumes steady-state H-mode conditions without the
+        massive, transient convective bursts caused by ELMs, which momentarily
+        flatten the entire midplane profile.
+
+        References
+        ----------
+        [1] T. Eich et al., “Scaling of the tokamak near the scrape-off layer H-mode
+        power width and implications for ITER,” Nuclear Fusion, vol. 53, no. 9,
+        p. 093031, Aug. 2013, doi: 10.1088/0029-5515/53/9/093031.
+
+        """
+        if np.any(r < (rmajor + rminor)):
+            raise ValueError(
+                f"Radial position r={r} must be greater than or equal to the plasma "
+                f"edge (rmajor + rminor)={rmajor + rminor}."
+            )
+
+        return pflux_plasma_outboard_sol_parallel_mw * np.exp(
+            -(r - (rmajor + rminor)) / len_plasma_sol_power_decay
+        )