XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Public Member Functions | Public Attributes | Private Member Functions | List of all members
sheath_module Module Reference

Module for variables and routines related to the logical sheath boundary condition. More...

Public Member Functions

subroutine init_sheath_mod (grid, sheath_mode)
 Initializes sheath module, i.e. memory allocation initialization of variables, and for sheath_mode==2 identification of vertex groups among the wall vertices. More...
 
subroutine finalize_sheath_mod (sheath_mode)
 Finalize sheath module: deallocation of sheath variables. More...
 
subroutine sheath_ptl_mem_alloc (nptl)
 Allocates memory for the particle based data required for sheath_mode==2. More...
 
subroutine sheath_ptl_mem_dealloc
 Deallocates particle memory used in sheath calculation with sheath_mode==2. More...
 
subroutine sheath_calculation (grid, psn, sp, iptl, type, itrout, pout, ith, sheath_mode)
 High level sheath routine to simplify calls Depending on the desired operating mode, the appropriate subroutine is called. More...
 
subroutine sheath_adjust (grid, sp, sheath_mode)
 High level routine for the sheath adjustment. More...
 
subroutine sheath_calculation_mode2 (grid, psn, sp, iptl, type, itrout, pout, ith)
 Updated sheath routine At this stage, we only gather a PDF in v_|| of particles hitting the wall so that we can determing the required sheath potential more accurately. Electrons will be picked for removal afterwards in a separate routine –> sheath_adjust_mode2. More...
 
subroutine sheath_adjust_mode2 (grid, sp)
 This new version of the logical sheath adustment does three things: 1) compute use energy distribution of incoming electrons to find the new sheath potential and 2) loop over electrons and decide which ones penetrate the sheath. More...
 
subroutine compute_elost_en (grid, sp, iep0, dep_inv)
 Computes the energy distribution of incoming electrons and the lost electron charge for the current value of the sheath potential sheath_pot. More...
 
subroutine find_sheath_pot (grid, dep, loss_shift, iep0, loss, q_minmax, solved)
 Finds sheath potential at which the electron loss is (approx.) equal to the ion loss plus a given shift (loss_shift) More...
 
subroutine sheath_diag_output (grid, i)
 Outputs some diagnostic data that helps to check the correctness of the algorithm that determines the sheath potential when sheath_mode==2. More...
 
subroutine sheath_calculation_mode1 (grid, psn, sp, iptl, stype, itrout, pout, ith)
 Old sheath routine (sheath_mode==1) Decides whether an electron is reflected in the sheath with a given value of the sheath potential and runs diagnostics. More...
 
subroutine sheath_adjust_mode1 (grid)
 Old sheath calculation (sheath_mode==1) Based on the electron and ion losses at the wall nodes a new value of the sheath potential is determined by adjusting the previous value by a certain amount –> feedback control loop. More...
 

Public Attributes

integer sheath_nep =100
 number of energy grid bins for sheath calculation More...
 
real(kind=8), dimension(:,:,:),
allocatable, target 
sheath_en_max
 Parallel velo cutoff (normalized to sqrt(T_e)) for sheath potential calculation. More...
 
real(kind=8), dimension(:,:),
allocatable 
sheath_dep
 Parallel energy grid bin size dep = en_max/nep. More...
 
real(kind=8), dimension(:,:),
allocatable 
sheath_dep_inv
 Inverse of grid bin size. More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
sheath_elost_en
 energy distribution of electron charge loss More...
 
real(kind=8), dimension(:),
allocatable 
sheath_te_ev
 Initial wall temperature for normalization. More...
 
integer sheath_nphi
 Number of poloidal planes for sheath potential. More...
 
real(8), dimension(:,:),
allocatable, target 
sheath_pot
 Sheath potential on the wall nodes. More...
 
real(8), dimension(:,:,:),
allocatable 
sheath_lost
 Total charge lost to the wall. More...
 
real(8), dimension(:,:,:),
allocatable 
sheath_ilost
 Ion charge lost to the wall. More...
 
integer sheath_ngroups
 Number of vertex groups for sheath calculation. More...
 
integer sheath_group_maxlen
 Max. number of members of a sheath vertex group. More...
 
integer, dimension(:), allocatable sheath_igroup
 Index of group to which a wall vertex belongs. More...
 
integer, dimension(:), allocatable sheath_group_len
 Number of vertices that belong to a group. More...
 
integer, dimension(:,:),
allocatable 
sheath_group_idx
 Indices of wall vertices that belong to a group. More...
 
integer, dimension(:),
allocatable, target 
sheath_ptl_widx
 Index of wall vertex where an electron crossed the wall. More...
 
integer, dimension(:),
allocatable, target 
sheath_ptl_iphi
 Index of toroidal sector where an electron crossed the wall. More...
 
real(8), dimension(:),
allocatable, target 
sheath_ptl_en_para
 Parallel kinetic energy of electrons that crossed the wall. More...
 
real(8), dimension(:),
allocatable, target 
sheath_ptl_en_perp
 Perpendicular kinetic energy of electrons that crossed the wall. More...
 
