Skip to content

profiles

Module for plasma profile definitions and utilities.

This module provides abstract base classes and implementations for creating and managing plasma profiles such as temperature and density distributions.

PlasmaProfileShapeType

Bases: IntEnum

Enum for i_plasma_pedestal method types

Source code in process/models/physics/profiles.py
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
@unique
class PlasmaProfileShapeType(IntEnum):
    """Enum for i_plasma_pedestal method types"""

    PARABOLIC_PROFILE = (0, "Parabolic Profile (L-mode)")
    PEDESTAL_PROFILE = (1, "Pedestal Profile (H-mode)")

    def __new__(cls, value: int, description: str):
        """Create a new PlasmaProfileShapeType instance."""
        obj = int.__new__(cls, value)
        obj._value_ = value
        obj._description_ = description
        return obj

    @DynamicClassAttribute
    def description(self):
        """The description of the plasma profile shape."""
        return self._description_

PARABOLIC_PROFILE = (0, 'Parabolic Profile (L-mode)') class-attribute instance-attribute

PEDESTAL_PROFILE = (1, 'Pedestal Profile (H-mode)') class-attribute instance-attribute

description()

The description of the plasma profile shape.

Source code in process/models/physics/profiles.py
36
37
38
39
@DynamicClassAttribute
def description(self):
    """The description of the plasma profile shape."""
    return self._description_

Profile

Bases: Model, ABC

Abstract base class used to create and hold profiles (temperature, density)

Source code in process/models/physics/profiles.py
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
class Profile(Model, ABC):
    """Abstract base class used to create and hold profiles (temperature, density)"""

    def __init__(self):
        """
        Initialize a Profiles object.

        Attributes
        ----------
        - profile_size (int): The size of the profile.
        - profile_x (ndarray): An array of values ranging from 0 to profile_size-1.
        - profile_y (ndarray): An array of zeros with length profile_size.
        - profile_integ (int): The integral of the profile_y array.
        - profile_dx (int): The step size between consecutive values in profile_x.
        """
        self.profile_integ = 0
        self.profile_dx = 0

    def run(self):
        """Initialise profile_x and profile_y"""
        self.profile_x = np.arange(
            self.data.physics.n_plasma_profile_elements, dtype=float
        )
        self.profile_y = np.zeros(self.data.physics.n_plasma_profile_elements)

    def output(self):
        """Profile model doesn't have any output"""

    def normalise_profile_x(self) -> None:
        """Normalizes the x-dimension of the profile.

        This method divides the values in the `profile_x` attribute by the maximum value
        in the `profile_x` array, resulting in a normalized version of the x-dimension.

        Example:
            If `profile_x` is [1, 2, 3, 4, 5], after normalization it will become
            [0.2, 0.4, 0.6, 0.8, 1.0].

        Note:
            This method modifies the `profile_x` attribute in-place.


        """
        self.profile_x /= max(self.profile_x)

    def calculate_profile_dx(self) -> None:
        """Calculates the differential between points in the profile.

        This method calculates the differential between points in the profile by
        dividing the maximum x value in the profile by the difference in size between
        the points. The result is stored in the `profile_dx` attribute.
        """
        self.profile_dx = max(self.profile_x) / (
            self.data.physics.n_plasma_profile_elements - 1
        )

    @abstractmethod
    def calculate_profile_y(self) -> None:
        """Use a profile function to act on self.profile_x to calculate and set the
        values of self.profile_y.
        """

    def integrate_profile_y(self) -> None:
        """Integrate profile_y values using scipy.integrate.simpson() function.

        This method calculates the integral of the profile_y values using the Simpson's
        rule provided by the scipy.integrate.simpson() function. The integral is stored
        in the `profile_integ` attribute.
        """
        self.profile_integ = sp.integrate.simpson(
            self.profile_y, x=self.profile_x, dx=self.profile_dx
        )

profile_integ = 0 instance-attribute

profile_dx = 0 instance-attribute

run()

Initialise profile_x and profile_y

Source code in process/models/physics/profiles.py
60
61
62
63
64
65
def run(self):
    """Initialise profile_x and profile_y"""
    self.profile_x = np.arange(
        self.data.physics.n_plasma_profile_elements, dtype=float
    )
    self.profile_y = np.zeros(self.data.physics.n_plasma_profile_elements)

output()

Profile model doesn't have any output

Source code in process/models/physics/profiles.py
67
68
def output(self):
    """Profile model doesn't have any output"""

normalise_profile_x()

Normalizes the x-dimension of the profile.

This method divides the values in the profile_x attribute by the maximum value in the profile_x array, resulting in a normalized version of the x-dimension.

Example: If profile_x is [1, 2, 3, 4, 5], after normalization it will become [0.2, 0.4, 0.6, 0.8, 1.0].

Note: This method modifies the profile_x attribute in-place.

Source code in process/models/physics/profiles.py
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
def normalise_profile_x(self) -> None:
    """Normalizes the x-dimension of the profile.

    This method divides the values in the `profile_x` attribute by the maximum value
    in the `profile_x` array, resulting in a normalized version of the x-dimension.

    Example:
        If `profile_x` is [1, 2, 3, 4, 5], after normalization it will become
        [0.2, 0.4, 0.6, 0.8, 1.0].

    Note:
        This method modifies the `profile_x` attribute in-place.


    """
    self.profile_x /= max(self.profile_x)

calculate_profile_dx()

Calculates the differential between points in the profile.

This method calculates the differential between points in the profile by dividing the maximum x value in the profile by the difference in size between the points. The result is stored in the profile_dx attribute.

Source code in process/models/physics/profiles.py
87
88
89
90
91
92
93
94
95
96
def calculate_profile_dx(self) -> None:
    """Calculates the differential between points in the profile.

    This method calculates the differential between points in the profile by
    dividing the maximum x value in the profile by the difference in size between
    the points. The result is stored in the `profile_dx` attribute.
    """
    self.profile_dx = max(self.profile_x) / (
        self.data.physics.n_plasma_profile_elements - 1
    )

