From c62cbfad3f7453821e298d3923f9fcbaf0d80344 Mon Sep 17 00:00:00 2001 From: mn3981 Date: Fri, 17 Oct 2025 13:59:42 +0100 Subject: [PATCH 1/3] :bug: Fix error in volume averaged thermal plasma pressure to be profile independant --- process/plasma_profiles.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/process/plasma_profiles.py b/process/plasma_profiles.py index 2953afd120..8f9e4afe8f 100644 --- a/process/plasma_profiles.py +++ b/process/plasma_profiles.py @@ -324,9 +324,11 @@ def calculate_profile_factors(self) -> None: # Shall assume that the pressure profile is parabolic. Can find volume average from # profile index and core value the same as for density and temperature physics_variables.pres_plasma_thermal_vol_avg = ( - physics_variables.pres_plasma_thermal_on_axis - / (physics_variables.alphap + 1) - ) + physics_variables.nd_plasma_electrons_vol_avg + * physics_variables.temp_plasma_electron_density_weighted_kev + + physics_variables.nd_ions_total + * physics_variables.temp_plasma_ion_density_weighted_kev + ) * constants.KILOELECTRON_VOLT # Central plasma current density (A/m^2) # Assumes a parabolic profile for the current density From 2d4bba9442c756c17e66cbe15e967762a6280ac5 Mon Sep 17 00:00:00 2001 From: mn3981 Date: Fri, 17 Oct 2025 18:07:09 +0100 Subject: [PATCH 2/3] :notes: Refine plasma pressure calculations to use density-weighted temperatures and clarify assumptions for parabolic profiles --- .../profiles/plasma_profiles.md | 48 +++++++++++-------- process/plasma_profiles.py | 22 ++++----- 2 files changed, 38 insertions(+), 32 deletions(-) diff --git a/documentation/physics-models/profiles/plasma_profiles.md b/documentation/physics-models/profiles/plasma_profiles.md index 93ca6922e9..ac0cf6ab29 100644 --- a/documentation/physics-models/profiles/plasma_profiles.md +++ b/documentation/physics-models/profiles/plasma_profiles.md @@ -271,7 +271,8 @@ The density and temperature profile runner function [`TeProfile/NeProfile.run()` Ratio of density-weighted to volume-averaged temperature factor is calculated: $$ -\texttt{pcoef} = \frac{(1+\alpha_n)(1+\alpha_T)}{(1+\alpha_T + \alpha_n)} +\texttt{f_temp_plasma_electron_density_vol_avg} =\frac{\langle T_{\text{e}} \rangle_{\text{n}}}{\underbrace{\langle T_{\text{e}} \rangle_{\text{V}}}_{\texttt{temp_plasma_electrons_vol_avg}}}\\ + = \frac{(1+\alpha_n)(1+\alpha_T)}{(1+\alpha_T + \alpha_n)} $$ -------- @@ -310,8 +311,10 @@ $\blacksquare$ The density weighted temperatures are set: $$\begin{aligned} -\texttt{temp_plasma_electron_density_weighted_kev} = \texttt{pcoef} \times T_\text{e} \\ -\texttt{temp_plasma_ion_density_weighted_kev} = \texttt{pcoef}\times T_\text{i} +\overbrace{\langle T_{\text{e}} \rangle_{\text{n}}}^{\texttt{temp_plasma_electron_density_weighted_kev}} = \\ +\texttt{f_temp_plasma_electron_density_vol_avg} \times \langle T_\text{e} \rangle \\ +\overbrace{\langle T_{\text{i}} \rangle_{\text{n}}}^{\texttt{temp_plasma_ion_density_weighted_kev}} = \\ +\texttt{f_temp_plasma_electron_density_vol_avg}\times \langle T_\text{i} \rangle \end{aligned}$$ ------ @@ -392,34 +395,43 @@ $$\begin{aligned} ##### `calculate_profile_factors()` -The central plasma pressure is calculated from the ideal gas law. +The central thermal plasma pressure is calculated from the ideal gas law. $$ -p_0 = (\texttt{nd_plasma_electron_on_axis} \times \texttt{temp_plasma_electron_on_axis_kev}\\ -+\texttt{nd_plasma_ions_on_axis}\times \texttt{temp_plasma_ion_on_axis_kev})\times (1000 \times \text{e}) +p_{\text{thermal,0}} = (\overbrace{n_{\text{e,0}}}^{\texttt{nd_plasma_electron_on_axis}} \times \overbrace{T_{\text{e,0}}}^{\texttt{temp_plasma_electron_on_axis_kev}}\\ ++\underbrace{n_{n\text{i,0}}}_{\texttt{nd_plasma_ions_on_axis}}\times \underbrace{T_{\text{i,0}}}_{\texttt{temp_plasma_ion_on_axis_kev})}\times (\times \text{keV}) $$ With the coefficients used to turn the temperature from $\text{keV}$ back to Joules. -Since we multiply the density and temperature profiles together to get the pressure we can find the pressure profile factor by adding the profile factors for temperature and density. + +We calculate the volume averaged thermal plasma pressure using the volume averaged density and the density weighted volume averaged temperature for the electrons and ions. The density weighted temperature is used as $\langle nT \rangle \neq \langle n \rangle \langle T \rangle$. This can be seen also in the [equality constraint for thermal $\beta$](../plasma_beta/plasma_beta.md#beta-consistency) $$ -\alpha_p = \alpha_n + \alpha_T +\langle p_{\text{thermal}} \rangle = \overbrace{\langle n_{\text{e}} \rangle_{\text{V}}}^{\texttt{nd_plasma_electrons_vol_avg}} \times \overbrace{\langle T_{\text{e}} \rangle_{\text{n}}}^{\texttt{temp_plasma_electron_density_weighted_kev}}\\ + + \underbrace{\langle n_{\text{i}} \rangle_{\text{V}}}_{\texttt{nd_plasma_ions_total_vol_avg}} \times \underbrace{\langle T_{\text{i}} \rangle_{\text{n}}}_{\texttt{temp_plasma_ion_density_weighted_kev}} $$ -The volume averaged pressure can then be set if we assume the pressure also has a parabolic profile. Using the standard relation used for both density and temeprature we can set the volume averaged pressure as: +???+ note "Parabolic pressure profile factor" -$$ -\langle p \rangle = \frac{p_0}{\alpha_p+1} -$$ + For a prabolic profile, since we multiply the density and temperature profiles together to get the pressure we can find the pressure profile factor by adding the profile factors for temperature and density. -???+ note "Pressure profile factor" + $$ + \alpha_p = \alpha_n + \alpha_T + $$ + + The volume averaged pressure can then be set if we assume the pressure also has a parabolic profile. Using the standard relation used for both density and temeprature we can set the volume averaged pressure as: + + $$ + \langle p \rangle = \frac{p_0}{\alpha_p+1} + $$ The calculation of $\alpha_p$ is only valid assuming a parabolic profile case. The calculatio of $p_0$ is still true as the core values are calculated independantly for each profile type. $p_0$ is NOT equal to $\langle p \rangle \times (1 + \alpha_p)$, but $p(\rho) = n(\rho)T(\rho)$ and $\langle p \rangle = \langle n \rangle$ $T_n$, where $T_n$ is the density-weighted temperature. + The on-axis plasma current density ($J_0$) is then calculated by assuming a parabolic profile and using the calculated plasma current ($I_{\text{plasma}}$) and the plasmas poloidal cross-sectional area ($A_{\text{pol}}$). $$ @@ -588,10 +600,8 @@ The density and temperature profile runner function [`TeProfile/NeProfile.run()` Density weighted temperatures are thus set as such $$ -\texttt{temp_plasma_electron_density_weighted_kev} = \\ -\frac{\int_0^1{\rho \ n(\rho) T(\rho) \ d\rho}}{\int_0^1{\rho \ n(\rho) \ d\rho}} \\ -\texttt{temp_plasma_ion_density_weighted_kev} = \\ -\texttt{temp_plasma_electron_density_weighted_kev}\times \frac{T_\text{i}}{T_\text{e}} + \overbrace{\langle T_{\text{e}} \rangle_{\text{n}}}^{\texttt{temp_plasma_electron_density_weighted_kev}} = \frac{\int_0^1{\rho \ n_{\text{e}}(\rho) T_{\text{e}}(\rho) \ d\rho}}{\int_0^1{\rho \ n_{\text{e}}(\rho) \ d\rho}} \\ +\overbrace{\langle T_{\text{i}} \rangle_{\text{n}}}^{\texttt{temp_plasma_ion_density_weighted_kev}} = \langle T_{\text{e}} \rangle_{\text{n}}\times \frac{\langle T_{\text{i}} \rangle_{\text{V}}}{\langle T_{\text{e}} \rangle_{\text{V}}} $$ The above is done numerically with [scipy.integrate.simpson](https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.simpson.html) integration. @@ -599,13 +609,13 @@ The above is done numerically with [scipy.integrate.simpson](https://docs.scipy. Set profile factor which is the ratio of density-weighted to volume-averaged temperature $$ -\texttt{pcoef} = \frac{\texttt{temp_plasma_electron_density_weighted_kev}}{T_\text{e}} +\texttt{f_temp_plasma_electron_density_vol_avg} =\frac{\langle T_{\text{e}} \rangle_{\text{n}}}{\underbrace{\langle T_{\text{e}} \rangle_{\text{V}}}_{\texttt{temp_plasma_electrons_vol_avg}}} $$ Calculate the line averaged electron density by integrating the normalised profile using the class [`integrate_profile_y()`](./plasma_profiles_abstract_class.md#calculate-the-profile-integral-value-integrate_profile_y) function $$ -\texttt{nd_plasma_electron_line} = \int_0^1{n(\rho) \ d\rho} + \overbrace{\bar{n_{\text{e}}}}^{\texttt{nd_plasma_electron_line}} = \int_0^1{n(\rho) \ d\rho} $$ A divertor variable `prn1` is set to be equal to the separatrix density over the mean density: diff --git a/process/plasma_profiles.py b/process/plasma_profiles.py index 8f9e4afe8f..c41221eccd 100644 --- a/process/plasma_profiles.py +++ b/process/plasma_profiles.py @@ -273,15 +273,11 @@ def calculate_profile_factors(self) -> None: # Central pressure (Pa), from ideal gas law : p = nkT physics_variables.pres_plasma_thermal_on_axis = ( - ( - physics_variables.nd_plasma_electron_on_axis - * physics_variables.temp_plasma_electron_on_axis_kev - + physics_variables.nd_plasma_ions_on_axis - * physics_variables.temp_plasma_ion_on_axis_kev - ) - * 1.0e3 - * constants.ELECTRON_CHARGE - ) + physics_variables.nd_plasma_electron_on_axis + * physics_variables.temp_plasma_electron_on_axis_kev + + physics_variables.nd_plasma_ions_on_axis + * physics_variables.temp_plasma_ion_on_axis_kev + ) * constants.KILOELECTRON_VOLT # Electron pressure profile (Pa) physics_variables.pres_plasma_electron_profile = self.neprofile.profile_y * ( @@ -314,19 +310,19 @@ def calculate_profile_factors(self) -> None: * physics_variables.f_temp_plasma_ion_electron ) - # Pressure profile index (N.B. no pedestal effects included here) + # Pressure profile index (only true for a parabolic profile) # N.B. pres_plasma_thermal_on_axis is NOT equal to

