XGCa
|
#include "petsc_version_defs.h"
#include <petsc/finclude/petsc.h>
Data Types | |
module | f90moduleinterfaces |
Explicit interfaces to somve PETSc function used by the FSA solver. More... | |
interface | f90moduleinterfaces::pcsetapplicationcontext |
interface | f90moduleinterfaces::pcgetapplicationcontext |
Macros | |
#define | DPOT(inode, iphi) psn%dpot(inode,1) |
#define | IDEN(inode, iphi) psn%idensity0(inode) |
#define | EDEN(inode, iphi) psn%edensity0(inode) |
Functions/Subroutines | |
subroutine | init_poisson (grid, psn, iflag_dummy) |
subroutine | psn_init_poisson_private (psn, grid, solver_data, nrhs) |
subroutine | pcshermanmorrisonapply (pc, xin, yout, ierr) |
Applies a Sherman-Morrison preconditioner. More... | |
subroutine | create_2field_solver (grid, bd, solver, solver_data) |
Sets up a 2x2 block solver for the Poisson equation. The first row is the Poisson equation, the second row is a constraint equation for the flux-surface averaged potential. More... | |
subroutine | make_12mat (grid, Bmat, BIdentMat, FSAMass_Te, solver, bd, scale, xgc_petsc) |
subroutine | make_21mat (grid, Cmat, Cio, Bmat, CConstBCVec, bd, xgc_petsc, petsc_xgc_bd) |
subroutine | solve_poisson (grid, psn, n0_only, iflag) |
High level Poisson solver routine, offers choice of outer iterative solver, two-step solver, and FSA solver. More... | |
subroutine | set_rhs_poisson (grid, psn, iflag) |
Routine to set and process the rhs arrays for both axisymmetric and turbulent poisson solvers. Also sets boundary data for XGCA. More... | |
subroutine | zero_out_wall (var) |
subroutine | solve_axisymmetric_poisson_iter (grid, psn, iflag) |
Poisson solver that uses an iterative method to solve for the axisymmetric electrostatic potential: The Poisson equation for the axisymmetric mode is
\[ -\nabla_\perp \cdot \xi \nabla_\perp \overline{\phi}_{k+1} + \frac{e n_0/T_{e,0}} \overline{\phi}_{k+1} = e \left( \langle \overline{\delta n_i} \rangle_g - \overline{\delta n_{e,NA}} \right) + \frac{e n_0/T_{e,0}} \langle \phi_{k} \rangle_{fs} \] where \(k\) is the iteration index, \(xi\) is the electric susceptibility, \(e\) is the elementary charge, \(n_0\) and \(T_0\) are the background density and temperature, \(\langle \overline{\delta n_i} \rangle_g\) is the gyroaveraged ion (gyrocenter) charge density, \(\overline{\delta n_{e,NA}}\) is the non-adiabatic electron charge density, and \(\langle\dots\rangle\) is the the flux-surface average. \(\overline{\dots}\) is the toroidal average. More... | |
subroutine | post_process_axisymmetric_poisson (grid, psn, n0_only) |
Routine to post process the solution from the axisymmetric poisson solver. More... | |
subroutine | solve_turb_poisson (grid, psn) |
Routine to solve the non-axisymmetric Poisson equation:
\[ -\nabla_\perp \cdot \xi \nabla_\perp \delta\phi + \frac{e n_0/T_{e,0}} \delta\phi = e \left( \langle \delta n_i \rangle_g - \overline{\delta n_{e,NA}} \right) \] . More... | |
subroutine | post_process_turb_poisson (grid, psn) |
Routine to post process the solution from the non-axisymmetric Poisson equation. More... | |
subroutine | solve_poisson_private (grid, psn, iflag) |
subroutine | update_poisson_solver (grid, psn, solver_data) |
Updates the LHS matrix of the Poisson equation (density and electron temperature) More... | |
subroutine | ptb_3db_da_solver_init (grid) |
We solve perturbed Ampere's law in the following form: dB = - grad(phi).grad(dA_phi) = grad(phi).grad(dpsi) –> grad(phi).curl(dB) = - grad(phi).curl(grad(phi).grad(dA_phi)) = div.(grad(phi) x (grad(phi) x grad(dA_phi)) = -div.[g^{phi,phi} * (grad(dA_phi) - d(dA_phi)/dphi grad(phi))] = -(1/R) * [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] = mu_0 j.grad(phi) –> [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] = -mu_0 R j^phi = -mu_0 (B_phi/B) j_parallel –> We need to set up one Helmholtz solver with alpha=1/R and beta=0 –> The RHS is -mu_0*(B_phi/B)*j_parallel –> j_parallel = psnijpar + psnejpar. More... | |
subroutine | ptb_3db_div_clean_init (grid) |
Initializes a solver for divergence cleaning of the perturbed magnetic field This equation is being solved: div.(R*phi_b) - (ntor^2/R)*phi_b = R*div.delta_B. More... | |
subroutine | ptb_3db_create_1field_solver (grid, solver, prefix) |
create 1 field solver with A00 built (Amat), just setoperator & setup like 2 field More... | |
#define DPOT | ( | inode, | |
iphi | |||
) | psn%dpot(inode,1) |
#define EDEN | ( | inode, | |
iphi | |||
) | psn%edensity0(inode) |
#define IDEN | ( | inode, | |
iphi | |||
) | psn%idensity0(inode) |
subroutine create_2field_solver | ( | type(grid_type) | grid, |
type(boundary2_type) | bd, | ||
type(xgc_solver) | solver, | ||
type(solver_init_data), intent(in) | solver_data | ||
) |
Sets up a 2x2 block solver for the Poisson equation. The first row is the Poisson equation, the second row is a constraint equation for the flux-surface averaged potential.
[in] | grid | XGC grid object, type(grid_type) |
[in] | bd | Object containing the indices of the solver boundary vertices in the XGC mesh. |
[in,out] | solver | XGC solver object into which to put the new solver, type(xgc_solver) |
[in] | solver_data | Data needed for settug up the ssolver, type(solver_init_data) |
subroutine init_poisson | ( | type(grid_type) | grid, |
type(psn_type) | psn, | ||
integer | iflag_dummy | ||
) |
subroutine make_12mat | ( | type(grid_type) | grid, |
intent(inout) | Bmat, | ||
BIdentMat, | |||
intent(in) | FSAMass_Te, | ||
type(xgc_solver) | solver, | ||
type(boundary2_type), intent(in) | bd, | ||
real (8), intent(in) | scale, | ||
dimension(grid%nnode), intent(in) | xgc_petsc | ||
) |
subroutine make_21mat | ( | type(grid_type) | grid, |
Cmat, | |||
Cio, | |||
intent(in) | Bmat, | ||
CConstBCVec, | |||
type(boundary2_type), intent(in) | bd, | ||
dimension(grid%nnode), intent(in) | xgc_petsc, | ||
dimension(grid%nnode), intent(in) | petsc_xgc_bd | ||
) |
subroutine pcshermanmorrisonapply | ( | pc, | |
xin, | |||
yout, | |||
ierr | |||
) |
Applies a Sherman-Morrison preconditioner.
[in,out] | pc | PETCs PC object (manages preconditioners) |
xin | ???, Vec | |
yout | ???, Vec | |
[out] | ierr | PETc error code |
subroutine post_process_axisymmetric_poisson | ( | type(grid_type) | grid, |
type(psn_type) | psn, | ||
logical, intent(in) | n0_only | ||
) |
Routine to post process the solution from the axisymmetric poisson solver.
[in] | grid | XGC grid data object |
[in,out] | psn | XGC field data object |
[in] | n0_only | Flag for solving only n=0 mode and then exit (no function in XGCa), logical |
subroutine post_process_turb_poisson | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(inout) | psn | ||
) |
Routine to post process the solution from the non-axisymmetric Poisson equation.
[in] | grid | XGC grid data object |
[in,out] | psn | XGC field data object |
subroutine psn_init_poisson_private | ( | type(psn_type) | psn, |
type(grid_type) | grid, | ||
type(solver_init_data), intent(in) | solver_data, | ||
integer | nrhs | ||
) |
subroutine ptb_3db_create_1field_solver | ( | type(grid_type), intent(in) | grid, |
type(xgc_solver), intent(inout) | solver, | ||
character(*), intent(in) | prefix | ||
) |
create 1 field solver with A00 built (Amat), just setoperator & setup like 2 field
grid | (in) type(grid_type), grid info |
solver | (inout) type(xgc_solver_type), XGC data structure for the solver |
prefix | (in) char(*), prefix for the solver |
subroutine ptb_3db_da_solver_init | ( | type(grid_type), intent(in) | grid | ) |
We solve perturbed Ampere's law in the following form: dB = - grad(phi).grad(dA_phi) = grad(phi).grad(dpsi) –> grad(phi).curl(dB) = - grad(phi).curl(grad(phi).grad(dA_phi)) = div.(grad(phi) x (grad(phi) x grad(dA_phi)) = -div.[g^{phi,phi} * (grad(dA_phi) - d(dA_phi)/dphi grad(phi))] = -(1/R) * [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] = mu_0 j.grad(phi) –> [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] = -mu_0 R j^phi = -mu_0 (B_phi/B) j_parallel –> We need to set up one Helmholtz solver with alpha=1/R and beta=0 –> The RHS is -mu_0*(B_phi/B)*j_parallel –> j_parallel = psnijpar + psnejpar.
grid | (in) type(grid_type), grid information |
subroutine ptb_3db_div_clean_init | ( | type(grid_type), intent(in) | grid | ) |
Initializes a solver for divergence cleaning of the perturbed magnetic field This equation is being solved: div.(R*phi_b) - (ntor^2/R)*phi_b = R*div.delta_B.
grid | (in) type(grid_type), grid information |
ntor | (in) integer, toroidal mode number |
subroutine set_rhs_poisson | ( | type(grid_type) | grid, |
type(psn_type) | psn, | ||
integer, intent(in) | iflag | ||
) |
Routine to set and process the rhs arrays for both axisymmetric and turbulent poisson solvers. Also sets boundary data for XGCA.
[in] | grid | XGC grid data object |
[in,out] | psn | XGC field data object |
[in] | iflag | Flag to indicate whether adiabatic elec or full eq. is solved, integer |
subroutine solve_axisymmetric_poisson_iter | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(inout) | psn, | ||
integer, intent(in) | iflag | ||
) |
Poisson solver that uses an iterative method to solve for the axisymmetric electrostatic potential: The Poisson equation for the axisymmetric mode is
\[ -\nabla_\perp \cdot \xi \nabla_\perp \overline{\phi}_{k+1} + \frac{e n_0/T_{e,0}} \overline{\phi}_{k+1} = e \left( \langle \overline{\delta n_i} \rangle_g - \overline{\delta n_{e,NA}} \right) + \frac{e n_0/T_{e,0}} \langle \phi_{k} \rangle_{fs} \]
where \(k\) is the iteration index, \(xi\) is the electric susceptibility, \(e\) is the elementary charge, \(n_0\) and \(T_0\) are the background density and temperature, \(\langle \overline{\delta n_i} \rangle_g\) is the gyroaveraged ion (gyrocenter) charge density, \(\overline{\delta n_{e,NA}}\) is the non-adiabatic electron charge density, and \(\langle\dots\rangle\) is the the flux-surface average. \(\overline{\dots}\) is the toroidal average.
The RHS and LHS of the axisymmetric equation can be Fourier-filtered (poloidal modes). The LHS of the non-axisymmetric equation can be Fourier-filtered in the toroidal and poloidal direction. The solver in both cases is a linear finite-element solver on an unstructured triangular grid.
[in] | grid | Solver grid data, type(grid_type) |
[in,out] | psn | Field data; the potential is stored there, type(psn_type) |
[in] | iflag | Flag to indicate whether adiabatic elec or full eq. is solved, integer |
subroutine solve_poisson | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(inout) | psn, | ||
logical, intent(in) | n0_only, | ||
integer, intent(in) | iflag | ||
) |
High level Poisson solver routine, offers choice of outer iterative solver, two-step solver, and FSA solver.
[in] | grid | XGC grid data object |
[in,out] | psn | XGC field data object |
[in] | n0_only | Flag for solving only n=0 mode and then exit (no function in XGCa), logical |
[in] | iflag | Flag to indicate whether adiabatic elec or full eq. is solved, integer |
subroutine solve_poisson_private | ( | type(grid_type) | grid, |
type(psn_type) | psn, | ||
integer, intent(in) | iflag | ||
) |
subroutine solve_turb_poisson | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(inout) | psn | ||
) |
Routine to solve the non-axisymmetric Poisson equation:
\[ -\nabla_\perp \cdot \xi \nabla_\perp \delta\phi + \frac{e n_0/T_{e,0}} \delta\phi = e \left( \langle \delta n_i \rangle_g - \overline{\delta n_{e,NA}} \right) \]
.
[in] | grid | XGC grid data object |
[in,out] | psn | XGC field data object |
subroutine update_poisson_solver | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(inout) | psn, | ||
type(solver_init_data), intent(in) | solver_data | ||
) |
Updates the LHS matrix of the Poisson equation (density and electron temperature)
[in] | grid | XGC grid data object, type(grid_type) |
[in,out] | psn | XGC field data object, type(psn_type) |
[in] | solver_data | Object with profile and magnetic data for solvers |
subroutine zero_out_wall | ( | real (8), dimension(grid%nnode), intent(inout) | var | ) |