calculate_profile_y() abstractmethod

Use a profile function to act on self.profile_x to calculate and set the values of self.profile_y.

Source code in process/models/physics/profiles.py
 98
 99
100
101
102
@abstractmethod
def calculate_profile_y(self) -> None:
    """Use a profile function to act on self.profile_x to calculate and set the
    values of self.profile_y.
    """

integrate_profile_y()

Integrate profile_y values using scipy.integrate.simpson() function.

This method calculates the integral of the profile_y values using the Simpson's rule provided by the scipy.integrate.simpson() function. The integral is stored in the profile_integ attribute.

Source code in process/models/physics/profiles.py
104
105
106
107
108
109
110
111
112
113
def integrate_profile_y(self) -> None:
    """Integrate profile_y values using scipy.integrate.simpson() function.

    This method calculates the integral of the profile_y values using the Simpson's
    rule provided by the scipy.integrate.simpson() function. The integral is stored
    in the `profile_integ` attribute.
    """
    self.profile_integ = sp.integrate.simpson(
        self.profile_y, x=self.profile_x, dx=self.profile_dx
    )

DensityProfilePedestalType

Bases: IntEnum

Enum for i_nd_plasma_pedestal_separatrix types

Source code in process/models/physics/profiles.py
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
@unique
class DensityProfilePedestalType(IntEnum):
    """Enum for i_nd_plasma_pedestal_separatrix types"""

    USER_INPUT = (0, "User input direct values")
    GREENWALD_FRACTION = (1, "Fractions of the Greenwald limit")

    def __new__(cls, value: int, description: str):
        """Create a new DensityProfilePedestalType instance.

        Parameters
        ----------
            value: Integer value for the enum member.
            description: Human-readable description for the enum member.
        """
        obj = int.__new__(cls, value)
        obj._value_ = value
        obj._description_ = description
        return obj

    @DynamicClassAttribute
    def description(self):
        """The description of the plasma profile shape."""
        return self._description_

USER_INPUT = (0, 'User input direct values') class-attribute instance-attribute

GREENWALD_FRACTION = (1, 'Fractions of the Greenwald limit') class-attribute instance-attribute

description()

The description of the plasma profile shape.

Source code in process/models/physics/profiles.py
136
137
138
139
@DynamicClassAttribute
def description(self):
    """The description of the plasma profile shape."""
    return self._description_

ElectronDensityProfile

Bases: Profile

Electron density (nₑ) profile class. Contains a function to calculate the electron density profile and store the data.