* (1 + alphap), but p(rho) = n(rho)*T(rho) # and

= .T_n where <...> denotes volume-averages and T_n is the # density-weighted temperature physics_variables.alphap = physics_variables.alphan + physics_variables.alphat - # Shall assume that the pressure profile is parabolic. Can find volume average from - # profile index and core value the same as for density and temperature + # Calculate the volume averaged plasma thermal pressure from the density-weighted temperatures + # Density-weighted temperatures are used as != * physics_variables.pres_plasma_thermal_vol_avg = ( physics_variables.nd_plasma_electrons_vol_avg * physics_variables.temp_plasma_electron_density_weighted_kev - + physics_variables.nd_ions_total + + physics_variables.nd_plasma_ions_total_vol_avg * physics_variables.temp_plasma_ion_density_weighted_kev ) * constants.KILOELECTRON_VOLT From 5369e1046a0fbce523d62c3958d452f2fafc3542 Mon Sep 17 00:00:00 2001 From: mn3981 Date: Tue, 21 Oct 2025 13:43:54 +0100 Subject: [PATCH 3/3] Add volume-average thermal pressure line and annotations to plasma pressure profiles plot --- process/io/plot_proc.py | 28 ++++++++++++++++++++++++++++ 1 file changed, 28 insertions(+) diff --git a/process/io/plot_proc.py b/process/io/plot_proc.py index 55756ab44d..305e10695e 100644 --- a/process/io/plot_proc.py +++ b/process/io/plot_proc.py @@ -10567,6 +10567,18 @@ def plot_plasma_pressure_profiles(axis, mfile_data, scan): color="green", label="Total", ) + + # Plot horizontal line for volume-average thermal pressure (converted to kPa) + p_vol_kpa = mfile_data.data["pres_plasma_thermal_vol_avg"].get_scan(scan) / 1000.0 + axis.axhline( + p_vol_kpa, + color="black", + linestyle="--", + linewidth=1.2, + label="Volume avg", + zorder=5, + ) + axis.set_xlabel("$\\rho$ [r/a]") axis.set_ylabel("Thermal Pressure [kPa]") axis.minorticks_on() @@ -10577,6 +10589,22 @@ def plot_plasma_pressure_profiles(axis, mfile_data, scan): axis.set_ylim(bottom=0) axis.legend() + textstr_pressure = "\n".join(( + rf"$p_0$: {mfile_data.data['pres_plasma_thermal_on_axis'].get_scan(scan) / 1000:,.3f} kPa", + rf"$\langle p_{{\text{{total}}}} \rangle_\text{{V}}$: {mfile_data.data['pres_plasma_thermal_vol_avg'].get_scan(scan) / 1000:,.3f} kPa", + )) + + axis.text( + 0.5, + 1.2, + textstr_pressure, + transform=axis.transAxes, + fontsize=9, + verticalalignment="top", + horizontalalignment="center", + bbox={"boxstyle": "round", "facecolor": "wheat", "alpha": 0.5}, + ) + def plot_plasma_pressure_gradient_profiles(axis, mfile_data, scan): # Get the plasma pressure profiles