PB3D
[2.45]
Ideal linear high-n MHD stability in 3-D
|
Calculates the normal ( \(\cdot_n\)) or geodesic ( \(\cdot_g\)) components of the plasma perturbation \(\vec{\xi}\) or the magnetic perturbation \(\vec{Q} = \nabla \times \left(\vec{x} \times \vec{B}\right)\). More...
Public Member Functions | |
integer function | calc_xuq_ind (mds_X, mds_sol, grid_eq, grid_X, grid_sol, eq_1, eq_2, X, sol, X_id, XUQ_style, time, XUQ, deriv) |
(time) individual version More... | |
integer function | calc_xuq_arr (mds_X, mds_sol, grid_eq, grid_X, grid_sol, eq_1, eq_2, X, sol, X_id, XUQ_style, time, XUQ, deriv) |
(time) array version More... | |
Calculates the normal ( \(\cdot_n\)) or geodesic ( \(\cdot_g\)) components of the plasma perturbation \(\vec{\xi}\) or the magnetic perturbation \(\vec{Q} = \nabla \times \left(\vec{x} \times \vec{B}\right)\).
Which variables are plot depends on XUQ_style
[17] :
XUQ_style
= 1: \(X\) (supports parallel derivative)XUQ_style
= 2: \(U\) (supports parallel derivative)XUQ_style
= 3: \(Q_n\)XUQ_style
= 4: \(Q_g\)\(X\) and \(U\) support the calculation of the parallel derivative.
For \(Q_n\) and \(Q_g\), the metric variables have to be provided as well.
The perturbation grid is assumed to have the same angular coordinates as the equilibrium grid, and the normal coordinates correspond to either the equilibrium grid (X_grid_style
1), the solution grid (X_grid_style
2) or the enriched equilibrium grid (X_grid_style
3).
The output grid, furthermore, has the angular part corresponding to the equilibrium grid, and the normal part to the solution grid.
Definition at line 56 of file sol_utilities.f90.
integer function sol_utilities::calc_xuq::calc_xuq_arr | ( | type(modes_type), intent(in), target | mds_X, |
type(modes_type), intent(in), target | mds_sol, | ||
type(grid_type), intent(in) | grid_eq, | ||
type(grid_type), intent(in) | grid_X, | ||
type(grid_type), intent(in) | grid_sol, | ||
type(eq_1_type), intent(in) | eq_1, | ||
type(eq_2_type), intent(in) | eq_2, | ||
type(x_1_type), intent(in), target | X, | ||
type(sol_type), intent(in) | sol, | ||
integer, intent(in) | X_id, | ||
integer, intent(in) | XUQ_style, | ||
real(dp), dimension(:), intent(in) | time, | ||
complex(dp), dimension(:,:,:,:), intent(inout) | XUQ, | ||
logical, intent(in), optional | deriv | ||
) |
(time) array version
[in] | mds_x | general modes variables in perturbation grid |
[in] | mds_sol | general modes variables in solution grid |
[in] | grid_eq | equilibrium grid |
[in] | grid_x | perturbation grid |
[in] | grid_sol | solution grid |
[in] | eq_1 | flux equilibrium |
[in] | eq_2 | metric equilibrium |
[in] | x | perturbation variables |
[in] | sol | solution variables |
[in] | x_id | nr. of Eigenvalue |
[in] | xuq_style | whether to calculate \(X\), \(U\), \(Q_n\) or \(Q_g\) |
[in] | time | time range |
[in,out] | xuq | \(X\), \(U\), \(Q_n\) or \(Q_g\) |
[in] | deriv | return parallel derivative |
Definition at line 67 of file sol_utilities.f90.
integer function sol_utilities::calc_xuq::calc_xuq_ind | ( | type(modes_type), intent(in) | mds_X, |
type(modes_type), intent(in) | mds_sol, | ||
type(grid_type), intent(in) | grid_eq, | ||
type(grid_type), intent(in) | grid_X, | ||
type(grid_type), intent(in) | grid_sol, | ||
type(eq_1_type), intent(in) | eq_1, | ||
type(eq_2_type), intent(in) | eq_2, | ||
type(x_1_type), intent(in) | X, | ||
type(sol_type), intent(in) | sol, | ||
integer, intent(in) | X_id, | ||
integer, intent(in) | XUQ_style, | ||
real(dp), intent(in) | time, | ||
complex(dp), dimension(:,:,:), intent(inout) | XUQ, | ||
logical, intent(in), optional | deriv | ||
) |
(time) individual version
[in] | mds_x | general modes variables in perturbation grid |
[in] | mds_sol | general modes variables in solution grid |
[in] | grid_eq | equilibrium grid |
[in] | grid_x | perturbation grid |
[in] | grid_sol | solution grid |
[in] | eq_1 | flux equilibrium |
[in] | eq_2 | metric equilibrium |
[in] | x | perturbation variables |
[in] | sol | solution variables |
[in] | x_id | nr. of Eigenvalue |
[in] | xuq_style | whether to calculate X, U, Qn or Qg |
[in,out] | xuq | normal component of perturbation |
[in] | deriv | return parallel derivative |
Definition at line 520 of file sol_utilities.f90.