Source code in process/models/physics/profiles.py
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
class ElectronDensityProfile(Profile):
    """Electron density (nₑ) profile class. Contains a function to calculate the electron
    density profile and store the data.
    """

    def run(self):
        """Subroutine which calls profile functions and stores neprofile data."""
        super().run()
        self.normalise_profile_x()
        self.calculate_profile_dx()
        self.set_physics_variables()
        self.calculate_profile_y(
            rho=self.profile_x,
            radius_plasma_pedestal_density_norm=self.data.physics.radius_plasma_pedestal_density_norm,
            nd_on_axis=self.data.physics.nd_plasma_electron_on_axis,
            nd_pedestal=self.data.physics.nd_plasma_pedestal_electron,
            nd_separatrix=self.data.physics.nd_plasma_separatrix_electron,
            alphan=self.data.physics.alphan,
        )
        self.integrate_profile_y()

    def calculate_profile_y(
        self,
        rho: np.array,
        radius_plasma_pedestal_density_norm: float,
        nd_on_axis: float,
        nd_pedestal: float,
        nd_separatrix: float,
        alphan: float,
    ) -> None:
        """Calculates the number density at each normalised minor radius (ρ) position.

        Parameters
        ----------
        rho :
            Normalised minor radius (ρ) vector.
        radius_plasma_pedestal_density_norm :
            Normalised minor radius pedestal position (ρₙ,pedestal).
        nd_on_axis :
            Central number density (n₀) [m⁻³].
        nd_pedestal :
            Pedestal density (n_pedestal) [m⁻³].
        nd_separatrix :
            Separatrix density (n_sep) [m⁻³].
        alphan :
            Density peaking parameter (αₙ).
        """  # noqa: RUF002
        if (
            PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
            == PlasmaProfileShapeType.PARABOLIC_PROFILE
        ):
            self.profile_y = nd_on_axis * (1 - rho**2) ** alphan

        # Input checks

        if nd_on_axis < nd_pedestal:
            logger.info(
                "NPROFILE: Pedestal density is higher than core density. %s, %s",
                nd_pedestal,
                nd_on_axis,
            )
        rho_index = rho <= radius_plasma_pedestal_density_norm
        self.profile_y[rho_index] = (
            nd_pedestal
            + (nd_on_axis - nd_pedestal)
            * (1 - (rho[rho_index] / radius_plasma_pedestal_density_norm) ** 2) ** alphan
        )
        # Invert the rho_index
        self.profile_y[~rho_index] = nd_separatrix + (nd_pedestal - nd_separatrix) * (
            1 - rho[~rho_index]
        ) / (1 - radius_plasma_pedestal_density_norm)

    @staticmethod
    def calculate_pedestal_profile_on_axis_density(
        radius_plasma_pedestal_density_norm: float,
        nd_pedestal: float,
        nd_separatrix: float,
        nd_vol_average: float,
        alphan: float,
    ) -> float:
        """Calculates the core density (n₀) of a pedestalised profile.

        Parameters
        ----------
        radius_plasma_pedestal_density_norm :
            Normalised minor radius pedestal position (ρₙ,pedestal).
        nd_pedestal: float,
            The pedestal density (n_pedestal) [m⁻³].
        nd_separatrix: float,
            The separatrix density (n_sep) [m⁻³].
        nd_vol_average: float,
            The volume averaged density (⟨n⟩) [m⁻³].
        alphan: float,
            The density peaking parameter (αₙ).

        Returns
        -------
        :
            The core on-axis density (n₀) [m⁻³].
        """
        nd_on_axis = (
            1
            / (3 * radius_plasma_pedestal_density_norm**2)
            * (
                3 * nd_vol_average * (1 + alphan)
                + nd_separatrix
                * (1 + alphan)
                * (
                    -2
                    + radius_plasma_pedestal_density_norm
                    + radius_plasma_pedestal_density_norm**2
                )
                - nd_pedestal
                * (
                    (1 + alphan) * (1 + radius_plasma_pedestal_density_norm)
                    + (alphan - 2) * radius_plasma_pedestal_density_norm**2
                )
            )
        )

        if nd_on_axis < 0.0:
            # Allows solver to continue and
            # warns the user to raise the lower bound on nd_plasma_electrons_vol_avg
            # if the run did not converge
            logger.error(
                "nd_on_axis is going negative when solving. Please raise the value of "
                "nd_plasma_electrons_vol_avg (⟨nₑ⟩) and or its lower limit."
            )
            nd_on_axis = 1.0e-6
        return nd_on_axis

    @staticmethod
    def calculate_parabolic_profile_on_axis_density(
        nd_vol_average: float,
        alphan: float,
    ) -> float:
        """Calculates the core density (n₀) of a parabolic profile.

        Parameters
        ----------
        nd_vol_average: float,
            The volume averaged density (⟨n⟩) [m⁻³].
        alphan: float,
            The density peaking parameter (αₙ).

        Returns
        -------
        :
            The core on-axis density (n₀) [m⁻³].
        """
        nd_on_axis = nd_vol_average * (1.0 + alphan)

        if nd_on_axis < 0.0:
            # Allows solver to continue and
            # warns the user to raise the lower bound on nd_plasma_electrons_vol_avg
            # if the run did not converge
            logger.error(
                "nd_on_axis is going negative when solving. Please raise the value of "
                "nd_plasma_electrons_vol_avg (⟨nₑ⟩) and or its lower limit."
            )
            nd_on_axis = 1.0e-6
        return nd_on_axis

    def set_pedestal_and_separatrix_values(self):
        """Sets the pedestal and separatrix density values based on the user input
        or greenwald fraction method.
        """
        i_nd_plasma_pedestal_separatrix = DensityProfilePedestalType(
            self.data.physics.i_nd_plasma_pedestal_separatrix
        )

        if i_nd_plasma_pedestal_separatrix == DensityProfilePedestalType.USER_INPUT:
            self.data.physics.f_nd_plasma_pedestal_greenwald = (
                self.data.physics.nd_plasma_pedestal_electron
                / (
                    PlasmaDensityLimit.calculate_greenwald_density_limit(
                        c_plasma=self.data.physics.plasma_current,
                        rminor=self.data.physics.rminor,
                    )
                )
            )

            self.data.physics.f_nd_plasma_separatrix_greenwald = (
                self.data.physics.nd_plasma_separatrix_electron
                / (
                    PlasmaDensityLimit.calculate_greenwald_density_limit(
                        c_plasma=self.data.physics.plasma_current,
                        rminor=self.data.physics.rminor,
                    )
                )
            )
        elif (
            i_nd_plasma_pedestal_separatrix
            == DensityProfilePedestalType.GREENWALD_FRACTION
        ):
            self.data.physics.nd_plasma_pedestal_electron = (
                self.data.physics.f_nd_plasma_pedestal_greenwald
                * PlasmaDensityLimit.calculate_greenwald_density_limit(
                    c_plasma=self.data.physics.plasma_current,
                    rminor=self.data.physics.rminor,
                )
            )
            self.data.physics.nd_plasma_separatrix_electron = (
                self.data.physics.f_nd_plasma_separatrix_greenwald
                * PlasmaDensityLimit.calculate_greenwald_density_limit(
                    c_plasma=self.data.physics.plasma_current,
                    rminor=self.data.physics.rminor,
                )
            )

    def set_physics_variables(self):
        """Calculates and sets physics variables required for the profile."""
        if (
            PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
            == PlasmaProfileShapeType.PARABOLIC_PROFILE
        ):
            self.data.physics.nd_plasma_electron_on_axis = (
                self.calculate_parabolic_profile_on_axis_density(
                    nd_vol_average=self.data.physics.nd_plasma_electrons_vol_avg,
                    alphan=self.data.physics.alphan,
                )
            )
        elif (
            PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
            == PlasmaProfileShapeType.PEDESTAL_PROFILE
        ):
            self.data.physics.nd_plasma_electron_on_axis = self.calculate_pedestal_profile_on_axis_density(  # noqa: E501
                radius_plasma_pedestal_density_norm=self.data.physics.radius_plasma_pedestal_density_norm,
                nd_pedestal=self.data.physics.nd_plasma_pedestal_electron,
                nd_separatrix=self.data.physics.nd_plasma_separatrix_electron,
                nd_vol_average=self.data.physics.nd_plasma_electrons_vol_avg,
                alphan=self.data.physics.alphan,
            )
        self.data.physics.nd_plasma_ions_on_axis = (
            self.data.physics.nd_plasma_ions_total_vol_avg
            / self.data.physics.nd_plasma_electrons_vol_avg
            * self.data.physics.nd_plasma_electron_on_axis
        )

run()

Subroutine which calls profile functions and stores neprofile data.

Source code in process/models/physics/profiles.py
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
def run(self):
    """Subroutine which calls profile functions and stores neprofile data."""
    super().run()
    self.normalise_profile_x()
    self.calculate_profile_dx()
    self.set_physics_variables()
    self.calculate_profile_y(
        rho=self.profile_x,
        radius_plasma_pedestal_density_norm=self.data.physics.radius_plasma_pedestal_density_norm,
        nd_on_axis=self.data.physics.nd_plasma_electron_on_axis,
        nd_pedestal=self.data.physics.nd_plasma_pedestal_electron,
        nd_separatrix=self.data.physics.nd_plasma_separatrix_electron,
        alphan=self.data.physics.alphan,
    )
    self.integrate_profile_y()

