Report Typos and Errors    
Semi-Lagrangian Library
Modular library for kinetic and gyrokinetic simulations of plasmas in fusion energy devices.
Derived types and interfaces | Functions/Subroutines
sll_m_mapping_3d Module Reference

Description

Module interfaces for coordinate transformation.

xi1, xi2, xi3 are the logical coordinates which are transformed to the physical coordinates x1, x2, x3 The Jacobian matrix is the derivate matrix of the transformation function the determinant of the Jacobian is needed for the transformation theorem the logical paramters are in the interval [0,1]

Author
Benedikt Perse

Derived types and interfaces

interface  sll_i_eval_function
 abstract interface for mapping functions More...
 
type  matrix_element
 
type  sll_t_mapping_3d
 type collecting functions for analytical coordinate mapping More...
 

Functions/Subroutines

real(kind=f64) function get_x1 (self, xi)
 Compute x1. More...
 
real(kind=f64) function get_x2 (self, xi)
 Compute x2. More...
 
real(kind=f64) function get_x3 (self, xi)
 Compute x3. More...
 
real(kind=f64) function, dimension(3) get_x (self, xi)
 Compute all three components of X. More...
 
real(kind=f64) function, dimension(3) get_xi (self, x)
 Compute the logical coordinate vector Xi. More...
 
real(kind=f64) function jacobian (self, xi)
 Compute the determinant of the Jacobi matrix. More...
 
real(kind=f64) function, dimension(3, 3) jacobian_matrix (self, xi)
 Compute the entries of the Jacobi matrix. More...
 
real(kind=f64) function, dimension(3, 3) jacobian_matrix_transposed (self, xi)
 Compute the entries of the transposed Jacobi matrix. More...
 
real(kind=f64) function, dimension(3, 3) jacobian_matrix_inverse (self, xi)
 Compute the entries of the inverse Jacobi matrix. More...
 
real(kind=f64) function, dimension(3, 3) jacobian_matrix_inverse_transposed (self, xi)
 Compute the entries of the transposed inverse Jacobi matrix. More...
 
real(kind=f64) function, dimension(3, 3) metric (self, xi)
 Compute the entries of the metric. More...
 
real(kind=f64) function metric_single (self, xi, component1, component2)
 Compute a single entriy of the metric. More...
 
real(kind=f64) function, dimension(3, 3) metric_inverse (self, xi)
 Compute the entries of the inverse metric. More...
 
real(kind=f64) function metric_inverse_single (self, xi, component1, component2)
 Compute a single entriy of the inverse metric. More...
 
real(kind=f64) function metric_single_jacobian (self, xi, component1, component2)
 Compute a single entriy of the metric divided by the jacobian. More...
 
real(kind=f64) function metric_inverse_single_jacobian (self, xi, component1, component2)
 Compute a single entriy of the inverse metric multiplied by the jacobian. More...
 
subroutine init (self, params, x1_func, x2_func, x3_func, jac11, jac12, jac13, jac21, jac22, jac23, jac31, jac32, jac33, jacob, xi1_func, xi2_func, xi3_func, flag2d, flag3d, n_cells, deg, Lx, volume)
 Initialize the coordinate transformation. More...
 
subroutine init_from_file (self, filename)
 Initialize the coordinate transformation from nml file. More...
 
subroutine free (self)
 Finalize the coordinate transformation. More...
 
subroutine calculate_interpolation_matrix_1d (n_cells, deg, xk, spline, matrix)
 Helper function for spline mapping. More...
 
subroutine calculate_interpolation_matrix_1d_periodic (n_cells, deg, spline, matrix)
 Helper function for spline mapping. More...
 
subroutine convert_x_to_xbox (self, position, xi, box)
 Helper function to compute xbox. More...
 

Function/Subroutine Documentation

◆ calculate_interpolation_matrix_1d()

subroutine sll_m_mapping_3d::calculate_interpolation_matrix_1d ( integer(kind=i32), intent(in)  n_cells,
integer(kind=i32), intent(in)  deg,
real(kind=f64), dimension(:), intent(in)  xk,
type( sll_t_spline_pp_1d), intent(in)  spline,
type(sll_t_matrix_csr), intent(out)  matrix 
)
private

Helper function for spline mapping.

Definition at line 938 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calculate_interpolation_matrix_1d_periodic()