real(8), dimension(:,:,:),
allocatable, target 
sheath_ptl_ph
 Buffer for particle coordinates (R,Z,phi) before and after the wall crossing. More...
 

Private Member Functions

subroutine, private sheath_adjust_init (nwall, nphi, nthreads)
 Initialize memory for sheath calculation in sheath_adjust_mode2. More...
 
subroutine, private sheath_adjust_finalize
 Free memory for sheath allocation. More...
 

Detailed Description

Module for variables and routines related to the logical sheath boundary condition.

Member Function/Subroutine Documentation

subroutine sheath_module::compute_elost_en ( type(grid_type), intent(in)  grid,
type(species_type), intent(in)  sp,
integer, dimension(grid%nwall,sheath_nphi), intent(in)  iep0,
real (kind=8), dimension(grid%nwall,sheath_nphi), intent(in)  dep_inv 
)

Computes the energy distribution of incoming electrons and the lost electron charge for the current value of the sheath potential sheath_pot.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in,out]spParticle data, type(species_type)
[in]iep0Bin on coarse energy grid in which to compute refined electron energy distribution, integer
[in]dep_invGrid size of energy grid, real(8)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::finalize_sheath_mod ( integer, intent(in)  sheath_mode)

Finalize sheath module: deallocation of sheath variables.

Parameters
[in]sheath_modeOperating mode of the sheath calculation

Here is the call graph for this function:

subroutine sheath_module::find_sheath_pot ( type(grid_type), intent(in)  grid,
real (kind=8), dimension(grid%nwall,sheath_nphi), intent(in)  dep,
real (kind=8), dimension(grid%nwall,sheath_nphi), intent(in)  loss_shift,
integer, dimension(grid%nwall,sheath_nphi), intent(inout)  iep0,
real (kind=8), dimension(grid%nwall,sheath_nphi), intent(out)  loss,
real (kind=8), dimension(grid%nwall,sheath_nphi,2), intent(out)  q_minmax,
logical, dimension(grid%nwall,sheath_nphi), intent(out)  solved 
)

Finds sheath potential at which the electron loss is (approx.) equal to the ion loss plus a given shift (loss_shift)

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in]depResolution of energy grid, real(8)
[in]loss_shiftShift in lost charge used for global loss balance, real(8)
[in,out]iep0Used for refinement of sheath potential on fine energy grid, integer
[out]lossApprox. charge loss determined on energy grid, real(8)
[out]q_minmaxMin. and max. of possible electron loss, real(8)
[out]solvedWhether local charge balance is successful on a node, logical

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::init_sheath_mod ( type(grid_type), intent(in)  grid,
integer, intent(in)  sheath_mode 
)

Initializes sheath module, i.e. memory allocation initialization of variables, and for sheath_mode==2 identification of vertex groups among the wall vertices.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in]sheath_modeOperating mode of the sheath calculation: 1) 1) feedback control adjustment 2) Instantaneous sheath adjustment based on energy distribution of incoming electrons

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::sheath_adjust ( type(grid_type), intent(in)  grid,
type(species_type), intent(inout)  sp,
integer, intent(in)  sheath_mode 
)

High level routine for the sheath adjustment.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in,out]spParticle data, type(species_type)
[in]sheath_modeOperating mode

Here is the call graph for this function:

subroutine, private sheath_module::sheath_adjust_finalize ( )
private

Free memory for sheath allocation.

Here is the caller graph for this function:

subroutine, private sheath_module::sheath_adjust_init ( integer, intent(in)  nwall,
integer, intent(in)  nphi,
integer, intent(in)  nthreads 
)
private

Initialize memory for sheath calculation in sheath_adjust_mode2.

Parameters
[in]nwallNumber of wall nodes, integer
[in]nphiNumber of toroidal segments, integer
[in]nthreadsNumber of OpenMP threads, integer

Here is the caller graph for this function:

subroutine sheath_module::sheath_adjust_mode1 ( type(grid_type)  grid)

Old sheath calculation (sheath_mode==1) Based on the electron and ion losses at the wall nodes a new value of the sheath potential is determined by adjusting the previous value by a certain amount –> feedback control loop.

Parameters
[in]gridUnstructured mesh data, type(grid_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::sheath_adjust_mode2 ( type(grid_type)  grid,
type(species_type), intent(inout)  sp 
)

This new version of the logical sheath adustment does three things: 1) compute use energy distribution of incoming electrons to find the new sheath potential and 2) loop over electrons and decide which ones penetrate the sheath.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in,out]spParticle data, type(species_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::sheath_calculation ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn,
type(species_type), intent(inout)  sp,
integer, intent(in)  iptl,
integer, intent(in)  type,
integer, intent(inout)  itrout,
real (kind=8), dimension(3), intent(inout)  pout,
integer, intent(in)  ith,
integer, intent(in)  sheath_mode 
)