calculate_profile_y(rho, radius_plasma_pedestal_density_norm, nd_on_axis, nd_pedestal, nd_separatrix, alphan)

Calculates the number density at each normalised minor radius (ρ) position.

Parameters:

Name Type Description Default
rho array

Normalised minor radius (ρ) vector.

required
radius_plasma_pedestal_density_norm float

Normalised minor radius pedestal position (ρₙ,pedestal).

required
nd_on_axis float

Central number density (n₀) [m⁻³].

required
nd_pedestal float

Pedestal density (n_pedestal) [m⁻³].

required
nd_separatrix float

Separatrix density (n_sep) [m⁻³].

required
alphan float

Density peaking parameter (αₙ).

required
Source code in process/models/physics/profiles.py
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
def calculate_profile_y(
    self,
    rho: np.array,
    radius_plasma_pedestal_density_norm: float,
    nd_on_axis: float,
    nd_pedestal: float,
    nd_separatrix: float,
    alphan: float,
) -> None:
    """Calculates the number density at each normalised minor radius (ρ) position.

    Parameters
    ----------
    rho :
        Normalised minor radius (ρ) vector.
    radius_plasma_pedestal_density_norm :
        Normalised minor radius pedestal position (ρₙ,pedestal).
    nd_on_axis :
        Central number density (n₀) [m⁻³].
    nd_pedestal :
        Pedestal density (n_pedestal) [m⁻³].
    nd_separatrix :
        Separatrix density (n_sep) [m⁻³].
    alphan :
        Density peaking parameter (αₙ).
    """  # noqa: RUF002
    if (
        PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
        == PlasmaProfileShapeType.PARABOLIC_PROFILE
    ):
        self.profile_y = nd_on_axis * (1 - rho**2) ** alphan

    # Input checks

    if nd_on_axis < nd_pedestal:
        logger.info(
            "NPROFILE: Pedestal density is higher than core density. %s, %s",
            nd_pedestal,
            nd_on_axis,
        )
    rho_index = rho <= radius_plasma_pedestal_density_norm
    self.profile_y[rho_index] = (
        nd_pedestal
        + (nd_on_axis - nd_pedestal)
        * (1 - (rho[rho_index] / radius_plasma_pedestal_density_norm) ** 2) ** alphan
    )
    # Invert the rho_index
    self.profile_y[~rho_index] = nd_separatrix + (nd_pedestal - nd_separatrix) * (
        1 - rho[~rho_index]
    ) / (1 - radius_plasma_pedestal_density_norm)

calculate_pedestal_profile_on_axis_density(radius_plasma_pedestal_density_norm, nd_pedestal, nd_separatrix, nd_vol_average, alphan) staticmethod

Calculates the core density (n₀) of a pedestalised profile.

Parameters:

Name Type Description Default
radius_plasma_pedestal_density_norm float

Normalised minor radius pedestal position (ρₙ,pedestal).

required
nd_pedestal float

The pedestal density (n_pedestal) [m⁻³].

required
nd_separatrix float

The separatrix density (n_sep) [m⁻³].

required
nd_vol_average float

The volume averaged density (⟨n⟩) [m⁻³].

required
alphan float

The density peaking parameter (αₙ).

required

Returns:

Type Description
float

The core on-axis density (n₀) [m⁻³].

Source code in process/models/physics/profiles.py
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
@staticmethod
def calculate_pedestal_profile_on_axis_density(
    radius_plasma_pedestal_density_norm: float,
    nd_pedestal: float,
    nd_separatrix: float,
    nd_vol_average: float,
    alphan: float,
) -> float:
    """Calculates the core density (n₀) of a pedestalised profile.

    Parameters
    ----------
    radius_plasma_pedestal_density_norm :
        Normalised minor radius pedestal position (ρₙ,pedestal).
    nd_pedestal: float,
        The pedestal density (n_pedestal) [m⁻³].
    nd_separatrix: float,
        The separatrix density (n_sep) [m⁻³].
    nd_vol_average: float,
        The volume averaged density (⟨n⟩) [m⁻³].
    alphan: float,
        The density peaking parameter (αₙ).

    Returns
    -------
    :
        The core on-axis density (n₀) [m⁻³].
    """
    nd_on_axis = (
        1
        / (3 * radius_plasma_pedestal_density_norm**2)
        * (
            3 * nd_vol_average * (1 + alphan)
            + nd_separatrix
            * (1 + alphan)
            * (
                -2
                + radius_plasma_pedestal_density_norm
                + radius_plasma_pedestal_density_norm**2
            )
            - nd_pedestal
            * (
                (1 + alphan) * (1 + radius_plasma_pedestal_density_norm)
                + (alphan - 2) * radius_plasma_pedestal_density_norm**2
            )
        )
    )

    if nd_on_axis < 0.0:
        # Allows solver to continue and
        # warns the user to raise the lower bound on nd_plasma_electrons_vol_avg
        # if the run did not converge
        logger.error(
            "nd_on_axis is going negative when solving. Please raise the value of "
            "nd_plasma_electrons_vol_avg (⟨nₑ⟩) and or its lower limit."
        )
        nd_on_axis = 1.0e-6
    return nd_on_axis

calculate_parabolic_profile_on_axis_density(nd_vol_average, alphan) staticmethod

Calculates the core density (n₀) of a parabolic profile.

Parameters:

Name Type Description Default
nd_vol_average float

The volume averaged density (⟨n⟩) [m⁻³].

required
alphan float

The density peaking parameter (αₙ).

required

Returns:

Type Description
float

The core on-axis density (n₀) [m⁻³].