subroutine sll_m_mapping_3d::calculate_interpolation_matrix_1d_periodic ( integer(kind=i32), intent(in)  n_cells,
integer(kind=i32), intent(in)  deg,
type( sll_t_spline_pp_1d), intent(in)  spline,
type(sll_t_matrix_csr), intent(out)  matrix 
)
private

Helper function for spline mapping.

Definition at line 1013 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ convert_x_to_xbox()

subroutine sll_m_mapping_3d::convert_x_to_xbox ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  position,
real(kind=f64), dimension(3), intent(out)  xi,
integer(kind=i32), dimension(3), intent(out)  box 
)
private

Helper function to compute xbox.

Parameters
[in,out]selfkernel smoother object
[in]positionPosition of the particle
[out]xiPosition of the particle
[out]boxPosition of the particle

Definition at line 1088 of file sll_m_mapping_3d.F90.

Here is the caller graph for this function:

◆ free()

subroutine sll_m_mapping_3d::free ( class(sll_t_mapping_3d), intent(inout)  self)
private

Finalize the coordinate transformation.

Definition at line 924 of file sll_m_mapping_3d.F90.

◆ get_x()

real(kind=f64) function, dimension(3) sll_m_mapping_3d::get_x ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute all three components of X.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
physical coordinates

Definition at line 200 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:

◆ get_x1()

real(kind=f64) function sll_m_mapping_3d::get_x1 ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute x1.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
physical coordinate

Definition at line 143 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:

◆ get_x2()

real(kind=f64) function sll_m_mapping_3d::get_x2 ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute x2.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
physical coordinate

Definition at line 162 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:

◆ get_x3()

real(kind=f64) function sll_m_mapping_3d::get_x3 ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute x3.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
physical coordinate

Definition at line 181 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:

◆ get_xi()

real(kind=f64) function, dimension(3) sll_m_mapping_3d::get_xi ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  x 
)
private

Compute the logical coordinate vector Xi.

Parameters
[in,out]selfcoordinate transformation
[in]xphysical coordinates
Returns
logical coordinates

Definition at line 223 of file sll_m_mapping_3d.F90.

◆ init()

