From d5b3db24051309491ccd62a994ed31ae8164e5c5 Mon Sep 17 00:00:00 2001 From: ym1906 Date: Fri, 14 Feb 2025 13:28:28 +0000 Subject: [PATCH 1/4] Updated tprofile to teprofile and nprofile to neprofile. --- process/impurity_radiation.py | 42 +++++++++++++------------- process/io/plot_proc.py | 8 ++--- process/plasma_profiles.py | 16 +++++----- process/profiles.py | 16 +++++----- tests/unit/test_plasma_profiles.py | 48 +++++++++++++++--------------- 5 files changed, 65 insertions(+), 65 deletions(-) diff --git a/process/impurity_radiation.py b/process/impurity_radiation.py index 370bbeb284..beb962f55e 100644 --- a/process/impurity_radiation.py +++ b/process/impurity_radiation.py @@ -340,19 +340,19 @@ def fradcore(rho, coreradius, coreradiationfraction): return fradcore -def zav_of_te(imp_element_index, tprofile): +def zav_of_te(imp_element_index, teprofile): """Calculates electron temperature dependent average atomic number :param imp_element_index: Impurity element index :type imp_element_index: int - :param tprofile: temperature profile - :type tprofile: numpy.array + :param teprofile: temperature profile + :type teprofile: numpy.array :return: zav_of_te - electron temperature dependent average atomic number :rtype: numpy.array """ - # less_than_imp_temp_mask = tprofile values less than impurity temperature. greater_than_imp_temp_mask = tprofile values higher than impurity temperature. + # less_than_imp_temp_mask = teprofile values less than impurity temperature. greater_than_imp_temp_mask = teprofile values higher than impurity temperature. bins = impurity_radiation_module.impurity_arr_temp_kev[imp_element_index] - indices = np.digitize(tprofile, bins) + indices = np.digitize(teprofile, bins) indices[indices >= bins.shape[0]] = bins.shape[0] - 1 indices[indices < 0] = 0 yi = impurity_radiation_module.impurity_arr_zav[imp_element_index, indices - 1] @@ -367,16 +367,16 @@ def zav_of_te(imp_element_index, tprofile): ) - xi ) - zav_of_te = yi + c * (np.log(tprofile) - xi) + zav_of_te = yi + c * (np.log(teprofile) - xi) less_than_imp_temp_mask = ( - tprofile + teprofile <= impurity_radiation_module.impurity_arr_temp_kev[imp_element_index, 0] ) zav_of_te[less_than_imp_temp_mask] = impurity_radiation_module.impurity_arr_zav[ imp_element_index, 0 ] greater_than_imp_temp_mask = ( - tprofile + teprofile >= impurity_radiation_module.impurity_arr_temp_kev[ imp_element_index, (impurity_radiation_module.impurity_arr_len_tab[imp_element_index]) - 1, @@ -390,21 +390,21 @@ def zav_of_te(imp_element_index, tprofile): return zav_of_te -def pimpden(imp_element_index, nprofile, tprofile): +def pimpden(imp_element_index, neprofile, teprofile): """Calculates the impurity radiation density (W/m3) :param imp_element_index: Impurity element index :type imp_element_index: int - :param nprofile: density profile - :type nprofile: numpy.array - :param tprofile: temperature profile - :type tprofile: numpy.array + :param neprofile: density profile + :type neprofile: numpy.array + :param teprofile: temperature profile + :type teprofile: numpy.array :return: pimpden - total impurity radiation density (W/m3) :rtype: numpy.array """ - # less_than_imp_temp_mask = tprofile values less than impurity temperature. greater_than_imp_temp_mask = tprofile values higher than impurity temperature. + # less_than_imp_temp_mask = teprofile values less than impurity temperature. greater_than_imp_temp_mask = teprofile values higher than impurity temperature. bins = impurity_radiation_module.impurity_arr_temp_kev[imp_element_index] - indices = np.digitize(tprofile, bins) + indices = np.digitize(teprofile, bins) indices[indices >= bins.shape[0]] = bins.shape[0] - 1 indices[indices < 0] = 0 @@ -425,17 +425,17 @@ def pimpden(imp_element_index, nprofile, tprofile): ) - xi ) - pimpden = np.exp(yi + c * (np.log(tprofile) - xi)) + pimpden = np.exp(yi + c * (np.log(teprofile) - xi)) pimpden = ( impurity_radiation_module.impurity_arr_frac[imp_element_index] - * nprofile - * nprofile + * neprofile + * neprofile * pimpden ) less_than_imp_temp_mask = ( - tprofile + teprofile <= impurity_radiation_module.impurity_arr_temp_kev[imp_element_index, 0] ) pimpden[less_than_imp_temp_mask] = impurity_radiation_module.impurity_arr_lz_wm3[ @@ -443,11 +443,11 @@ def pimpden(imp_element_index, nprofile, tprofile): ] # if not impurity_radiation_module.toolow: # Only print warning once during a run # impurity_radiation_module.toolow = True - # error_handling.fdiags[0] = tprofile + # error_handling.fdiags[0] = teprofile # error_handling.report_error(35) greater_than_imp_temp_mask = ( - tprofile + teprofile >= impurity_radiation_module.impurity_arr_temp_kev[ imp_element_index, impurity_radiation_module.impurity_arr_len_tab[imp_element_index] - 1, diff --git a/process/io/plot_proc.py b/process/io/plot_proc.py index fae3962962..f9f1cb6b3f 100644 --- a/process/io/plot_proc.py +++ b/process/io/plot_proc.py @@ -747,7 +747,7 @@ def arc_fill(axis, r1, r2, color="pink"): axis.add_patch(patch) -def plot_nprofile(prof, demo_ranges): +def plot_neprofile(prof, demo_ranges): """Function to plot density profile Arguments: prof --> axis object to add plot to @@ -882,7 +882,7 @@ def plot_jprofile(prof): ) -def plot_tprofile(prof, demo_ranges): +def plot_teprofile(prof, demo_ranges): """Function to plot temperature profile Arguments: prof --> axis object to add plot to @@ -3384,12 +3384,12 @@ def main_plot( # Plot density profiles plot_4 = fig2.add_subplot(231) # , aspect= 0.05) plot_4.set_position([0.075, 0.55, 0.25, 0.4]) - plot_nprofile(plot_4, demo_ranges) + plot_neprofile(plot_4, demo_ranges) # Plot temperature profiles plot_5 = fig2.add_subplot(232) plot_5.set_position([0.375, 0.55, 0.25, 0.4]) - plot_tprofile(plot_5, demo_ranges) + plot_teprofile(plot_5, demo_ranges) # Plot impurity profiles plot_8 = fig2.add_subplot(233) diff --git a/process/plasma_profiles.py b/process/plasma_profiles.py index a298c9f31d..698de11d48 100644 --- a/process/plasma_profiles.py +++ b/process/plasma_profiles.py @@ -15,14 +15,14 @@ class PlasmaProfile: """ - Plasma profile class. Initiates the density and temperature profiles and + Plasma profile class. Initiates the electron density and electron temperature profiles and handles the required physics variables. Attributes: profile_size (int): The size of the plasma profile. outfile (str): The output file path. - neprofile (NProfile): An instance of the NProfile class. - teprofile (TProfile): An instance of the TProfile class. + neprofile (NeProfile): An instance of the NeProfile class. + teprofile (TeProfile): An instance of the TeProfile class. Methods: run(): Subroutine to execute PlasmaProfile functions. @@ -40,14 +40,14 @@ def __init__(self) -> None: Args: profile_size (int): The size of the plasma profile. outfile (str): The output file path. - neprofile (NProfile): An instance of the NProfile class. - teprofile (TProfile): An instance of the TProfile class. + neprofile (NeProfile): An instance of the NeProfile class. + teprofile (TeProfile): An instance of the TeProfile class. """ # Default profile_size = 501, but it's possible to experiment with this value. self.profile_size = 501 self.outfile = constants.nout - self.neprofile = profiles.NProfile(self.profile_size) - self.teprofile = profiles.TProfile(self.profile_size) + self.neprofile = profiles.NeProfile(self.profile_size) + self.teprofile = profiles.TeProfile(self.profile_size) def run(self) -> None: """ @@ -186,7 +186,7 @@ def pedestal_parameterisation(self) -> None: Returns: None """ - # Run TProfile and NProfile class methods: + # Run TeProfile and NeProfile class methods: # Re-caluclate core and profile values self.teprofile.run() diff --git a/process/profiles.py b/process/profiles.py index 7e6813d5fe..758ee1930f 100644 --- a/process/profiles.py +++ b/process/profiles.py @@ -90,7 +90,7 @@ def integrate_profile_y(self) -> None: ) -class NProfile(Profile): +class NeProfile(Profile): """Electron density profile class. Contains a function to calculate the electron density profile and store the data. @@ -98,14 +98,14 @@ class NProfile(Profile): Inherits attributes from the base class `Profile`. Methods: - run(): Subroutine which calls functions and stores nprofile data. + run(): Subroutine which calls functions and stores neprofile data. calculate_profile_y(rho, rhopedn, n0, nped, nsep, alphan): Calculates the density at each normalised minor radius position. ncore(rhopedn, nped, nsep, nav, alphan): Calculates the core density of a pedestalised profile. set_physics_variables(): Calculates and sets physics variables required for the profile. """ def run(self) -> None: - """Subroutine which calls profile functions and stores nprofile data.""" + """Subroutine which calls profile functions and stores neprofile data.""" self.normalise_profile_x() self.calculate_profile_dx() self.set_physics_variables() @@ -130,7 +130,7 @@ def calculate_profile_y( ) -> None: """ This routine calculates the density at each normalised minor radius position - rho for a HELIOS-type density pedestal profile (nprofile). + rho for a HELIOS-type density pedestal profile (neprofile). Authors: R Kemp, CCFE, Culham Science Centre @@ -237,11 +237,11 @@ def set_physics_variables(self) -> None: ) -class TProfile(Profile): - """Temperature profile class. Contains a function to calculate the temperature profile and store the data.""" +class TeProfile(Profile): + """Electron temperature profile class. Contains a function to calculate the temperature profile and store the data.""" def run(self) -> None: - """Subroutine to initialise nprofile and execute calculations.""" + """Subroutine to initialise neprofile and execute calculations.""" self.normalise_profile_x() self.calculate_profile_dx() self.set_physics_variables() @@ -267,7 +267,7 @@ def calculate_profile_y( tbeta: float, ) -> None: """ - Calculates the temperature at a normalised minor radius position rho for a pedestalised profile (tprofile). + Calculates the temperature at a normalised minor radius position rho for a pedestalised profile (teprofile). If ipedestal = 0 the original parabolic profile form is used instead. References: Jean, J. (2011). HELIOS: A Zero-Dimensional Tool for Next Step and Reactor Studies. Fusion Science and Technology, 59(2), 308-349. https://doi.org/10.13182/FST11-A11650 diff --git a/tests/unit/test_plasma_profiles.py b/tests/unit/test_plasma_profiles.py index 4d3ab9f291..27aa89cbb3 100644 --- a/tests/unit/test_plasma_profiles.py +++ b/tests/unit/test_plasma_profiles.py @@ -5,7 +5,7 @@ from process.fortran import divertor_variables, physics_variables from process.plasma_profiles import PlasmaProfile -from process.profiles import NProfile, TProfile +from process.profiles import NeProfile, TeProfile class ProfileParam(NamedTuple): @@ -13,27 +13,27 @@ class ProfileParam(NamedTuple): expected_profile: np.array -class NProfileParam(NamedTuple): +class NeProfileParam(NamedTuple): nesep: float = 0.0 ipedestal: float = 0.0 ne0: float = 0.0 neped: float = 0.0 rhopedn: float = 0.0 alphan: float = 0.0 - expected_nprofile: list | None = None + expected_neprofile: list | None = None @pytest.mark.parametrize( - "nprofileparam", + "neprofileparam", ( - NProfileParam( + NeProfileParam( nesep=3.6421334486704804e19, ipedestal=1, ne0=0.0, neped=6.1916268627398164e19, rhopedn=0.94000000000000006, alphan=1, - expected_nprofile=[ + expected_neprofile=[ 1.12500000e20, 1.12275191e20, 1.11587869e20, @@ -49,45 +49,45 @@ class NProfileParam(NamedTuple): ), ids=["baseline_2018"], ) -def test_nprofile(nprofileparam: ProfileParam, monkeypatch): - nprofile = NProfile(10) - nprofile.run() - assert nprofile.profile_y == pytest.approx(nprofileparam.expected_nprofile) +def test_neprofile(neprofileparam: ProfileParam, monkeypatch): + neprofile = NeProfile(10) + neprofile.run() + assert neprofile.profile_y == pytest.approx(neprofileparam.expected_neprofile) def test_ncore(): - nprofile = NProfile(10) + neprofile = NeProfile(10) rhopedn = 0.94 nped = 5.8300851381352219e19 nsep = 3.4294618459618943e19 nav = 7.4321e19 alphan = 1.0 - assert nprofile.ncore(rhopedn, nped, nsep, nav, alphan) == pytest.approx( + assert neprofile.ncore(rhopedn, nped, nsep, nav, alphan) == pytest.approx( 9.7756974320342041e19 ) -class TProfileParam(NamedTuple): +class TeProfileParam(NamedTuple): rhopedt: float = 0.0 tbeta: float = 0.0 tesep: float = 0.0 ipedestal: float = 0.0 alphat: float = 0.0 teped: float = 0.0 - expected_tprofile: Any = np.array + expected_teprofile: Any = np.array @pytest.mark.parametrize( - "tprofileparam", + "teprofileparam", ( - TProfileParam( + TeProfileParam( rhopedt=0.94000000000000006, tbeta=2, tesep=0.10000000000000001, ipedestal=1, alphat=1.45, teped=5.5, - expected_tprofile=[ + expected_teprofile=[ 18.85, 18.739472621498333, 18.403679368327712, @@ -103,15 +103,15 @@ class TProfileParam(NamedTuple): ), ids=["baseline_2018"], ) -def test_tprofile(tprofileparam: ProfileParam, monkeypatch): - monkeypatch.setattr(physics_variables, "ipedestal", tprofileparam.ipedestal) - tprofile = TProfile(10) - tprofile.run() - assert tprofile.profile_y == pytest.approx(tprofileparam.expected_tprofile) +def test_teprofile(teprofileparam: ProfileParam, monkeypatch): + monkeypatch.setattr(physics_variables, "ipedestal", teprofileparam.ipedestal) + teprofile = TeProfile(10) + teprofile.run() + assert teprofile.profile_y == pytest.approx(teprofileparam.expected_teprofile) def test_tcore(): - tprofile = TProfile(10) + teprofile = TeProfile(10) rhopedt = 0.94 tped = 3.7775374842470044 tsep = 0.1 @@ -119,7 +119,7 @@ def test_tcore(): alphat = 1.45 tbeta = 2.0 - assert tprofile.tcore(rhopedt, tped, tsep, tav, alphat, tbeta) == pytest.approx( + assert teprofile.tcore(rhopedt, tped, tsep, tav, alphat, tbeta) == pytest.approx( 28.09093632260765 ) From 43c4a695da5d0ccfd19867341532e18f89736562 Mon Sep 17 00:00:00 2001 From: ym1906 Date: Fri, 14 Feb 2025 14:30:18 +0000 Subject: [PATCH 2/4] "Updated docs with teprofile and neprofile change" --- .../RF/culham_electron_cyclotron.md | 14 ++-- .../RF/culham_lower_hybrid.md | 9 +-- .../profiles/plasma_density_profile.md | 41 +++++------- .../profiles/plasma_profiles.md | 67 +++++++------------ .../profiles/plasma_temperature_profile.md | 38 +++++------ 5 files changed, 67 insertions(+), 102 deletions(-) diff --git a/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_electron_cyclotron.md b/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_electron_cyclotron.md index 204f6fd053..b0246bff67 100644 --- a/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_electron_cyclotron.md +++ b/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_electron_cyclotron.md @@ -1,13 +1,11 @@ # Culham Electron Cyclotron Model | `culecd()` - -- `iefrf/iefrffix` = 7 +- `iefrf/iefrffix` = 7 This routine calculates the current drive parameters for a electron cyclotron system, based on the AEA FUS 172 model[^1] - -1. Local electron temperature $(\mathtt{tlocal})$ is calculated using the `tprofile` method -2. Local electron density $(\mathtt{dlocal})$ is calculated using the `nprofile` method +1. Local electron temperature $(\mathtt{tlocal})$ is calculated using the `teprofile` method +2. Local electron density $(\mathtt{dlocal})$ is calculated using the `neprofile` method 3. Calculate the inverse aspect ratio `epsloc`. @@ -47,7 +45,7 @@ Uses the `eccdef` model found [here](ec_overview.md) $$ \mathtt{ecgam} = 0.25(\mathtt{ecgam1} + \mathtt{ecgam2} +\mathtt{ecgam3} + \mathtt{ecgam4}) - $$ + $$ 7. Calculate the current drive efficiency by dividing `ecgam` by `(dlocal * physics_variables.rmajor)`. @@ -57,8 +55,6 @@ $$ Note: The `eccdef` method is called to calculate the current drive efficiency at each poloidal angle. - - [^1]: Hender, T.C., Bevir, M.K., Cox, M., Hastie, R.J., Knight, P.J., Lashmore-Davies, C.N., Lloyd, B., Maddison, G.P., Morris, A.W., O’Brien, M.R. and Turner, M.F., 1992. *"Physics assessment for the European reactor study."* AEA FUS, 172. -[^2]: Wesson, J. and Campbell, D.J., 2011. *"Tokamaks"* (Vol. 149). Oxford university press. \ No newline at end of file +[^2]: Wesson, J. and Campbell, D.J., 2011. *"Tokamaks"* (Vol. 149). Oxford university press. diff --git a/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_lower_hybrid.md b/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_lower_hybrid.md index 60325886d5..e9be2ae5c9 100644 --- a/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_lower_hybrid.md +++ b/documentation/proc-pages/eng-models/heating_and_current_drive/RF/culham_lower_hybrid.md @@ -2,18 +2,15 @@ - `iefrf/iefrffix` = 6 - - This routine calculates the current drive parameters for a lower hybrid system, based on the AEA FUS 172 model. AEA FUS 251: A User's Guide to the PROCESS Systems Code AEA FUS 172: Physics Assessment for the European Reactor Study[^1] - 1. Call the [`lhrad()`](#lower-hybrid-wave-absorption-radius--lhrad) method to calculate the lower hybrid wave absorption radius, `rratio`. 2. Calculate the penetration radius, `rpenet`, by multiplying `rratio` with the minor radius of the plasma. -3. Calculate the local density, `dlocal`, using the `nprofile()` function from the `profiles_module` module. This function takes into account various plasma parameters such as the density profile, electron density at the edge, pedestal density, separatrix density, and the value of the parameter `alphan`. -4. Similarly, calculate the local temperature, `tlocal`, using the `tprofile()` function from the `profiles_module` module. This function considers parameters such as the temperature profile, electron temperature at the edge, pedestal temperature, separatrix temperature, `alphat`, and `tbeta`. +3. Calculate the local density, `dlocal`, using the `neprofile()` function from the `profiles_module` module. This function takes into account various plasma parameters such as the density profile, electron density at the edge, pedestal density, separatrix density, and the value of the parameter `alphan`. +4. Similarly, calculate the local temperature, `tlocal`, using the `teprofile()` function from the `profiles_module` module. This function considers parameters such as the temperature profile, electron temperature at the edge, pedestal temperature, separatrix temperature, `alphat`, and `tbeta`. 5. Calculate the local toroidal magnetic field, `blocal`, using the formula `bt * rmajor / (rmajor - rpenet)`. Here, `bt` is the toroidal magnetic field at the magnetic axis, and `rmajor` is the major radius of the plasma. 6. Calculate the parallel refractive index, `nplacc`, which is needed for plasma access. It uses the local density `dlocal` and the local magnetic field `blocal` to calculate a fraction `frac`. `nplacc` is then obtained by adding `frac` to the square root of `1.0 + frac * frac`. 7. Calculate the local inverse aspect ratio, `epslh`, by dividing `rpenet` by `rmajor`. @@ -46,5 +43,3 @@ Iterate over the following steps to find the minor radius ratio: - Update $\mathtt{rat0}$ with the new approximation, $\mathtt{rat1}$. 3. If the loop completes all 100 iterations without finding a satisfactory solution, report an error and set $\mathtt{rat0}$ to 0.8. 4. Return the final value of $\mathtt{rat0}$. - -[^1]: T. C. Hender, M. K. Bevir, M. Cox, R. J. Hastie, P. J. Knight, C. N. Lashmore-Davies, B. Lloyd, G. P. Maddison, A. W. Morris, M. R. O'Brien, M.F. Turner abd H. R. Wilson, *"Physics Assessment for the European Reactor Study"*, AEA Fusion Report AEA FUS 172 (1992) \ No newline at end of file diff --git a/documentation/proc-pages/physics-models/profiles/plasma_density_profile.md b/documentation/proc-pages/physics-models/profiles/plasma_density_profile.md index 41c904211a..21ce9c6070 100644 --- a/documentation/proc-pages/physics-models/profiles/plasma_density_profile.md +++ b/documentation/proc-pages/physics-models/profiles/plasma_density_profile.md @@ -1,4 +1,4 @@ -# Density Profile | `NProfile(Profile)` +# Density Profile | `NeProfile(Profile)` The desnity profile class is organised around a central runner function that is called each time the plasma is parameterised by the parent [`PlasmaProfile()`](./plasma_profiles.md) class. It is called by [`pedestal_parameterisation()`](plasma_profiles.md#pedestal_parameterisation) and [`parabolic parameterisation()`](./plasma_profiles.md#parabolic_paramterisation). The sequence of the runner function can be seen below along with explanation of the following calculations. @@ -8,10 +8,9 @@ The desnity profile class is organised around a central runner function that is 2. The steps between the normalized points is then done by [`calculate_profile_dx()`](./plasma_profiles_abstract_class.md#calculate-the-profile-steps-in-x--calculate_profile_dx) which divides the max x-dimension by the number of points. -3. The core electron and ion temperatures are claculated via [`set_physics_variables()`]() +3. The core electron and ion temperatures are claculated via [`set_physics_variables()`]() - - ### Calculate core values | `set_physics_variables()` + ### Calculate core values | `set_physics_variables()` The core electron density is calculated using the [`ncore`](plasma_density_profile.md#electron-core-density-of-a-pedestalised-profile--ncore) method. The core ion density is then set from $n_{\text{i}}$ (`nd_ions_total`) which is the total ion density such as: @@ -19,10 +18,10 @@ The desnity profile class is organised around a central runner function that is $$ n_{\text{i0}} = \left(\frac{n_\text{i}}{n_\text{e}}\right)n_{\text{e0}} $$ - - #### Electron core density of a pedestalised profile | `ncore()` - This function calculates the core electron density for a pedestalsied profile (`ipedestal == 1`). It takes in values of + #### Electron core density of a pedestalised profile | `ncore()` + + This function calculates the core electron density for a pedestalsied profile (`ipedestal == 1`). It takes in values of | Profile parameter / Input | Density | |----------------------------------|-----------| @@ -32,7 +31,6 @@ The desnity profile class is organised around a central runner function that is | Average density | `dene`, $\langle n \rangle$ | | Profile index/ peaking parameter | `alphan`, $\alpha_n$ | - $$ n_0 = \frac{1}{3\rho_{\text{ped,n}}^2}\left[3\langle n_{\text{e}} \rangle (1+\alpha_n) + n_{\text{sep}} (1+\alpha_n) \left(-2 + \rho_{\text{ped,n}} + \rho_{\text{ped,n}}^2\right) \\ @@ -42,7 +40,7 @@ The desnity profile class is organised around a central runner function that is If `ncore` is returned as being less than 0, it is forced into a state of `ncore = 1E-6` in order to help convergence. This will also give a warning to the user to raise the lower bound of the average electron density `dene`. - ##### Derivation + ##### Derivation We calculate the volume integrated profile and then divide by the volume of integration to get the volume average density $\langle n_{\text{e}} \rangle$. If we assume the plasma to be a torus of circular cross-section then we can use spherical coordinates. We can simplify the problem by representing the torus as a cylinder of height equal to the circumference of the torus equal to $2\pi R$ where $R$ is the major radius of the torus, and $a$ is the plasma minor radius in the poloidal plane. @@ -71,26 +69,26 @@ The desnity profile class is organised around a central runner function that is \frac{\rho^2}{\rho_{\text{ped},n}^2}\right)^{\alpha_n}\right) \ d\rho \\ +\int^1_{\rho_{\text{ped,n}}} \rho\left(n_{\text{sep}} + (n_{\text{ped}} - n_{\text{sep}})\left( \frac{1- \rho}{1-\rho_{\text{ped},n}}\right)\right)\right] \ d\rho $$ - + Integrating each part within its bounds: - + $$ 4\pi^2R\left[ \frac{\left(n_{\text{ped}} {\alpha}_{n} + n_{0}\right) {\rho}_{\text{ped,n}}^{2}}{2{\alpha}_{n} + 2} \\ - +\frac{\left(1-{\rho}_{\text{ped,n}}\right) \left(\left(n_{\text{sep}} + 2n_{\text{ped}}\right) {\rho}_{\text{ped,n}} + 2n_{\text{sep}} + n_{\text{ped}}\right)}{6}\right] + +\frac{\left(1-{\rho}_{\text{ped,n}}\right) \left(\left(n_{\text{sep}} + 2n_{\text{ped}}\right) {\rho}_{\text{ped,n}} + 2n_{\text{sep}} + n_{\text{ped}}\right)}{6}\right] $$ In the form of volume average density where the volume integrated density function has to be divided by the volume of the cylinder / torus, within the volume bounded by that pedestal position we get: $$ \langle n_{\text{e}} \rangle = 4\pi^2R\left[ \frac{\frac{\left(n_{\text{ped}} {\alpha}_{n} + n_{0}\right) {\rho}_{\text{ped,n}}^{2}}{2{\alpha}_{n} + 2} - +\frac{\left(1-{\rho}_{\text{ped,n}}\right) \left(\left(n_{\text{sep}} + 2n_{\text{ped}}\right) {\rho}_{\text{ped,n}} + 2n_{\text{sep}} + n_{\text{ped}}\right)}{6}}{2\pi^2 R \rho^2}\right] + +\frac{\left(1-{\rho}_{\text{ped,n}}\right) \left(\left(n_{\text{sep}} + 2n_{\text{ped}}\right) {\rho}_{\text{ped,n}} + 2n_{\text{sep}} + n_{\text{ped}}\right)}{6}}{2\pi^2 R \rho^2}\right] $$ In this case, the value of $\rho$ is equal to 1 as we integrated over the full profile. $$ \langle n_{\text{e}} \rangle = 2\left[\frac{\left(n_{\text{ped}} {\alpha}_{n} + n_{0}\right) {\rho}_{\text{ped,n}}^{2}}{2{\alpha}_{n} + 2} \\ - +\frac{\left(1-{\rho}_{\text{ped,n}}\right) \left(\left(n_{\text{sep}} + 2n_{\text{ped}}\right) {\rho}_{\text{ped,n}} + 2n_{\text{sep}} + n_{\text{ped}}\right)}{6}\right] + +\frac{\left(1-{\rho}_{\text{ped,n}}\right) \left(\left(n_{\text{sep}} + 2n_{\text{ped}}\right) {\rho}_{\text{ped,n}} + 2n_{\text{sep}} + n_{\text{ped}}\right)}{6}\right] $$ $$ @@ -105,16 +103,16 @@ The desnity profile class is organised around a central runner function that is - n_{\text{ped}}\left( (1 + \alpha_n)(1+ \rho_{\text{ped,n}}) + (\alpha_n -2) \rho_{\text{ped,n}}^2\right)\right] $$ - + $\blacksquare$ ------ 4. The y profile is then calculated using [`calculate_profile_y()`](plasma_density_profile.md#calculate-desnity-at-each-radius-position-calculate_profile_y). This routine calculates the density at each normalised minor radius position $\rho$ for a HELIOS-type density pedestal profile[^1] - ### Calculate desnity at each radius position | `calculate_profile_y()` + ### Calculate desnity at each radius position | `calculate_profile_y()` - A table of the input variables can be found below + A table of the input variables can be found below | Profile parameter / Input | Density | |----------------------------------|-----------| @@ -128,7 +126,7 @@ The desnity profile class is organised around a central runner function that is If `ipedestal == 0` then the original parabolic profile form is used $$ - n(\rho) = n_0(1 - \rho^2)^{\alpha_n} + n(\rho) = n_0(1 - \rho^2)^{\alpha_n} $$ The central density ($n_0$) is then checked to make sure it is not less than the pedestal density, $n_{\text{ped}}$. @@ -136,10 +134,8 @@ The desnity profile class is organised around a central runner function that is Values of the profile density are then assigned based on the density function below across bounds from 0 to `rhopedn` and `rhopedn` to 1. - - $$\begin{aligned} - \mbox{Density:} \ n(\rho) = \left\{ + \mbox{Density:} \ n(\rho) = \left\{ \begin{aligned} & n_{\text{ped}} + (n_0 - n_{\text{ped}}) \left( 1 - \frac{\rho^2}{\rho_{\text{ped,n}}^2}\right)^{\alpha_n} @@ -149,8 +145,7 @@ The desnity profile class is organised around a central runner function that is \end{aligned} \right. \end{aligned}$$ - 5. Profile is then integrated with `integrate_profile_y()` using Simpsons integration from the profile abstract base class -[^1]: Jean, J. (2011). *HELIOS: A Zero-Dimensional Tool for Next Step and Reactor Studies*. Fusion Science and Technology, 59(2), 308–349. https://doi.org/10.13182/FST11-A11650 \ No newline at end of file +[^1]: Jean, J. (2011). *HELIOS: A Zero-Dimensional Tool for Next Step and Reactor Studies*. Fusion Science and Technology, 59(2), 308–349. diff --git a/documentation/proc-pages/physics-models/profiles/plasma_profiles.md b/documentation/proc-pages/physics-models/profiles/plasma_profiles.md index 82dd4b638e..14e683d016 100644 --- a/documentation/proc-pages/physics-models/profiles/plasma_profiles.md +++ b/documentation/proc-pages/physics-models/profiles/plasma_profiles.md @@ -4,14 +4,13 @@ In `PROCESS` the density, temperature and current profiles of the plasma for electrons and ions can take two forms depending on the switch value for `ipedestal`. Either without a [pedestal](http://fusionwiki.ciemat.es/wiki/Pedestal), `ipedestal == 0` or with a pedestal `ipedestal == 1`. `ipedestal == 0` is better suited for modelling L-mode plasmas, while `ipedestal == 1` is better suited for modelling [H-mode](https://en.wikipedia.org/wiki/High-confinement_mode) plasmas. -The files responsible for calculating and storing the profiles are `plasma_profiles.py` and `profiles.py`. A central plasma profile object is created from the [`PlasmaProfile`](plasma_profiles.md#plasma-profile-class-plasmaprofile) class that contains attributes for the plasma density and temperature. The density and temperature profiles are in themselves objects of the [`Profile`](./plasma_profiles_abstract_class.md) abstract base class. [`Profile`](./plasma_profiles_abstract_class.md), [`NProfile`](plasma_density_profile.md) and [`TProfile`](./plasma_temperature_profile.md) are all defined in `profiles.py`. [`PlasmaProfile`](plasma_profiles.md#plasma-profile-class-plasmaprofile) is exclusively in `plasma_profiles.py` +The files responsible for calculating and storing the profiles are `plasma_profiles.py` and `profiles.py`. A central plasma profile object is created from the [`PlasmaProfile`](plasma_profiles.md#plasma-profile-class-plasmaprofile) class that contains attributes for the plasma density and temperature. The density and temperature profiles are in themselves objects of the [`Profile`](./plasma_profiles_abstract_class.md) abstract base class. [`Profile`](./plasma_profiles_abstract_class.md), [`NeProfile`](plasma_density_profile.md) and [`TeProfile`](./plasma_temperature_profile.md) are all defined in `profiles.py`. [`PlasmaProfile`](plasma_profiles.md#plasma-profile-class-plasmaprofile) is exclusively in `plasma_profiles.py`
![UML of profiles](./uml_classes_PlasmaProfile.png){height="1000px"}
Figure 1: UML class breakdown of the plasma profiles
- ## Parabolic Profile | L-mode If `ipedestal == 0` then no pedestal is present and the function describing the profiles is given by: @@ -40,7 +39,6 @@ be described as 1/2-D. The relevant profile index variables are While PROCESS assumes a standard parabolic profile to be the shape of the current profile as per the 1989 ITER physics guidelines[^2] , it does not calculate its shape or values apart from the core value. The profile peaking factor `alphaj` is calculated in the plasma current calculation relating to `iprofile` found [here](../plasma_current.md).The on-axis current density is analytically calculated in [`calculate_profile_factors()`](#calculate_profile_factors) Only the temeprature and density profiles are calculated fully withing the `PlasmaProfiles` class. - The graph below is for a standard parabolic profile. You can vary the core value (`n0`) and the profile index (`alphan`) to see how the function behaves @@ -113,9 +111,8 @@ The graph below is for a standard parabolic profile. You can vary the core value If `ipedestal == 1` there is now a pedestal present in the profile and they follow the form shown below: - $$\begin{aligned} -\mbox{Density:} \qquad n(\rho) = \left\{ +\mbox{Density:} \qquad n(\rho) = \left\{ \begin{aligned} & \text{n}_{\text{ped}} + (n_0 - \text{n}_{\text{ped}}) \left( 1 - \frac{\rho^2}{\rho_{\text{ped},n}^2}\right)^{\alpha_n} @@ -126,9 +123,8 @@ $$\begin{aligned} \right. \end{aligned}$$ - $$\begin{aligned} -\mbox{Temperature:} \ T(\rho) = \left\{ +\mbox{Temperature:} \ T(\rho) = \left\{ \begin{aligned} & \text{T}_{\text{ped}} + (T_0 - \text{T}_{\text{ped}}) \left( 1 - \frac{\rho^{\beta_T}} {\rho_{\text{ped},T}^{\beta_T}}\right)^{\alpha_T} & 0 \leq \rho \leq \rho_{\text{ped},T} \\ @@ -136,7 +132,7 @@ $$\begin{aligned} & \rho_{\text{ped},T} < \rho \leq 1 \end{aligned} \right. -\end{aligned}$$ +\end{aligned}$$ Subscripts $0$, $\text{ped}$ and $\text{sep}$, denote values at the centre ($\rho = 0$), the pedestal ($\rho = \rho_{\text{ped}}$) and the separatrix ($\rho=1$), @@ -152,21 +148,19 @@ scaled from the electron values by the ratio of the volume-averaged values). $\beta_T$ can have an important impact not on just the plasma profile shape but also synchrotron radiation calculations. For more info on its effect, visit the radiation section [here](../plasma_radiation.md). - !!! warning " Pedestal setting" - If `ipedestal == 1` then the pedestal density `neped` is set as a fraction `fgwped` of the - Greenwald density (providing `fgwped` >= 0). The default value of `fgwped` is 0.8[^2]. + If `ipedestal == 1` then the pedestal density `neped` is set as a fraction `fgwped` of the + Greenwald density (providing `fgwped` >= 0). The default value of `fgwped` is 0.8[^2]. A table of the the associated variables can be seen below - -| Profile parameter | Density | Temperature | +| Profile parameter | Density | Temperature | |----------------------------------|-----------|-------------| | Pedestal radius (r/a) | `rhopedn`, $\rho_{\text{ped},n}$ | `rhopedt`, $\rho_{\text{ped},T}$ | -| Plasma centre value | `ne0`, $n_0$ | `te0`, $T_0$ | -| Pedestal value | `neped`, $n_{\text{ped}}$ | `teped`, $T_{\text{ped}}$ | -| Separatrix value | `nesep`, $n_{\text{sep}}$ | `tesep`, $T_{\text{sep}}$ | -| Profile index/ peaking parameter | `alphan`, $\alpha_n$ | `alphat`, $\alpha_T$ | +| Plasma centre value | `ne0`, $n_0$ | `te0`, $T_0$ | +| Pedestal value | `neped`, $n_{\text{ped}}$ | `teped`, $T_{\text{ped}}$ | +| Separatrix value | `nesep`, $n_{\text{sep}}$ | `tesep`, $T_{\text{sep}}$ | +| Profile index/ peaking parameter | `alphan`, $\alpha_n$ | `alphat`, $\alpha_T$ | | Profile index $\beta$ | | `tbeta`, $\beta_T$ | The graph below is for a standard pedestal profile. You can vary its attributes given in the table above to see how the function behaves. @@ -237,18 +231,18 @@ The graph below is for a standard pedestal profile. You can vary its attributes !!! warning " Un-realistic profiles" - If setting `ipedestal == 1` it is highly recommended to make sure constraint equation 81 (icc=81) is active. This enforces solutions in which $n_0$ has to be greater than $n_{\text{ped}}$. + If setting `ipedestal == 1` it is highly recommended to make sure constraint equation 81 (icc=81) is active. This enforces solutions in which $n_0$ has to be greater than $n_{\text{ped}}$. Negative $n_0$ values can also arise during iteration, so it is important to be weary on how low the lower bound for $n_{\text{e}} (\mathtt{dene})$ is set. More info can be found [here](plasma_profiles.md#pedestal-density-upper-limit) -------- -## Plasma Profile Class | `PlasmaProfile` +## Plasma Profile Class | `PlasmaProfile` ### Initialization | `__init__()` -The parent plasma profile class is `PlasmaProfile`. Initialization sets the profile class size and `neprofile` and `teprofile` to [`NProfile`](plasma_density_profile.md) & [`TProfile`](plasma_temperature_profile.md) objects from `profiles.py` +The parent plasma profile class is `PlasmaProfile`. Initialization sets the profile class size and `neprofile` and `teprofile` to [`NeProfile`](plasma_density_profile.md) & [`TeProfile`](plasma_temperature_profile.md) objects from `profiles.py` ???+ Note - Profile sizes are set to 501 point by default. this can be varied in the `__init__` of `PlasmaProfile`. Changing this will affect the values when doing Simpsons rule integration on the profiles. + Profile sizes are set to 501 point by default. this can be varied in the `__init__` of `PlasmaProfile`. Changing this will affect the values when doing Simpsons rule integration on the profiles. ---------- @@ -272,8 +266,7 @@ Depending on the value of `ipedestal` different functions will be ran, the diffe If pedestal profile values are set they are reset to have values that agree with the original form of the parabolic profiles. Such that $\rho_{\text{ped}} = 1$ and that pedestal and separatrix densities and temepratures are zero. This will then warn the user in the terminal. -The density and temperature profile runner function [`TProfile/NProfile.run()`](plasma_density_profile.md#runner-function-run) is then called to re-calculate the profile and core values. - +The density and temperature profile runner function [`TeProfile/NeProfile.run()`](plasma_density_profile.md#runner-function-run) is then called to re-calculate the profile and core values. Ratio of density-weighted to volume-averaged temperature factor is calculated: @@ -304,9 +297,8 @@ $$ $\Gamma$ is the [gamma function](https://en.wikipedia.org/wiki/Gamma_function) and is calculated in the code with [scipy.special.gamma()](https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.gamma.html) - $$ -\mathtt{dnla} / \bar{n_{\text{e}}} = \frac{n_0}{2}\frac{\Gamma(1/2)\Gamma(\alpha_n+1)}{\Gamma(\alpha_n+3/2)} +\mathtt{dnla} / \bar{n_{\text{e}}} = \frac{n_0}{2}\frac{\Gamma(1/2)\Gamma(\alpha_n+1)}{\Gamma(\alpha_n+3/2)} $$ This is in agreement with the derivation from the ITER Physics Design 1989 [^2] @@ -317,7 +309,6 @@ $\blacksquare$ The density weighted temperatures are set: - $$\begin{aligned} \mathtt{ten} = \mathtt{pcoef} \times T_\text{e} \\ \mathtt{tin} = \mathtt{pcoef}\times T_\text{i} @@ -381,7 +372,6 @@ $$ \langle n \rangle = \frac{2\pi \int^1_0 \rho \left(n_0(1-\rho^2)^{\alpha_n}\right) \ d\rho}{\pi\rho^2} $$ - $\blacksquare$ ------ @@ -400,7 +390,6 @@ $$\begin{aligned} ----- - ##### `calculate_profile_factors()` The central plasma pressure is calculated from the ideal gas law. @@ -425,7 +414,7 @@ $$ ???+ note "Pressure profile factor" - 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. + 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. @@ -502,7 +491,7 @@ This solution for the profile gradient only holds true if $\alpha_T \ge 1$. In the region $0 \le \alpha_T \le 1$ when substituting the roots of the second derivative into the first derivative the function diverges into the complex number solution space. To overcome this we can assume the second derivative root to be a value. In this case we assume a default value of $\mathtt{rho\_te\_max}$ = 0.9. -Then we just substitute this value into the first derivative to +Then we just substitute this value into the first derivative to $$ -2T_0\alpha_T(0.9)(1-(0.9)^2)^{-1+\alpha_T} @@ -575,9 +564,6 @@ You can use the slider in the graph below to experiment with the value of $\math - - - We can now set the normalised gradient length: $$ @@ -594,7 +580,7 @@ $$ ##### `pedestal_parameterisation()` -The density and temperature profile runner function [`TProfile/NProfile.run()`](plasma_density_profile.md#runner-function-run) is firstly called to re-calculate the profile and core values. +The density and temperature profile runner function [`TeProfile/NeProfile.run()`](plasma_density_profile.md#runner-function-run) is firstly called to re-calculate the profile and core values. Perform integrations to calculate ratio of density-weighted to volume-averaged temperature, etc. Density-weighted temperature = $\frac{\int{nT \ dV}}{\int{n \ dV}}$, which is approximately equal to the ratio $\frac{\int{\rho \ n(\rho) T(\rho) \ d\rho}}{\int{\rho \ n(\rho) \ d\rho}}$ @@ -637,9 +623,9 @@ The same function is run from the `ipedestal == 0 ` profile case, found [here](p -------- -### Setting pedestal values as fractions of the Greenwald limit +### Setting pedestal values as fractions of the Greenwald limit -By default, the values of $n_{\text{ped}}$ and $n_{\text{sep}}$ are set as fractions of the [Greenwald](https://wiki.fusion.ciemat.es/wiki/Greenwald_limit) limit such as: +By default, the values of $n_{\text{ped}}$ and $n_{\text{sep}}$ are set as fractions of the [Greenwald](https://wiki.fusion.ciemat.es/wiki/Greenwald_limit) limit such as: $$ n_{\text{ped}} = \mathtt{fgwped} \times \frac{I_p [\text{A}]}{\pi a^2 [\text{m}^2]} \times 10^{14} @@ -668,10 +654,10 @@ To prevent unrealistic profiles when iterating the values of `ne0, nesep, neped` This constraint can be activated by stating `icc = 76` in the input file -Eich et.al has given a usable formula for the critical separatrix density in terms of fundamental plasma parameters[^4]. A function for the ballooning parameter at the separatrix $(\alpha_{\text{sep}}^{\text{crit}})$ is seen to increase about -linearly with the separatrix density normalised to Greenwald density. This is seen in JET and ASDEX based on Thomson-scattering measurements, a clear correlation of the density limit of the tokamak H-mode high-confinement regime with -the approach to the ideal ballooning instability threshold at the periphery of the plasma. - +Eich et.al has given a usable formula for the critical separatrix density in terms of fundamental plasma parameters[^4]. A function for the ballooning parameter at the separatrix $(\alpha_{\text{sep}}^{\text{crit}})$ is seen to increase about +linearly with the separatrix density normalised to Greenwald density. This is seen in JET and ASDEX based on Thomson-scattering measurements, a clear correlation of the density limit of the tokamak H-mode high-confinement regime with +the approach to the ideal ballooning instability threshold at the periphery of the plasma. + $$ \alpha_{\text{sep}}^{\text{crit}} \approx \kappa^{1.2}(1+1.5\delta) \approx 2.5 $$ @@ -690,7 +676,6 @@ Where $A$ is the plasma aspect ratio, $P_{\text{sep}}$ is the power crossing the The value of `nesep` is check to make sure that it does not go higher than $n_{\text{sep}}^{\text{crit}}$. This can be scaled with `fnesep` - [^1]: M. Bernert et al. Plasma Phys. Control. Fus. **57** (2015) 014038 [^2]: N.A. Uckan and ITER Physics Group, 'ITER Physics Design Guidelines: 1989', [^3]: Johner Jean (2011) HELIOS: A Zero-Dimensional Tool for Next Step and Reactor Studies, Fusion Science and Technology, 59:2, 308-349, DOI: 10.13182/FST11-A11650 diff --git a/documentation/proc-pages/physics-models/profiles/plasma_temperature_profile.md b/documentation/proc-pages/physics-models/profiles/plasma_temperature_profile.md index 642ba8924e..1920fb4e50 100644 --- a/documentation/proc-pages/physics-models/profiles/plasma_temperature_profile.md +++ b/documentation/proc-pages/physics-models/profiles/plasma_temperature_profile.md @@ -1,4 +1,4 @@ -# Temperature Profile | `TProfile(Profile)` +# Temperature Profile | `TeProfile(Profile)` The temperature profile class is organised around a central runner function that is called each time the plasma is parameterised by the parent [`PlasmaProfile()`](plasma_profiles.md#plasma-profile-class-plasmaprofile) class. It is called by [`pedestal_parameterisation()`](plasma_profiles.md#pedestal_parameterisation) and [`parabolic parameterisation()`](plasma_profiles.md#parabolic_paramterisation). The sequence of the runner function can be seen below along with explanation of the following calculations. @@ -8,14 +8,13 @@ The temperature profile class is organised around a central runner function that 2. The steps between the normalized points is then done by [`calculate_profile_dx()`](./plasma_profiles_abstract_class.md#calculate-the-profile-steps-in-x--calculate_profile_dx) which divides the max x-dimension by the number of points. -3. The core electron and ion temperatures are claculated via [`set_physics_variables()`]() +3. The core electron and ion temperatures are claculated via [`set_physics_variables()`]() - - ### Calculate core values | `set_physics_variables()` + ### Calculate core values | `set_physics_variables()` The core electron temperature is calculated using the [`tcore`](plasma_temperature_profile.md#electron-core-density-of-a-pedestalised-profile--tcore) method. - - #### Electron core density of a pedestalised profile | `tcore()` + + #### Electron core density of a pedestalised profile | `tcore()` This function calculates the core electron density for a pedestalsied profile in $\text{keV}$. The inclusion of a new $\beta_T$ exponent term allows a more accurate description of temperature profiles with a triangular shape or a strong gradient near the pedestal (characteristic of regimes with an [internal transport barrier](https://wiki.fusion.ciemat.es/wiki/Internal_Transport_Barrier)). @@ -30,7 +29,6 @@ The temperature profile class is organised around a central runner function that | Profile index/ peaking parameter | `alphan`, $\alpha_T$ | | Profile index/ peaking parameter | `tbeta`, $\beta_T$ | - $$ T_0 = T_{\text{ped}}+\frac{\beta_T(3\langle T_{\text{e}} \rangle +T_{\text{sep}}(-2+\rho_{\text{ped}} +\rho_{\text{ped}}^2)-T_{\text{ped}}(1+\rho_{\text{ped}}+\rho_{\text{ped}}^2))}{6\rho_{\text{ped}}^2\text{B}\left[1+\alpha_T,\frac{2}{\beta_T}\right]} $$ @@ -39,7 +37,7 @@ The temperature profile class is organised around a central runner function that --------- - ##### Derivation + ##### Derivation We calculate the volume integrated profile and then divide by the volume of integration to get the volume average density $\langle T_{\text{e}} \rangle$. If we assume the plasma to be a torus of circular cross-section then we can use spherical coordinates. We can simplify the problem by representing the torus as a cylinder of height equal to the circumference of the torus equal to $2\pi R$ where $R$ is the major radius of the torus, and $a$ is the plasma minor radius in the poloidal plane. @@ -72,23 +70,23 @@ The temperature profile class is organised around a central runner function that In the form of volume average temperature where the volume integrated temperature function has to be divided by the volume of the cylinder / torus, within the volume bounded by that pedestal position we get: $$ - \langle T_{\text{e}} \rangle = 4\pi^2R\left[ \frac{\frac{\left(T_{\text{ped}}\beta_T+(2T_0-2T_{\text{ped}})B\left(\alpha_T+1,\frac{2}{\beta_T}\right)\right)\rho_{\text{ped},T}^2}{2\beta_T}+\frac{(1-\rho_{\text{ped},T})\left((T_{\text{sep}}+2T_{\text{ped}})\rho_\text{ped}+2T_{\text{sep}}+T_{\text{ped}}\right)}{6}}{2\pi^2 R \rho^2}\right] + \langle T_{\text{e}} \rangle = 4\pi^2R\left[ \frac{\frac{\left(T_{\text{ped}}\beta_T+(2T_0-2T_{\text{ped}})B\left(\alpha_T+1,\frac{2}{\beta_T}\right)\right)\rho_{\text{ped},T}^2}{2\beta_T}+\frac{(1-\rho_{\text{ped},T})\left((T_{\text{sep}}+2T_{\text{ped}})\rho_\text{ped}+2T_{\text{sep}}+T_{\text{ped}}\right)}{6}}{2\pi^2 R \rho^2}\right] $$ In this case, the value of $\rho$ is equal to 1 as we integrated over the full profile. $$ \langle T_{\text{e}} \rangle = 2\left[ \frac{\left(T_{\text{ped}}\beta_T+(2T_0-2T_{\text{ped}})\text{B}\left(\alpha_T+1,\frac{2}{\beta_T}\right)\right)\rho_{\text{ped},T}^2}{2\beta_T} \\ - +\frac{(1-\rho_{\text{ped},T})\left((T_{\text{sep}}+2T_{\text{ped}})\rho_\text{ped}+2T_{\text{sep}}+T_{\text{ped}}\right)}{6}\right] + +\frac{(1-\rho_{\text{ped},T})\left((T_{\text{sep}}+2T_{\text{ped}})\rho_\text{ped}+2T_{\text{sep}}+T_{\text{ped}}\right)}{6}\right] $$ $$ \langle T_{\text{e}} \rangle = \frac{\left(T_{\text{ped}}\beta_T+(2T_0-2T_{\text{ped}})\text{B}\left[\alpha_T+1,\frac{2}{\beta_T}\right]\right)\rho_{\text{ped},T}^2}{\beta_T} \\ - +\frac{(1-\rho_{\text{ped},T})\left((T_{\text{sep}}+2T_{\text{ped}})\rho_\text{ped}+2T_{\text{sep}}+T_{\text{ped}}\right)}{3} + +\frac{(1-\rho_{\text{ped},T})\left((T_{\text{sep}}+2T_{\text{ped}})\rho_\text{ped}+2T_{\text{sep}}+T_{\text{ped}}\right)}{3} $$ - + Where $\text{B}$ is the [Beta function](https://en.wikipedia.org/wiki/Beta_function) - + Re-arranging to get $T_0$ we get: $$ @@ -105,12 +103,11 @@ The temperature profile class is organised around a central runner function that T_{\text{i0}} = \left(\frac{T_{\text{i}}}{T_{\text{e}}}\right)T_{\text{e0}} $$ - 4. The y profile is then calculated using [`calculate_profile_y()`](plasma_temperature_profile.md#calculate-temperature-at-each-radius-position-calculate_profile_y). This routine calculates the temperature at each normalised minor radius position $\rho$ for a HELIOS-type temperature pedestal profile[^1] ------- - ### Calculate temperature at each radius position | `calculate_profile_y()` + ### Calculate temperature at each radius position | `calculate_profile_y()` A table of the input variables can be found below @@ -127,7 +124,7 @@ The temperature profile class is organised around a central runner function that If `ipedestal == 0` then the original parabolic profile form is used $$ - T(\rho) = T_0(1 - \rho^2)^{\alpha_T} + T(\rho) = T_0(1 - \rho^2)^{\alpha_T} $$ The central temperature ($T_0$) is then checked to make sure it is not less than the pedestal temperature, $T_{\text{ped}}$. @@ -135,10 +132,8 @@ The temperature profile class is organised around a central runner function that Values of the profile temperature are then assigned based on the density function below across bounds from 0 to `rhopedn` and `rhopedn` to 1. - - $$\begin{aligned} - \mbox{Temperature:} \ \ T(\rho) = \left\{ + \mbox{Temperature:} \ \ T(\rho) = \left\{ \begin{aligned} & \text{T}_{\text{ped}} + (T_0 - \text{T}_{\text{ped}}) \left( 1 - \frac{\rho^{\beta_T}} {\rho_{\text{ped},T}^{\beta_T}}\right)^{\alpha_T} \ 0 \leq \rho \leq \rho_{\text{ped},T} \\ @@ -146,9 +141,8 @@ The temperature profile class is organised around a central runner function that \ \rho_{\text{ped},T} < \rho \leq 1 \end{aligned} \right. - \end{aligned}$$ - + \end{aligned}$$ 5. Profile is then integrated with `integrate_profile_y()` using Simpsons integration from the profile abstract base class -[^1]: Jean, J. (2011). *HELIOS: A Zero-Dimensional Tool for Next Step and Reactor Studies*. Fusion Science and Technology, 59(2), 308–349. https://doi.org/10.13182/FST11-A11650 \ No newline at end of file +[^1]: Jean, J. (2011). *HELIOS: A Zero-Dimensional Tool for Next Step and Reactor Studies*. Fusion Science and Technology, 59(2), 308–349. From d8c2ec84f49252829d856c2cbb69a020772c3821 Mon Sep 17 00:00:00 2001 From: ym1906 Date: Mon, 17 Feb 2025 10:12:04 +0000 Subject: [PATCH 3/4] update profile names --- process/io/plot_proc.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/process/io/plot_proc.py b/process/io/plot_proc.py index f9f1cb6b3f..f3600db105 100644 --- a/process/io/plot_proc.py +++ b/process/io/plot_proc.py @@ -747,7 +747,7 @@ def arc_fill(axis, r1, r2, color="pink"): axis.add_patch(patch) -def plot_neprofile(prof, demo_ranges): +def plot_n_profiles(prof, demo_ranges): """Function to plot density profile Arguments: prof --> axis object to add plot to @@ -882,7 +882,7 @@ def plot_jprofile(prof): ) -def plot_teprofile(prof, demo_ranges): +def plot_t_profiles(prof, demo_ranges): """Function to plot temperature profile Arguments: prof --> axis object to add plot to @@ -3384,12 +3384,12 @@ def main_plot( # Plot density profiles plot_4 = fig2.add_subplot(231) # , aspect= 0.05) plot_4.set_position([0.075, 0.55, 0.25, 0.4]) - plot_neprofile(plot_4, demo_ranges) + plot_n_profiles(plot_4, demo_ranges) # Plot temperature profiles plot_5 = fig2.add_subplot(232) plot_5.set_position([0.375, 0.55, 0.25, 0.4]) - plot_teprofile(plot_5, demo_ranges) + plot_t_profiles(plot_5, demo_ranges) # Plot impurity profiles plot_8 = fig2.add_subplot(233) From 23c71a44f36caddc0ee6a6c79ae58066cf104b62 Mon Sep 17 00:00:00 2001 From: ym1906 Date: Mon, 17 Feb 2025 10:13:07 +0000 Subject: [PATCH 4/4] added explicity reference to electron --- process/impurity_radiation.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/process/impurity_radiation.py b/process/impurity_radiation.py index beb962f55e..28b56b1547 100644 --- a/process/impurity_radiation.py +++ b/process/impurity_radiation.py @@ -395,9 +395,9 @@ def pimpden(imp_element_index, neprofile, teprofile): :param imp_element_index: Impurity element index :type imp_element_index: int - :param neprofile: density profile + :param neprofile: electron density profile :type neprofile: numpy.array - :param teprofile: temperature profile + :param teprofile: electron temperature profile :type teprofile: numpy.array :return: pimpden - total impurity radiation density (W/m3) :rtype: numpy.array