Source code in process/models/physics/profiles.py
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
@staticmethod
def calculate_parabolic_profile_on_axis_density(
    nd_vol_average: float,
    alphan: float,
) -> float:
    """Calculates the core density (n₀) of a parabolic profile.

    Parameters
    ----------
    nd_vol_average: float,
        The volume averaged density (⟨n⟩) [m⁻³].
    alphan: float,
        The density peaking parameter (αₙ).

    Returns
    -------
    :
        The core on-axis density (n₀) [m⁻³].
    """
    nd_on_axis = nd_vol_average * (1.0 + alphan)

    if nd_on_axis < 0.0:
        # Allows solver to continue and
        # warns the user to raise the lower bound on nd_plasma_electrons_vol_avg
        # if the run did not converge
        logger.error(
            "nd_on_axis is going negative when solving. Please raise the value of "
            "nd_plasma_electrons_vol_avg (⟨nₑ⟩) and or its lower limit."
        )
        nd_on_axis = 1.0e-6
    return nd_on_axis

set_pedestal_and_separatrix_values()

Sets the pedestal and separatrix density values based on the user input or greenwald fraction method.

Source code in process/models/physics/profiles.py
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
def set_pedestal_and_separatrix_values(self):
    """Sets the pedestal and separatrix density values based on the user input
    or greenwald fraction method.
    """
    i_nd_plasma_pedestal_separatrix = DensityProfilePedestalType(
        self.data.physics.i_nd_plasma_pedestal_separatrix
    )

    if i_nd_plasma_pedestal_separatrix == DensityProfilePedestalType.USER_INPUT:
        self.data.physics.f_nd_plasma_pedestal_greenwald = (
            self.data.physics.nd_plasma_pedestal_electron
            / (
                PlasmaDensityLimit.calculate_greenwald_density_limit(
                    c_plasma=self.data.physics.plasma_current,
                    rminor=self.data.physics.rminor,
                )
            )
        )

        self.data.physics.f_nd_plasma_separatrix_greenwald = (
            self.data.physics.nd_plasma_separatrix_electron
            / (
                PlasmaDensityLimit.calculate_greenwald_density_limit(
                    c_plasma=self.data.physics.plasma_current,
                    rminor=self.data.physics.rminor,
                )
            )
        )
    elif (
        i_nd_plasma_pedestal_separatrix
        == DensityProfilePedestalType.GREENWALD_FRACTION
    ):
        self.data.physics.nd_plasma_pedestal_electron = (
            self.data.physics.f_nd_plasma_pedestal_greenwald
            * PlasmaDensityLimit.calculate_greenwald_density_limit(
                c_plasma=self.data.physics.plasma_current,
                rminor=self.data.physics.rminor,
            )
        )
        self.data.physics.nd_plasma_separatrix_electron = (
            self.data.physics.f_nd_plasma_separatrix_greenwald
            * PlasmaDensityLimit.calculate_greenwald_density_limit(
                c_plasma=self.data.physics.plasma_current,
                rminor=self.data.physics.rminor,
            )
        )

set_physics_variables()

Calculates and sets physics variables required for the profile.

Source code in process/models/physics/profiles.py
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
def set_physics_variables(self):
    """Calculates and sets physics variables required for the profile."""
    if (
        PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
        == PlasmaProfileShapeType.PARABOLIC_PROFILE
    ):
        self.data.physics.nd_plasma_electron_on_axis = (
            self.calculate_parabolic_profile_on_axis_density(
                nd_vol_average=self.data.physics.nd_plasma_electrons_vol_avg,
                alphan=self.data.physics.alphan,
            )
        )
    elif (
        PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
        == PlasmaProfileShapeType.PEDESTAL_PROFILE
    ):
        self.data.physics.nd_plasma_electron_on_axis = self.calculate_pedestal_profile_on_axis_density(  # noqa: E501
            radius_plasma_pedestal_density_norm=self.data.physics.radius_plasma_pedestal_density_norm,
            nd_pedestal=self.data.physics.nd_plasma_pedestal_electron,
            nd_separatrix=self.data.physics.nd_plasma_separatrix_electron,
            nd_vol_average=self.data.physics.nd_plasma_electrons_vol_avg,
            alphan=self.data.physics.alphan,
        )
    self.data.physics.nd_plasma_ions_on_axis = (
        self.data.physics.nd_plasma_ions_total_vol_avg
        / self.data.physics.nd_plasma_electrons_vol_avg
        * self.data.physics.nd_plasma_electron_on_axis
    )

ElectronTemperatureProfile

Bases: Profile

Electron temperature (Tₑ) profile class. Contains a function to calculate the temperature profile and store the data.