High level sheath routine to simplify calls Depending on the desired operating mode, the appropriate subroutine is called.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in,out]spParticle data, type(species_type)
[in]iptlIndex of particle to process, integer
[in]typeSpecies type of particle, integer
[in,out]itroutOutput triangle index, integer
[in,out]poutTriangle interpolation weights, real(8)
[in]ithOMP thread index, integer
[in]sheath_modeOperating mode

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::sheath_calculation_mode1 ( type(grid_type)  grid,
type(psn_type), intent(inout)  psn,
type(species_type)  sp,
integer, intent(in)  iptl,
integer, intent(in)  stype,
integer, intent(out)  itrout,
real (kind=8), dimension(3), intent(out)  pout,
integer, intent(in)  ith 
)

Old sheath routine (sheath_mode==1) Decides whether an electron is reflected in the sheath with a given value of the sheath potential and runs diagnostics.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in,out]spParticle data, type(species_type)
[in]iptlIndex of particle to process, integer
[in]typeSpecies type of particle, integer
[in,out]itroutOutput triangle index, integer
[in,out]poutTriangle interpolation weights, real(8)
[in]ithOMP thread index, integer

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::sheath_calculation_mode2 ( type(grid_type)  grid,
type(psn_type)  psn,
type(species_type)  sp,
integer, intent(in)  iptl,
integer, intent(in)  type,
integer  itrout,
real (kind=8), dimension(3)  pout,
integer, intent(in)  ith 
)

Updated sheath routine At this stage, we only gather a PDF in v_|| of particles hitting the wall so that we can determing the required sheath potential more accurately. Electrons will be picked for removal afterwards in a separate routine –> sheath_adjust_mode2.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in,out]spParticle data, type(species_type)
[in]iptlIndex of particle to process, integer
[in]typeSpecies type of particle, integer
[in,out]itroutOutput triangle index, integer
[in,out]poutTriangle interpolation weights, real(8)
[in]ithOMP thread index, integer

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine sheath_module::sheath_diag_output ( type(grid_type), intent(in)  grid,
integer, intent(in)  i 
)

Outputs some diagnostic data that helps to check the correctness of the algorithm that determines the sheath potential when sheath_mode==2.

Parameters
[in]gridUnstructured mesh data, type(grid_type)
[in]iFlag indicating whether data should be written, integer
subroutine sheath_module::sheath_ptl_mem_alloc ( integer, intent(in)  nptl)

Allocates memory for the particle based data required for sheath_mode==2.

Parameters
[in]nptlMax. number of particles/rank

Here is the caller graph for this function:

subroutine sheath_module::sheath_ptl_mem_dealloc ( )

Deallocates particle memory used in sheath calculation with sheath_mode==2.

Here is the caller graph for this function:

Member Data Documentation

real (kind=8), dimension(:,:), allocatable sheath_module::sheath_dep

Parallel energy grid bin size dep = en_max/nep.

real (kind=8), dimension(:,:), allocatable sheath_module::sheath_dep_inv

Inverse of grid bin size.

real (kind=8), dimension(:,:,:,:), allocatable sheath_module::sheath_elost_en

energy distribution of electron charge loss

real (kind=8), dimension(:,:,:), allocatable, target sheath_module::sheath_en_max

Parallel velo cutoff (normalized to sqrt(T_e)) for sheath potential calculation.

integer, dimension(:,:), allocatable sheath_module::sheath_group_idx

Indices of wall vertices that belong to a group.

integer, dimension(:), allocatable sheath_module::sheath_group_len

Number of vertices that belong to a group.

integer sheath_module::sheath_group_maxlen

Max. number of members of a sheath vertex group.

integer, dimension(:), allocatable sheath_module::sheath_igroup

Index of group to which a wall vertex belongs.

real (8), dimension(:,:,:), allocatable sheath_module::sheath_ilost

Ion charge lost to the wall.

real (8), dimension(:,:,:), allocatable sheath_module::sheath_lost

Total charge lost to the wall.

integer sheath_module::sheath_nep =100

number of energy grid bins for sheath calculation

integer sheath_module::sheath_ngroups

Number of vertex groups for sheath calculation.

integer sheath_module::sheath_nphi

Number of poloidal planes for sheath potential.

real (8), dimension(:,:), allocatable, target sheath_module::sheath_pot

Sheath potential on the wall nodes.

real (8), dimension(:), allocatable, target sheath_module::sheath_ptl_en_para

Parallel kinetic energy of electrons that crossed the wall.

real (8), dimension(:), allocatable, target sheath_module::sheath_ptl_en_perp

Perpendicular kinetic energy of electrons that crossed the wall.

integer, dimension(:), allocatable, target sheath_module::sheath_ptl_iphi

Index of toroidal sector where an electron crossed the wall.

real (8), dimension(:,:,:), allocatable, target sheath_module::sheath_ptl_ph

Buffer for particle coordinates (R,Z,phi) before and after the wall crossing.

integer, dimension(:), allocatable, target sheath_module::sheath_ptl_widx

Index of wall vertex where an electron crossed the wall.

real (kind=8), dimension(:), allocatable sheath_module::sheath_te_ev

Initial wall temperature for normalization.


The documentation for this module was generated from the following file: