Skip to content

Guiding Center Hamiltonian

dexter.GCHamiltonian()

Guiding center Hamiltonian.

The guiding center Hamiltonian is defined as:

$$ \mathcal{H}(\theta, \psi; P_\zeta, \mu) = \dfrac{(P_\zeta+\psi_p)^2B^2}{2g^2} + \mu B. $$

This object has no internal state, it only provides evaluation methods of the Hamiltonian and its derived quantities.

eval(machine: Machine, theta: float, flux: MagneticFlux, pzeta: float | MagneticFlux, mu: float) -> Array1 classmethod

Evaluates the Hamiltonian $\mathcal{H}(\theta, \psi; P_\zeta, \mu)$.

Parameters:

  • machine (Machine) –

    The machine on which the Hamiltonian is defined.

  • theta (float) –

    The $\theta$ value.

  • flux (MagneticFlux) –

    The magnetic flux $\psi$.

  • pzeta (float | MagneticFlux) –

    The $P_\zeta$ constant of motion.

  • mu (float) –

    The $\mu$ constant of motion.

Example
>>> geometry = dex.LarGeometry(baxis=3, raxis=1.75, rlast=0.5)
>>> LCFS = geometry.psit_last
>>> machine = dex.Machine(
...     geometry=geometry,
...     qfactor=dex.ParabolicQfactor(qaxis=1.1, qlast=3.9, lcfs=LCFS),
...     current=dex.LarCurrent(),
...     bfield=dex.LarBfield(),
... )
>>>
>>> psit = dex.MagneticFlux.Toroidal(0.01)
>>> pzeta = -0.1 * machine.psip_last
>>> h = dex.GCHamiltonian.eval(machine, 3.14, psit, pzeta, 1e-4)

eval_array(machine: Machine, theta_array: Array1, flux_array: MagneticFluxArray, pzeta: float | MagneticFlux, mu: float) -> Array1 classmethod

Evaluates the Hamiltonian $\mathcal{H}(\theta, \psi; P_\zeta, \mu)$ over a 1D array.

Note

The produced array has the form:

$$| \mathcal{H}(\theta_0, \psi_0),\ \mathcal{H}(\theta_1, \psi_1),\ \dotsc,\ \mathcal{H}(\theta_N, \psi_N)\ |$$

This method uses the flux_array.kind flux coordinate to perform all evaluations. Since only $\psi_p$ appears explicitly in the Hamiltonian, if flux_array.kind == "Toroidal" then an extra evaluation is needed to calculate $\psi_p$.

Parameters:

  • machine (Machine) –

    The machine on which the Hamiltonian is defined.

  • theta_array (Array1) –

    The $\theta$ values.

  • flux_array (MagneticFluxArray) –

    The flux $\psi$ values.

  • pzeta (float | MagneticFlux) –

    The $P_\zeta$ constant of motion.

  • mu (float) –

    The $\mu$ constant of motion.

Example
>>> geometry = dex.LarGeometry(baxis=3, raxis=1.75, rlast=0.5)
>>> LCFS = geometry.psit_last
>>> machine = dex.Machine.FromNetcdf(path, "Akima", "Bicubic")
>>>
>>> pzeta, mu = -0.1 * machine.psip_last, 1e-4
>>> theta_array = np.linspace(0, 6.28, 10)
>>>
>>> # H(ψt, θ)
>>> psit_array = dex.MagneticFluxArray.Toroidal(np.linspace(0, machine.psit_last.value, 10))
>>> h = dex.GCHamiltonian.eval_array(machine, theta_array, psit_array, pzeta, mu)
>>>
>>> # H(ψp, θ)
>>> psip_array = dex.MagneticFluxArray.Poloidal(np.linspace(0, machine.psip_last.value, 10))
>>> h = dex.GCHamiltonian.eval_array(machine, theta_array, psip_array, pzeta, mu)

eval_grid(machine: Machine, theta_array: Array1, flux_array: MagneticFluxArray, pzeta: float | MagneticFlux, mu: float) -> Array2 classmethod