Source code in process/models/physics/profiles.py
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
class ElectronTemperatureProfile(Profile):
    """Electron temperature (Tₑ) profile class. Contains a function to calculate the
    temperature profile and store the data.
    """

    def run(self):
        """Subroutine to initialise neprofile and execute calculations."""
        super().run()
        self.normalise_profile_x()
        self.calculate_profile_dx()
        self.set_physics_variables()
        self.calculate_profile_y(
            rho=self.profile_x,
            radius_plasma_pedestal_temp_norm=self.data.physics.radius_plasma_pedestal_temp_norm,
            temp_on_axis_kev=self.data.physics.temp_plasma_electron_on_axis_kev,
            temp_pedestal_kev=self.data.physics.temp_plasma_pedestal_kev,
            temp_separatrix_kev=self.data.physics.temp_plasma_separatrix_kev,
            alphat=self.data.physics.alphat,
            tbeta=self.data.physics.tbeta,
        )
        self.integrate_profile_y()

    def calculate_profile_y(
        self,
        rho: np.array,
        radius_plasma_pedestal_temp_norm: float,
        temp_on_axis_kev: float,
        temp_pedestal_kev: float,
        temp_separatrix_kev: float,
        alphat: float,
        tbeta: float,
    ) -> None:
        """Calculates the temperature at each normalised minor radius (ρ) position.

        Parameters
        ----------
        rho :
            Normalised minor radius (ρ) vector
        radius_plasma_pedestal_temp_norm :
            Normalised minor radius pedestal position (ρₜ,pedestal).
        temp_on_axis_kev :
            Central on-axis temperature (T₀) [keV].
        temp_pedestal_kev :
            Pedestal temperature (T_pedestal) [keV].
        temp_separatrix_kev :
            Separatrix temperature (T_separatrix) [keV].
        alphat :
            Temperature peaking parameter (αₜ).
        tbeta :
            Second temperature exponent (βₜ).

        Raises
        ------
        ProcessValueError
            If negative temperature in plasma profile

        """  # noqa: RUF002
        if (
            PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
            == PlasmaProfileShapeType.PARABOLIC_PROFILE
        ):
            # profile values of 0 cause divide by 0 errors so ensure the profile value
            # is at least 1e-8
            # which is small enough that it won't make a difference to any calculations
            self.profile_y = np.maximum(temp_on_axis_kev * (1 - rho**2) ** alphat, 1e-8)
            return

        if temp_on_axis_kev < temp_pedestal_kev:
            logger.info(
                "TPROFILE: Pedestal temperature is higher than core temperature. %s, %s",
                temp_pedestal_kev,
                temp_on_axis_kev,
            )

        rho_index = rho <= radius_plasma_pedestal_temp_norm
        self.profile_y[rho_index] = (
            temp_pedestal_kev
            + (temp_on_axis_kev - temp_pedestal_kev)
            * (1 - (rho[rho_index] / radius_plasma_pedestal_temp_norm) ** tbeta)
            ** alphat
        )
        self.profile_y[~rho_index] = temp_separatrix_kev + (
            temp_pedestal_kev - temp_separatrix_kev
        ) * (1 - rho[~rho_index]) / (1 - radius_plasma_pedestal_temp_norm)

        # Check for any negative temperature in profile: always fatal in
        # later models eventually
        if (self.profile_y < 0).any():
            raise ProcessValueError("Negative temperature in plasma profile")

    @staticmethod
    def calculate_pedestal_profile_on_axis_temperature(
        radius_plasma_pedestal_temp_norm: float,
        temp_pedestal_kev: float,
        temp_separatrix_kev: float,
        temp_vol_avg_kev: float,
        alphat: float,
        tbeta: float,
    ) -> float:
        """Calculates the core on-axis temperature (T₀) of a pedestalised profile.

        Parameters
        ----------
        radius_plasma_pedestal_temp_norm :
            Normalised minor radius pedestal position (ρₜ,pedestal).
        temp_pedestal_kev :
            Pedestal temperature (T_pedestal) [keV].
        temp_separatrix_kev :
            Separatrix temperature (T_separatrix) [keV].
        temp_vol_avg_kev :
            Volume average temperature (⟨T⟩) [keV].
        alphat :
            Temperature peaking parameter (αₜ).
        tbeta :
            Second temperature exponent (βₜ).

        Returns
        -------
        :
            The core on-axis temperature (T₀) [keV]

        """
        #  Calculate core temperature

        return temp_pedestal_kev + (
            (
                tbeta
                * (
                    3 * temp_vol_avg_kev
                    + temp_separatrix_kev
                    * (
                        -2.0
                        + radius_plasma_pedestal_temp_norm
                        + radius_plasma_pedestal_temp_norm**2
                    )
                    - temp_pedestal_kev
                    * (
                        1
                        + radius_plasma_pedestal_temp_norm
                        + radius_plasma_pedestal_temp_norm**2
                    )
                )
            )
            / (
                6
                * radius_plasma_pedestal_temp_norm**2
                * sp.special.beta(1 + alphat, 2 / tbeta)
            )
        )

    @staticmethod
    def calculate_parabolic_profile_on_axis_temperature(
        temp_vol_avg_kev: float,
        alphat: float,
    ) -> float:
        """Calculates the core on-axis temperature (T₀) of a parabolic profile.

        Parameters
        ----------
        temp_vol_avg_kev :
            Volume average temperature (⟨T⟩) [keV].
        alphat :
            Temperature peaking parameter (αₜ).

        Returns
        -------
        :
            The core on-axis temperature (T₀) [keV]

        """
        #  Calculate core temperature

        return temp_vol_avg_kev * (1.0 + alphat)

    def set_physics_variables(self):
        """Calculates and sets physics variables required for the temperature profile."""
        if (
            PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
            == PlasmaProfileShapeType.PARABOLIC_PROFILE
        ):
            self.data.physics.temp_plasma_electron_on_axis_kev = (
                self.calculate_parabolic_profile_on_axis_temperature(
                    temp_vol_avg_kev=self.data.physics.temp_plasma_electron_vol_avg_kev,
                    alphat=self.data.physics.alphat,
                )
            )
        elif (
            PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
            == PlasmaProfileShapeType.PEDESTAL_PROFILE
        ):
            self.data.physics.temp_plasma_electron_on_axis_kev = self.calculate_pedestal_profile_on_axis_temperature(  # noqa: E501
                radius_plasma_pedestal_temp_norm=self.data.physics.radius_plasma_pedestal_temp_norm,
                temp_pedestal_kev=self.data.physics.temp_plasma_pedestal_kev,
                temp_separatrix_kev=self.data.physics.temp_plasma_separatrix_kev,
                temp_vol_avg_kev=self.data.physics.temp_plasma_electron_vol_avg_kev,
                alphat=self.data.physics.alphat,
                tbeta=self.data.physics.tbeta,
            )

        self.data.physics.temp_plasma_ion_on_axis_kev = (
            self.data.physics.temp_plasma_ion_vol_avg_kev
            / self.data.physics.temp_plasma_electron_vol_avg_kev
            * self.data.physics.temp_plasma_electron_on_axis_kev
        )

run()

Subroutine to initialise neprofile and execute calculations.

