Skip to content

adam_prism_common_object

Source: src/app/prism/common/adam_prism_common_object.F90

Dependencies

Contents

Derived Types

prism_common_object

Maxwell equations system class definition, common data to all backends.

Inheritance

Extends: realm_object

Components

NameTypeAttributesDescription
iotype(io_object)IO handler.
amrtype(amr_object)AMR marker handler.
slicestype(slices_object)Slices handler.
blanesmoantype(blanesmoan_object)Blanes-Moan integrator.
cfmtype(cfm_object)Commutator-Free Magnus integrator.
leapfrogtype(leapfrog_object)Leapfrog integrator.
flailtype(flail_object)Linear algebra methods handler.
wenotype(weno_object)WENO reconstructor.
ibtype(ib_object)Immersed boundary.
rktype(rk_object)Runge-Kutta integrator.
adamtype(adam_object)ADAM, grid + tree + field + maps orchestrator.
fdv_schemecharacter(len=:)allocatableFDV scheme, fd/fv.
fdv_orderinteger(kind=I4P)Order of finite difference/volume schemes, general order.
fdv_half_stencilinteger(kind=I4P)Half stencil length of finite difference/volume schemes.
fdv_half_stencilsinteger(kind=I4P)Half stencil length of fdv schemes for each derivative up to 6.
ngcinteger(kind=I4P)pointerNumber of ghost cells.
niinteger(kind=I4P)pointerNumber of cells in i direction.
njinteger(kind=I4P)pointerNumber of cells in j direction.
nkinteger(kind=I4P)pointerNumber of cells in k direction.
nbinteger(kind=I4P)pointerTotal blocks number for MPI.
blocks_numberinteger(kind=I4P)pointerActual blocks number.
nvinteger(kind=I4P)pointerNumber of variables in q vector.
realm_indexinteger(kind=I4P)Inter-realm seam coupling — currently-published stage buffer key.
stage_activeinteger(kind=I4P)Procedure pointer TBPs for FDV operators (set at initialization by backend).
compute_block_total_variationprocedure(compute_block_total_variation_interface)pass(self), pointerCompute TV.
compute_curlprocedure(compute_curl_interface)pass(self), pointerCompute curl.
compute_derivative1procedure(compute_derivative1_interface)pass(self), pointerCompute deriv1.
compute_derivative2procedure(compute_derivative2_interface)pass(self), pointerCompute deriv2.
compute_derivative4procedure(compute_derivative4_interface)pass(self), pointerCompute deriv4.
compute_divergenceprocedure(compute_divergence_interface)pass(self), pointerCompute dive.
compute_gradientprocedure(compute_gradient_interface)pass(self), pointerCompute grad.
compute_laplacianprocedure(compute_laplacian_interface)pass(self), pointerCompute laplac.
nv_cinteger(kind=I4P)pointerNumber of conservative variables in q vector.
nv_sinteger(kind=I4P)pointerNumber of source variables in q vector.
nv_clinteger(kind=I4P)pointerNumber of divergence cleaning variables in q vector.
qreal(kind=R8P)allocatableConservative cell centered variables.
dqreal(kind=R8P)allocatableResiduals right hand side.
q_picreal(kind=R8P)allocatablePIC variables.
pic_fieldsreal(kind=R8P)allocatableFields value at particle locations.
curlreal(kind=R8P)allocatableCurl fields.
divergencereal(kind=R8P)allocatableDivergence fields.
q_nametype(string)allocatableFields names [1:nv].
dq_nametype(string)allocatableResiduals names [1:nv].
q_pic_nametype(string)allocatablePIC Fields names.
curl_nametype(string)allocatableCurl fields names.
div_nametype(string)allocatableDivergence fields names.
energy_Dreal(kind=R8P)allocatableEnergy of field D, time history.
energy_Breal(kind=R8P)allocatableEnergy of field B, time history.
coil_powerreal(kind=R8P)allocatablePower of coils, time history.
Poynting_fluxreal(kind=R8P)allocatableTotal Poynting flux from boundary, time history.
rms_energy_error_Dreal(kind=R8P)RMS energy error of D field.
rms_energy_error_Breal(kind=R8P)RMS energy error of B field.
max_divergence_Dreal(kind=R8P)Maximum of divergence of D field.
max_divergence_Breal(kind=R8P)Maximum of divergence of B field.
max_divergence_Jreal(kind=R8P)Maximum of divergence of J field.
bctype(prism_bc_object)Boundary conditions.
coiltype(prism_coil_object)Coil source term.
external_fieldstype(prism_external_fields_object)External fields.
fWLayertype(prism_fWLayer_object)Far-field weighted layer.
ictype(prism_ic_object)Initial conditions.
leapfrog_pictype(prism_leapfrog_pic_object)Leapfrog PIC integrator.
numericstype(prism_numerics_object)Numerics configuration.
particle_injectiontype(prism_particle_injection_object)Particle injection.
physicstype(prism_physics_object)Physics configuration.
pictype(prism_pic_object)PIC state.
rk_bctype(prism_rk_bc_object)Runge-Kutta BC handler.
rk_pictype(prism_rk_pic_object)Runge-Kutta PIC handler.
timetype(prism_time_object)Time integration state.