subroutine sll_m_mapping_3d::init ( class(sll_t_mapping_3d), intent(out)  self,
real(kind=f64), dimension(:), intent(in)  params,
procedure(sll_i_eval_function x1_func,
procedure(sll_i_eval_function x2_func,
procedure(sll_i_eval_function x3_func,
procedure(sll_i_eval_function), optional  jac11,
procedure(sll_i_eval_function), optional  jac12,
procedure(sll_i_eval_function), optional  jac13,
procedure(sll_i_eval_function), optional  jac21,
procedure(sll_i_eval_function), optional  jac22,
procedure(sll_i_eval_function), optional  jac23,
procedure(sll_i_eval_function), optional  jac31,
procedure(sll_i_eval_function), optional  jac32,
procedure(sll_i_eval_function), optional  jac33,
procedure(sll_i_eval_function), optional  jacob,
procedure(sll_i_eval_function), optional  xi1_func,
procedure(sll_i_eval_function), optional  xi2_func,
procedure(sll_i_eval_function), optional  xi3_func,
logical, optional  flag2d,
logical, optional  flag3d,
integer(kind=i32), dimension(3), optional  n_cells,
integer(kind=i32), dimension(3), optional  deg,
real(kind=f64), dimension(3), optional  Lx,
real(kind=f64), optional  volume 
)
private

Initialize the coordinate transformation.

Parameters
[out]selfcoordinate transformation
[in]paramstransformation parameters
x1_functransformation function
x2_functransformation function
x3_functransformation function
jac11entry of jacobian matrix
jac12entry of jacobian matrix
jac13entry of jacobian matrix
jac21entry of jacobian matrix
jac22entry of jacobian matrix
jac23entry of jacobian matrix
jac31entry of jacobian matrix
jac32entry of jacobian matrix
jac33entry of jacobian matrix
jacobdeterminant of jacobian matrix
xi1_functransformation function
xi2_functransformation function
xi3_functransformation function
flag2dlogical flag for mappings mixing the first two logical coordiantes
flag3dlogical flag for mappings mixing all three logical coordinates
n_cellsnumber of cells (and grid points)
degspline deg
lxlength of the physical domain
volumevolume of the physical domain

Definition at line 526 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:

◆ init_from_file()

subroutine sll_m_mapping_3d::init_from_file ( class(sll_t_mapping_3d), intent(out)  self,
character(len=*), intent(in)  filename 
)
private

Initialize the coordinate transformation from nml file.

Parameters
[out]selfcoordinate transformation

mapping cases defined in sll_m_3d_coordinate_transformations

Definition at line 738 of file sll_m_mapping_3d.F90.

◆ jacobian()

real(kind=f64) function sll_m_mapping_3d::jacobian ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the determinant of the Jacobi matrix.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates

Definition at line 235 of file sll_m_mapping_3d.F90.

◆ jacobian_matrix()

real(kind=f64) function, dimension(3,3) sll_m_mapping_3d::jacobian_matrix ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the entries of the Jacobi matrix.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
jacobian matrix

Definition at line 262 of file sll_m_mapping_3d.F90.

Here is the call graph for this function:

◆ jacobian_matrix_inverse()

real(kind=f64) function, dimension(3,3) sll_m_mapping_3d::jacobian_matrix_inverse ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the entries of the inverse Jacobi matrix.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
inverse jacobian matrix

Definition at line 325 of file sll_m_mapping_3d.F90.

◆ jacobian_matrix_inverse_transposed()

real(kind=f64) function, dimension(3,3) sll_m_mapping_3d::jacobian_matrix_inverse_transposed ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the entries of the transposed inverse Jacobi matrix.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
jacobian matrix

Definition at line 366 of file sll_m_mapping_3d.F90.

◆ jacobian_matrix_transposed()

real(kind=f64) function, dimension(3,3) sll_m_mapping_3d::jacobian_matrix_transposed ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the entries of the transposed Jacobi matrix.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
transposed jacobian matrix

Definition at line 303 of file sll_m_mapping_3d.F90.

◆ metric()

real(kind=f64) function, dimension(3,3) sll_m_mapping_3d::metric ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the entries of the metric.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
metric

Definition at line 389 of file sll_m_mapping_3d.F90.

◆ metric_inverse()

real(kind=f64) function, dimension(3,3) sll_m_mapping_3d::metric_inverse ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi 
)
private

Compute the entries of the inverse metric.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
Returns
inverse metric

Definition at line 427 of file sll_m_mapping_3d.F90.

◆ metric_inverse_single()

real(kind=f64) function sll_m_mapping_3d::metric_inverse_single ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi,
integer(kind=i32), intent(in)  component1,
integer(kind=i32), intent(in)  component2 
)
private

Compute a single entriy of the inverse metric.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
[in]component2components of the wished entry
Returns
single entry of the inverse metric

Definition at line 448 of file sll_m_mapping_3d.F90.

◆ metric_inverse_single_jacobian()

real(kind=f64) function sll_m_mapping_3d::metric_inverse_single_jacobian ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi,
integer(kind=i32), intent(in)  component1,
integer(kind=i32), intent(in)  component2 
)
private

Compute a single entriy of the inverse metric multiplied by the jacobian.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
[in]component2components of the wished entry
Returns
single entry of the inverse metric multiplied by the jacobian

Definition at line 493 of file sll_m_mapping_3d.F90.

◆ metric_single()

real(kind=f64) function sll_m_mapping_3d::metric_single ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi,
integer(kind=i32), intent(in)  component1,
integer(kind=i32), intent(in)  component2 
)
private

Compute a single entriy of the metric.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
[in]component2components of the wished entry
Returns
single entry of the metric

Definition at line 410 of file sll_m_mapping_3d.F90.

◆ metric_single_jacobian()

real(kind=f64) function sll_m_mapping_3d::metric_single_jacobian ( class(sll_t_mapping_3d), intent(inout)  self,
real(kind=f64), dimension(3), intent(in)  xi,
integer(kind=i32), intent(in)  component1,
integer(kind=i32), intent(in)  component2 
)
private

Compute a single entriy of the metric divided by the jacobian.

Parameters
[in,out]selfcoordinate transformation
[in]xilogical coordinates
[in]component2components of the wished entry
Returns
single entry of the metric divided by the jacobian

Definition at line 464 of file sll_m_mapping_3d.F90.

    Report Typos and Errors