Source code in process/models/physics/profiles.py
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
def run(self):
    """Subroutine to initialise neprofile and execute calculations."""
    super().run()
    self.normalise_profile_x()
    self.calculate_profile_dx()
    self.set_physics_variables()
    self.calculate_profile_y(
        rho=self.profile_x,
        radius_plasma_pedestal_temp_norm=self.data.physics.radius_plasma_pedestal_temp_norm,
        temp_on_axis_kev=self.data.physics.temp_plasma_electron_on_axis_kev,
        temp_pedestal_kev=self.data.physics.temp_plasma_pedestal_kev,
        temp_separatrix_kev=self.data.physics.temp_plasma_separatrix_kev,
        alphat=self.data.physics.alphat,
        tbeta=self.data.physics.tbeta,
    )
    self.integrate_profile_y()

calculate_profile_y(rho, radius_plasma_pedestal_temp_norm, temp_on_axis_kev, temp_pedestal_kev, temp_separatrix_kev, alphat, tbeta)

Calculates the temperature at each normalised minor radius (ρ) position.

Parameters:

Name Type Description Default
rho array

Normalised minor radius (ρ) vector

required
radius_plasma_pedestal_temp_norm float

Normalised minor radius pedestal position (ρₜ,pedestal).

required
temp_on_axis_kev float

Central on-axis temperature (T₀) [keV].

required
temp_pedestal_kev float

Pedestal temperature (T_pedestal) [keV].

required
temp_separatrix_kev float

Separatrix temperature (T_separatrix) [keV].

required
alphat float

Temperature peaking parameter (αₜ).

required
tbeta float

Second temperature exponent (βₜ).

required

Raises:

Type Description
ProcessValueError

If negative temperature in plasma profile

Source code in process/models/physics/profiles.py
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
def calculate_profile_y(
    self,
    rho: np.array,
    radius_plasma_pedestal_temp_norm: float,
    temp_on_axis_kev: float,
    temp_pedestal_kev: float,
    temp_separatrix_kev: float,
    alphat: float,
    tbeta: float,
) -> None:
    """Calculates the temperature at each normalised minor radius (ρ) position.

    Parameters
    ----------
    rho :
        Normalised minor radius (ρ) vector
    radius_plasma_pedestal_temp_norm :
        Normalised minor radius pedestal position (ρₜ,pedestal).
    temp_on_axis_kev :
        Central on-axis temperature (T₀) [keV].
    temp_pedestal_kev :
        Pedestal temperature (T_pedestal) [keV].
    temp_separatrix_kev :
        Separatrix temperature (T_separatrix) [keV].
    alphat :
        Temperature peaking parameter (αₜ).
    tbeta :
        Second temperature exponent (βₜ).

    Raises
    ------
    ProcessValueError
        If negative temperature in plasma profile

    """  # noqa: RUF002
    if (
        PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
        == PlasmaProfileShapeType.PARABOLIC_PROFILE
    ):
        # profile values of 0 cause divide by 0 errors so ensure the profile value
        # is at least 1e-8
        # which is small enough that it won't make a difference to any calculations
        self.profile_y = np.maximum(temp_on_axis_kev * (1 - rho**2) ** alphat, 1e-8)
        return

    if temp_on_axis_kev < temp_pedestal_kev:
        logger.info(
            "TPROFILE: Pedestal temperature is higher than core temperature. %s, %s",
            temp_pedestal_kev,
            temp_on_axis_kev,
        )

    rho_index = rho <= radius_plasma_pedestal_temp_norm
    self.profile_y[rho_index] = (
        temp_pedestal_kev
        + (temp_on_axis_kev - temp_pedestal_kev)
        * (1 - (rho[rho_index] / radius_plasma_pedestal_temp_norm) ** tbeta)
        ** alphat
    )
    self.profile_y[~rho_index] = temp_separatrix_kev + (
        temp_pedestal_kev - temp_separatrix_kev
    ) * (1 - rho[~rho_index]) / (1 - radius_plasma_pedestal_temp_norm)

    # Check for any negative temperature in profile: always fatal in
    # later models eventually
    if (self.profile_y < 0).any():
        raise ProcessValueError("Negative temperature in plasma profile")

calculate_pedestal_profile_on_axis_temperature(radius_plasma_pedestal_temp_norm, temp_pedestal_kev, temp_separatrix_kev, temp_vol_avg_kev, alphat, tbeta) staticmethod

Calculates the core on-axis temperature (T₀) of a pedestalised profile.

Parameters:

Name Type Description Default
radius_plasma_pedestal_temp_norm float

Normalised minor radius pedestal position (ρₜ,pedestal).

required
temp_pedestal_kev float

Pedestal temperature (T_pedestal) [keV].

required
temp_separatrix_kev float

Separatrix temperature (T_separatrix) [keV].

required
temp_vol_avg_kev float

Volume average temperature (⟨T⟩) [keV].

required
alphat float

Temperature peaking parameter (αₜ).

required
tbeta float

Second temperature exponent (βₜ).

required

Returns:

Type Description
float

The core on-axis temperature (T₀) [keV]