Evaluates the Hamiltonian $\mathcal{H}(\theta, \psi; P_\zeta, \mu)$ over a 2D grid.

Note

The $\theta$ and $\psi$ arrays do not need to have the same length. The produced grid has the form:

$$|\ \mathcal{H}(\theta_0, \psi_0), \ \mathcal{H}(\theta_1, \psi_0), \ \dotsc, \ \mathcal{H}(\theta_M, \psi_0)\ | \\ |\ \mathcal{H}(\theta_0, \psi_1), \ \mathcal{H}(\theta_1, \psi_1), \ \dotsc, \ \mathcal{H}(\theta_M, \psi_1)\ | \\ ... \\ |\mathcal{H}(\theta_0, \psi_N), \mathcal{H}(\theta_1, \psi_N), \dotsc, \mathcal{H}(\theta_M, \psi_N) | $$

This method uses the flux_array.kind flux coordinate to perform all evaluations. Since only $\psi_p$ appears explicitly in the Hamiltonian, if flux_array.kind == "Toroidal" then an extra evaluation is needed to calculate $\psi_p$.

Parameters:

  • machine (Machine) –

    The machine on which the Hamiltonian is defined.

  • theta_array (Array1) –

    The $\theta$ values.

  • flux_array (MagneticFluxArray) –

    The flux $\psi$ values.

  • pzeta (float | MagneticFlux) –

    The $P_\zeta$ constant of motion.

  • mu (float) –

    The $\mu$ constant of motion.

Example
>>> geometry = dex.LarGeometry(baxis=3, raxis=1.75, rlast=0.5)
>>> LCFS = geometry.psit_last
>>> machine = dex.Machine.FromNetcdf(path, "Akima", "Bicubic")
>>>
>>> pzeta, mu = -0.1 * machine.psip_last, 1e-4
>>> theta_array = np.linspace(0, 6.28, 4)
>>>
>>> # H(ψt, θ)
>>> psit_array = dex.MagneticFluxArray.Toroidal(np.linspace(0, machine.psit_last.value, 5))
>>> h = dex.GCHamiltonian.eval_grid(machine, theta_array, psit_array, pzeta, mu)
>>>
>>> # H(ψp, θ)
>>> psip_array = dex.MagneticFluxArray.Poloidal(np.linspace(0, machine.psip_last.value, 5))
>>> h = dex.GCHamiltonian.eval_grid(machine, theta_array, psip_array, pzeta, mu)

eval_rho_array(machine: Machine, theta_array: Array1, flux_array: MagneticFluxArray, energy: float, mu: float) -> Array1 classmethod

Calculates $\rho_{||}$ from a given energy and $\mu$ value by solving the Hamiltonian w.r.t. $\rho_{||}$.

This method is useful for instantiating particles in Boozer coordinates with equal energy.

Note

  • The computed $\rho_{||}$ values are always non-negative, as they are the result of a square root. To instantiate Counter-Passing particles, the values must be negated manually.

  • If no energetically allowed orbits exist for the given $\psi$, $\theta$ and $E$, the corresponding array element is set to nan.

  • Perturbations are ignored.

Parameters:

  • machine (Machine) –

    The machine on which the Hamiltonian is defined.

  • theta_array (Array1) –

    The $\theta$ values.

  • flux_array (MagneticFluxArray) –

    The flux $\psi$ values.

  • energy (float) –

    The target energy value.

  • mu (float) –

    The $\mu$ constant of motion.