Type-Bound Procedures

NameAttributesDescription
load_fdv_from_filepass(self)Load FDV config from file.
initialize_forestpass(self)Invoked by forest%initialize per realm at startup.
compute_local_dt_forestpass(self)Invoked by forest%compute_global_dt during the min reduction.
advance_one_step_forestpass(self)Invoked by forest%evolve_one_step per realm per timestep.
stages_per_step_forestpass(self)Number of integrator stages this realm exposes per step.
open_step_forestpass(self)Per-step prologue (multi-realm path).
begin_stage_forestpass(self)Begin one integrator stage on this realm (multi-realm path).
end_stage_forestpass(self)End the stage: residuals + assignment (multi-realm path).
close_step_forestpass(self)Per-step epilogue (multi-realm path).
post_step_forestpass(self)Invoked by forest%post_step per realm per timestep.
is_done_forestpass(self)Invoked by forest%is_done during the termination reduction.
finalize_forestpass(self)Invoked by forest%finalize per realm at shutdown.
finalize_mpi_forestpass(self)Process-global MPI finalize; forest calls it ONCE after all.
fill_seam_from_peer_forestpass(self)Receive-side roundtrip: copy peer's interior into self's ghosts.
after_topology_build_forestpass(self)Backend hook invoked once after the forest builds the seam maps.
apply_reflux_to_stage_forestpass(self)Apply Berger-Colella reflux to self's stage buffer.
close_block_xh5fnopassClose XH5F file block.
close_file_xh5fnopassClose XH5F file.
open_block_xh5fpass(self)Open block file XH5F.
open_file_xh5fpass(self)Open file XH5F.
save_q_xh5fpass(self)Save in XH5F (XDMF/HDF5) format.
allocate_commonpass(self)Allocate common data.
compute_auxiliary_fieldspass(self)Compute auxiliary fields.
initializepass(self)Initialize the equation common data.
load_restart_filespass(self)Load restart files.
save_energy_errorpass(self)Save energy error history.
save_energy_historypass(self)Save energy history.
save_divergence_historypass(self)Save divergence history.
save_restart_filespass(self)Save restart files.
save_xh5fpass(self)Save simulation data in XH5F format.
initialize_coilspass(self)Initialize coils.
set_rectangular_coil_xpass(self)Subroutine to set a rectangular coil source with +-x normal
set_rectangular_coil_ypass(self)Subroutine to set a rectangular coil source with +-y normal
set_rectangular_coil_zpass(self)Subroutine to set a rectangular coil source with +-z normal
set_circular_coil_xpass(self)Subroutine to set a circular coil source with +-x normal
set_circular_coil_ypass(self)Subroutine to set a circular coil source with +-y normal
set_circular_coil_zpass(self)Subroutine to set a circular coil source with +-z normal
set_solenoid_xpass(self)Subroutine to set a solenoid source with +-x normal
set_solenoid_ypass(self)Subroutine to set a solenoid source with +-y normal
set_solenoid_zpass(self)Subroutine to set a solenoid source with +-z normal
set_helicon_coilpass(self)Subroutine to set a helicon coil source.
coupling_descriptor_forestpass(self)Report (scheme_time, rk_scheme, nv) for β admissibility.

Subroutines

allocate_common

Allocate common data.

fortran
subroutine allocate_common(self)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.

Call graph

compute_auxiliary_fields

Compute auxiliary fields.

fortran
subroutine compute_auxiliary_fields(self)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.

Call graph

initialize

Initialize the equation common data.

fortran
subroutine initialize(self, filename, memory_avail, nv, verbose, L0)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inouttargetThe equation.
filenamecharacter(len=*)inInput file name.
memory_availreal(kind=R8P)invalueMemory available for single MPI process.
nvinteger(kind=I4P)inoptionalNumber of field variables.
verboselogicalinoptionalTrigger verbose output.
L0real(kind=R8P)inoptionalAdimensionalization parameter.

Call graph

load_restart_files