Source code in process/models/physics/profiles.py
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
@staticmethod
def calculate_pedestal_profile_on_axis_temperature(
    radius_plasma_pedestal_temp_norm: float,
    temp_pedestal_kev: float,
    temp_separatrix_kev: float,
    temp_vol_avg_kev: float,
    alphat: float,
    tbeta: float,
) -> float:
    """Calculates the core on-axis temperature (T₀) of a pedestalised profile.

    Parameters
    ----------
    radius_plasma_pedestal_temp_norm :
        Normalised minor radius pedestal position (ρₜ,pedestal).
    temp_pedestal_kev :
        Pedestal temperature (T_pedestal) [keV].
    temp_separatrix_kev :
        Separatrix temperature (T_separatrix) [keV].
    temp_vol_avg_kev :
        Volume average temperature (⟨T⟩) [keV].
    alphat :
        Temperature peaking parameter (αₜ).
    tbeta :
        Second temperature exponent (βₜ).

    Returns
    -------
    :
        The core on-axis temperature (T₀) [keV]

    """
    #  Calculate core temperature

    return temp_pedestal_kev + (
        (
            tbeta
            * (
                3 * temp_vol_avg_kev
                + temp_separatrix_kev
                * (
                    -2.0
                    + radius_plasma_pedestal_temp_norm
                    + radius_plasma_pedestal_temp_norm**2
                )
                - temp_pedestal_kev
                * (
                    1
                    + radius_plasma_pedestal_temp_norm
                    + radius_plasma_pedestal_temp_norm**2
                )
            )
        )
        / (
            6
            * radius_plasma_pedestal_temp_norm**2
            * sp.special.beta(1 + alphat, 2 / tbeta)
        )
    )

calculate_parabolic_profile_on_axis_temperature(temp_vol_avg_kev, alphat) staticmethod

Calculates the core on-axis temperature (T₀) of a parabolic profile.

Parameters:

Name Type Description Default
temp_vol_avg_kev float

Volume average temperature (⟨T⟩) [keV].

required
alphat float

Temperature peaking parameter (αₜ).

required

Returns:

Type Description
float

The core on-axis temperature (T₀) [keV]

Source code in process/models/physics/profiles.py
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
@staticmethod
def calculate_parabolic_profile_on_axis_temperature(
    temp_vol_avg_kev: float,
    alphat: float,
) -> float:
    """Calculates the core on-axis temperature (T₀) of a parabolic profile.

    Parameters
    ----------
    temp_vol_avg_kev :
        Volume average temperature (⟨T⟩) [keV].
    alphat :
        Temperature peaking parameter (αₜ).

    Returns
    -------
    :
        The core on-axis temperature (T₀) [keV]

    """
    #  Calculate core temperature

    return temp_vol_avg_kev * (1.0 + alphat)

set_physics_variables()

Calculates and sets physics variables required for the temperature profile.

Source code in process/models/physics/profiles.py
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
def set_physics_variables(self):
    """Calculates and sets physics variables required for the temperature profile."""
    if (
        PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
        == PlasmaProfileShapeType.PARABOLIC_PROFILE
    ):
        self.data.physics.temp_plasma_electron_on_axis_kev = (
            self.calculate_parabolic_profile_on_axis_temperature(
                temp_vol_avg_kev=self.data.physics.temp_plasma_electron_vol_avg_kev,
                alphat=self.data.physics.alphat,
            )
        )
    elif (
        PlasmaProfileShapeType(self.data.physics.i_plasma_pedestal)
        == PlasmaProfileShapeType.PEDESTAL_PROFILE
    ):
        self.data.physics.temp_plasma_electron_on_axis_kev = self.calculate_pedestal_profile_on_axis_temperature(  # noqa: E501
            radius_plasma_pedestal_temp_norm=self.data.physics.radius_plasma_pedestal_temp_norm,
            temp_pedestal_kev=self.data.physics.temp_plasma_pedestal_kev,
            temp_separatrix_kev=self.data.physics.temp_plasma_separatrix_kev,
            temp_vol_avg_kev=self.data.physics.temp_plasma_electron_vol_avg_kev,
            alphat=self.data.physics.alphat,
            tbeta=self.data.physics.tbeta,
        )

    self.data.physics.temp_plasma_ion_on_axis_kev = (
        self.data.physics.temp_plasma_ion_vol_avg_kev
        / self.data.physics.temp_plasma_electron_vol_avg_kev
        * self.data.physics.temp_plasma_electron_on_axis_kev
    )

calculate_vol_avg_of_profile(profile_x, profile_y, profile_dx=None)

Calculate the volume averaged value (⟨profile_y⟩) of a radially normalised profile.

Parameters:

Name Type Description Default
profile_x ndarray

The x-values of the profile.

required
profile_y ndarray

The y-values of the profile.

required
profile_dx float | None

The spacing between consecutive x-values in the profile.

None

Returns:

Type Description
float

The volume-averaged value (⟨profile_y⟩) of the profile.

Raises:

Type Description
ValueError

If profile_x is not a 1D array, contains fewer than 2 points, does not span from 0 to 1,or is not strictly increasing.

Notes
  • The 2 factor in the calculation arises using both sides of the profile and the radial normalisation of the profile.
Source code in process/models/physics/profiles.py
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
def calculate_vol_avg_of_profile(
    profile_x: np.ndarray, profile_y: np.ndarray, profile_dx: float | None = None
) -> float:
    """Calculate the volume averaged value (⟨profile_y⟩) of a radially normalised
    profile.

    Parameters
    ----------
    profile_x :
        The x-values of the profile.
    profile_y :
        The y-values of the profile.
    profile_dx :
        The spacing between consecutive x-values in the profile.

    Returns
    -------
    float
        The volume-averaged value (⟨profile_y⟩) of the profile.

    Raises
    ------
    ValueError
        If profile_x is not a 1D array, contains fewer than 2 points,
        does not span from 0 to 1,or is not strictly increasing.


    Notes
    -----
    - The 2 factor in the calculation arises using both sides of the profile and the
    radial normalisation of the profile.


    """
    if profile_x.ndim != 1:
        raise ValueError("profile_x must be a 1D array.")

    if profile_x.size < 2:
        raise ValueError("profile_x must contain at least 2 points.")

    if not np.isclose(profile_x[0], 0.0) or not np.isclose(profile_x[-1], 1.0):
        raise ValueError("profile_x must span from 0 to 1.")

    if np.any(np.diff(profile_x) <= 0):
        raise ValueError("profile_x must be strictly increasing.")

    return 2.0 * sp.integrate.simpson(
        profile_y * profile_x,
        x=profile_x,
        dx=profile_dx if profile_dx is not None else profile_x[1] - profile_x[0],
    )