Example
>>> geometry = dex.LarGeometry(1, 1.75, 0.5)
>>> LCFS = geometry.psit_last
>>> machine = dex.Machine(
... geometry=geometry,
... qfactor=dex.ParabolicQfactor(qaxis=1.1, qlast=3.5, lcfs=LCFS),
... current=dex.LarCurrent(),
... bfield=dex.LarBfield(),
... )
>>>
>>> num = 10
>>> psit0_values = np.linspace(1e-6, LCFS.value, num)
>>> psit0s = dex.MagneticFluxArray.Toroidal(psit0_values)
>>> theta0s = np.zeros(num)
>>>
>>> # ======================
>>>
>>> mu = 5e-5
>>> energy = machine.quantity(21, "keV").to("NormJoule").m
>>> rho0s = dex.GCHamiltonian.eval_rho_array(machine, theta0s, psit0s, energy, mu)
>>>
>>> initial_conditions = dex.QueueInitialConditions.Boozer(
... t0=np.zeros(num),
... flux0=psit0s,
... theta0=theta0s,
... zeta0=np.zeros(num),
... rho0=rho0s,
... mu0=np.full(num, mu),
... )
>>>
>>> queue = dex.Queue(initial_conditions)
>>> queue.classify(machine)
>>>
>>> from math import isclose
>>> for particle in queue:
...     energy = machine.quantity(particle.initial_energy, "NormJoule").to("keV")
...     assert isclose(energy.magnitude, 21) # keV

eval_theta_dot_array(machine: Machine, theta_array: Array1, flux_array: MagneticFluxArray, pzeta: float | MagneticFlux, mu: float) -> Array1 classmethod

Evaluates $\dot\theta = \partial\mathcal{H}(\theta, \psi; P_\zeta, \mu)/\partial P_\theta$ over a 1D array

This method uses the flux_array.kind flux coordinate to perform all evaluations. Since only $\psi_p$ appears explicitly in the equation (inside the $\rho_{||}$ term), if flux_array.kind == "Toroidal" then an extra evaluation is needed to calculate $\psi_p$.

Parameters:

  • machine (Machine) –

    The machine on which the Hamiltonian is defined.

  • theta_array (Array1) –

    The $\theta$ values.

  • flux_array (MagneticFluxArray) –

    The flux $\psi$ values.

  • pzeta (float | MagneticFlux) –

    The $P_\zeta$ constant of motion.

  • mu (float) –

    The $\mu$ constant of motion.

Example
>>> geometry = dex.LarGeometry(baxis=3, raxis=1.75, rlast=0.5)
>>> LCFS = geometry.psit_last
>>> machine = dex.Machine.FromNetcdf(path, "Akima", "Bicubic")
>>>
>>> pzeta, mu = -0.1 * machine.psip_last, 1e-4
>>> theta_array = np.linspace(0, 6.28, 10)
>>>
>>> # dot(θ(ψt, θ))
>>> psit_array = dex.MagneticFluxArray.Toroidal(np.linspace(0, machine.psit_last.value, 10))
>>> theta_dot = dex.GCHamiltonian.eval_theta_dot_array(
...     machine,
...     theta_array,
...     psit_array,
...     pzeta,
...     mu
... )
>>>
>>> # dot(θ(ψp, θ))
>>> psip_array = dex.MagneticFluxArray.Poloidal(np.linspace(0, machine.psip_last.value, 10))
>>> theta_dot = dex.GCHamiltonian.eval_theta_dot_array(
...     machine,
...     theta_array,
...     psip_array,
...     pzeta,
...     mu
... )

dexter.CriticalPoint

A critical point of the Hamiltonian.

A simple container type that is created through EnergyPzetaPlane.calculate_critical_points.

Attributes:

  • theta (float) –

    The point's $\theta$ position.

  • flux (MagneticFlux) –

    The point's $\psi$ position.

  • energy (float) –

    The point's energy.

  • pzeta (float) –

    The point's $P_\zeta$.

  • mu (float) –

    The point's $\mu$.

  • kind (CriticalPointKind) –

    The point's kind.

dexter.BoundaryExtremum

A boundary extremum of the Hamiltonian.

A simple container type that is created through EnergyPzetaPlane.calculate_boundary_extrema.

Attributes:

  • theta (float) –

    The point's $\theta$ position.

  • flux (MagneticFlux) –

    The point's $\psip$ position.

  • energy (float) –

    The point's energy.

  • pzeta (float) –

    The point's $P_\zeta$.

  • mu (float) –

    The point's $\mu$.

  • kind (BoundaryExtremumKind) –

    The point's kind.