Save restart files.

fortran
subroutine load_restart_files(self, t, time)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.
tinteger(kind=I4P)outTime iteration.
timereal(kind=R8P)outTime.

Call graph

save_divergence_history

Save divergence history.

fortran
subroutine save_divergence_history(self, is_to_open, is_to_close)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.
is_to_openlogicalinoptionalFlag to open file before first saving.
is_to_closelogicalinoptionalFlag to close file after last saving.

Call graph

save_energy_error

Save energy error history.

fortran
subroutine save_energy_error(self, is_to_open, is_to_close)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.
is_to_openlogicalinoptionalFlag to open file before first saving.
is_to_closelogicalinoptionalFlag to close file after last saving.

Call graph

save_energy_history

Save energy history.

fortran
subroutine save_energy_history(self, is_to_open, is_to_close)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.
is_to_openlogicalinoptionalFlag to open file before first saving.
is_to_closelogicalinoptionalFlag to close file after last saving.

Call graph

save_restart_files

Save restart files.

fortran
subroutine save_restart_files(self)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.

Call graph

save_xh5f

Save simulation data in HDF5 format.

fortran
subroutine save_xh5f(self, output_basename, with_ghost)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.
output_basenamecharacter(len=*)inoptionalOutput basename.
with_ghostlogicalinoptionalFlag to save ghost cells.

Call graph

initialize_coils

Initialize coils.

fortran
subroutine initialize_coils(self)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutThe equation.

Call graph

set_helicon_coil

Set helicon coil.

fortran
subroutine set_helicon_coil(self, n)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.

Call graph

compute_Helicon_coil_amplitude

Compute amplitude correction for helicon coil

fortran
subroutine compute_Helicon_coil_amplitude(I_target, rb, ra, r1, r2, sigma, n, amplitude)

Arguments

NameTypeIntentAttributesDescription
I_targetreal(kind=R8P)in
rbreal(kind=R8P)in
rareal(kind=R8P)in
r1real(kind=R8P)in
r2real(kind=R8P)in
sigmareal(kind=R8P)in
ninteger(kind=I4P)in
amplitudereal(kind=R8P)out

Call graph

cartesian_to_cylindrical

Convert cartesian coordinates to cylindrical coordinates with respect to the coil center and normal direction.

Attributes: pure

fortran
subroutine cartesian_to_cylindrical(x, y, z, x_c, y_c, z_c, normal, r, theta, axial)

Arguments

NameTypeIntentAttributesDescription
xreal(kind=R8P)inCartesian coordinates.
yreal(kind=R8P)inCartesian coordinates.
zreal(kind=R8P)inCartesian coordinates.
x_creal(kind=R8P)inCoil center coordinates.
y_creal(kind=R8P)inCoil center coordinates.
z_creal(kind=R8P)inCoil center coordinates.
normalcharacter(len=2)inNormal direction.
rreal(kind=R8P)outCylindrical coordinates.
thetareal(kind=R8P)outCylindrical coordinates.
axialreal(kind=R8P)outCylindrical coordinates.

Call graph

cylindrical_to_cartesian

Convert cylindrical coordinates to cartesian coordinates with respect to the coil center and normal direction.

Attributes: pure

fortran
subroutine cylindrical_to_cartesian(r, theta, axial, x_c, y_c, z_c, normal, x, y, z)

Arguments

NameTypeIntentAttributesDescription
rreal(kind=R8P)in
thetareal(kind=R8P)in
axialreal(kind=R8P)in
x_creal(kind=R8P)in
y_creal(kind=R8P)in
z_creal(kind=R8P)in
normalcharacter(len=2)in
xreal(kind=R8P)out
yreal(kind=R8P)out
zreal(kind=R8P)out

cylindrical_vector_to_cartesian

Convert cylindrical vector A components to cartesian vector components.

Attributes: pure

fortran
subroutine cylindrical_vector_to_cartesian(A_r, A_theta, A_axial, theta, normal, A_x, A_y, A_z)

Arguments

NameTypeIntentAttributesDescription
A_rreal(kind=R8P)in
A_thetareal(kind=R8P)in
A_axialreal(kind=R8P)in
thetareal(kind=R8P)in
normalcharacter(len=2)in
A_xreal(kind=R8P)out
A_yreal(kind=R8P)out
A_zreal(kind=R8P)out

Call graph

set_rectangular_coil_x

Set rectangular coil with normal direction parallel to x.

fortran
subroutine set_rectangular_coil_x(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inCoil normal direction, +1=+x, -1=-x.

Call graph

set_rectangular_coil_y

Set rectangular coil with normal direction parallel to y.

fortran
subroutine set_rectangular_coil_y(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inCoil normal direction, +1=+y, -1=-y.

Call graph

set_rectangular_coil_z

Set rectangular coil with normal direction parallel to z.

fortran
subroutine set_rectangular_coil_z(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inCoil normal direction, +1=+z, -1=-z.

Call graph

set_circular_coil_x

Set circular coil with normal direction parallel to x.

fortran
subroutine set_circular_coil_x(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inCoil normal direction, +1=+x, -1=-x.

Call graph

set_circular_coil_y

Set circular coil with normal direction parallel to y.

fortran
subroutine set_circular_coil_y(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inCoil normal direction, +1=+y, -1=-y.

Call graph

set_circular_coil_z

Set circular coil with normal direction parallel to z.

fortran
subroutine set_circular_coil_z(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inCoil normal direction, +1=+z, -1=-z.

Call graph

set_solenoid_x

Set solenoid with axis direction parallel to x.

fortran
subroutine set_solenoid_x(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inSolenoid axis direction, +1=+x, -1=-x.

Call graph

set_solenoid_y

Set solenoid with axis direction parallel to y.

fortran
subroutine set_solenoid_y(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inSolenoid axis direction, +1=+y, -1=-y.

Call graph

set_solenoid_z

Set solenoid with axis direction parallel to z.

fortran
subroutine set_solenoid_z(self, n, verse)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inoutCpu object.
ninteger(kind=I4P)inCoil number.
versereal(kind=R8P)inSolenoid axis direction, +1=+z, -1=-z.

Call graph

coupling_descriptor_forest

Report this realm's coupling descriptor for β admissibility.

Returns (numerics%scheme_time, rk%scheme, physics%nv). β admits a seam between two realms iff all three values agree. Spatial operator (numerics%scheme_space) is intentionally NOT in the descriptor — β exists precisely to support different spatial operators per realm.

Provided on prism_common_object so both prism_cpu_object and prism_fnl_object inherit a single implementation.

fortran
subroutine coupling_descriptor_forest(self, scheme_time, rk_scheme, nv)

Arguments

NameTypeIntentAttributesDescription
selfclass(prism_common_object)inThe realm.
scheme_timecharacter(len=:)outallocatableTime-integration family tag ("runge_kutta", "leapfrog", ...).
rk_schemecharacter(len=:)outallocatableWithin-family scheme tag ("runge-kutta-ssp-54", ...).
nvinteger(kind=I4P)outNumber of conserved variables on this realm.

Call graph

Functions

helicon_coil_A_at_point

Compute helicon/saddle coil auxiliary vector potential A at one point.

Returns: real(kind=R8P)

fortran
function helicon_coil_A_at_point(cell_coord, x_c, y_c, z_c, normal, Radius, N_points, s_map, csi_map, area_signed, eta, theta_c, sigma) result(Ap)

Arguments

NameTypeIntentAttributesDescription
cell_coordreal(kind=R8P)inCoordinate cartesiane del punto di valutazione.
x_creal(kind=R8P)inCoordinata x centro cilindro.
y_creal(kind=R8P)inCoordinata x centro cilindro.
z_creal(kind=R8P)inCoordinata z centro cilindro.
normalcharacter(len=2)inAsse cilindro.
Radiusreal(kind=R8P)inRaggio cilindro.
N_pointsinteger(kind=I4P)inNumero di punti della spline.
s_mapreal(kind=R8P)ins coordinates of the spline points for the helicon coil case.
csi_mapreal(kind=R8P)incsi coordinates of the spline points for the helicon coil case.
area_signedreal(kind=R8P)in
etareal(kind=R8P)in
theta_creal(kind=R8P)in
sigmareal(kind=R8P)inSmearing.

Call graph

wrap_to_pi

Wrap angle to [-pi, pi).

Attributes: pure, elemental

Returns: real(kind=R8P)

fortran
function wrap_to_pi(alpha) result(alpha_w)

Arguments

NameTypeIntentAttributesDescription
alphareal(kind=R8P)in

Call graph

compute_signed_area

Compute signed area of a a 2D polygon

Returns: real(kind=R8P)

fortran
function compute_signed_area(N_points, x_points, y_points) result(signed_area)

Arguments

NameTypeIntentAttributesDescription
N_pointsinteger(kind=I4P)in
x_pointsreal(kind=R8P)in
y_pointsreal(kind=R8P)in

Call graph

compute_windings_number

Returns: real(kind=R8P)

fortran
function compute_windings_number(N_points, x_points, y_points, coord) result(w_n)

Arguments

NameTypeIntentAttributesDescription
N_pointsinteger(kind=I4P)in
x_pointsreal(kind=R8P)in
y_pointsreal(kind=R8P)in
coordreal(kind=R8P)in

Call graph

rectangular_coil_A_at_point

Returns: real(kind=R8P)

fortran
function rectangular_coil_A_at_point(xp, xc, idir_n, idir_a, idir_b, a1, a2, b1, b2, sigma) result(Ap)

Arguments

NameTypeIntentAttributesDescription
xpreal(kind=R8P)inPunto di valutazione.
xcreal(kind=R8P)inCentro spira.
idir_ninteger(kind=I4P)inComponente normale: 1=x, 2=y, 3=z.
idir_ainteger(kind=I4P)inPrima coordinata locale nel piano spira.
idir_binteger(kind=I4P)inSeconda coordinata locale nel piano spira.
a1real(kind=R8P)inEstremi lato locale a.
a2real(kind=R8P)inEstremi lato locale a.
b1real(kind=R8P)inEstremi lato locale b.
b2real(kind=R8P)inEstremi lato locale b.
sigmareal(kind=R8P)inSmearing unico della spira.

Call graph

erf_function

Returns: real(kind=R8P)

fortran
function erf_function(s, mu, sigma) result(res)

Arguments

NameTypeIntentAttributesDescription
sreal(kind=R8P)in
mureal(kind=R8P)in
sigmareal(kind=R8P)in

Call graph

tangential_window

Returns: real(kind=R8P)

fortran
function tangential_window(s, smin, smax, sigma) result(res)

Arguments

NameTypeIntentAttributesDescription
sreal(kind=R8P)in
sminreal(kind=R8P)in
smaxreal(kind=R8P)in
sigmareal(kind=R8P)in

Call graph

erf_primitive_function

Primitive of erf((s-mu)/(sqrt(2)*sigma)).

Returns: real(kind=R8P)

fortran
function erf_primitive_function(s, mu, sigma) result(res)

Arguments

NameTypeIntentAttributesDescription
sreal(kind=R8P)in
mureal(kind=R8P)in
sigmareal(kind=R8P)in

Call graph

build_A_ghost_stencil

Build ghost cell A array to compute J_vec value through the curl, useful in multiblock application

Returns: real(kind=R8P)

fortran
function build_A_ghost_stencil(hs, dxyz, gc_coord, x_c, y_c, z_c, sigma, i_dir_n, i_dir_a, i_dir_b, a1, a2, b1, b2, N_Points, s_map, csi_map, area_signed, eta, theta_c, r_coil, normal, coil_type) result(A_gc)

Arguments

NameTypeIntentAttributesDescription
hsinteger(kind=I4P)inHalf-stencil size for the curl computation.
dxyzreal(kind=R8P)inSpace steps.
gc_coordreal(kind=R8P)inGhost cell coordinates.
x_creal(kind=R8P)inCenter coordinates.
y_creal(kind=R8P)inCenter coordinates.
z_creal(kind=R8P)inCenter coordinates.
sigmareal(kind=R8P)inWidth of the Gaussian profile.
i_dir_ninteger(kind=I4P)inDirection indices for the normal and tangential directions.
i_dir_ainteger(kind=I4P)inDirection indices for the normal and tangential directions.
i_dir_binteger(kind=I4P)inDirection indices for the normal and tangential directions.
a1real(kind=R8P)inoptionalLocal coordinates of the rectangular coil along the first
a2real(kind=R8P)inoptionalLocal coordinates of the rectangular coil along the first
b1real(kind=R8P)inoptionalLocal coordinates of the rectangular coil along the second
b2real(kind=R8P)inoptionalLocal coordinates of the rectangular coil along the second
N_Pointsinteger(kind=I4P)in
s_mapreal(kind=R8P)inoptionals coordinates of the spline points for the helicon coil case.
csi_mapreal(kind=R8P)inoptionalcsi coordinates of the spline points for the helicon coil case.
area_signedreal(kind=R8P)inoptional
etareal(kind=R8P)inoptional
theta_creal(kind=R8P)inoptional
r_coilreal(kind=R8P)inoptionalRadius of the helicon coil.
normalcharacter(len=2)inoptionalNormal vector for the helicon coil case.
coil_typereal(kind=R8P)inCoil type.

Call graph

crossproduct

Returns: real(kind=R8P)

fortran
function crossproduct(a, b) result(cross)

Arguments

NameTypeIntentAttributesDescription
areal(kind=R8P)inLeft hand side.
breal(kind=R8P)inLeft hand side.

Call graph