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_particle_mesh_coupling_spline_cl_1d Module Reference

Description

Kernel smoother for 2d with splines of arbitrary degree placed on a uniform mesh.

Author
Benedikt Perse, IPP

Spline with index i starts at point i

Derived types and interfaces

type  sll_t_particle_mesh_coupling_spline_cl_1d
 Spline kernel smoother in1d. More...
 

Functions/Subroutines

subroutine add_charge_single_spline_cl_1d (self, position, marker_charge, rho_dofs)
 Destructor. More...
 
subroutine add_particle_mass_spline_cl_1d (self, position, marker_charge, particle_mass)
 Add the contribution of one particle to the approximate mass matrix. More...
 
subroutine add_particle_mass_spline_cl_1d_full (self, position, marker_charge, particle_mass)
 Add the contribution of one particle to the approximate mass matrix. More...
 
subroutine evaluate_field_single_spline_cl_1d (self, position, field_dofs, field_value)
 Evaluate field at at position position. More...
 
subroutine evaluate_multiple_spline_cl_1d (self, position, components, field_dofs, field_value)
 Evaluate several fields at position position. More...
 
subroutine add_current_spline_cl_1d (self, position_old, position_new, marker_charge, j_dofs)
 Add current with integration over x. More...
 
subroutine update_j (self, lower, upper, index, marker_charge, sign, j_dofs)
 Helper function for add_current. More...
 
subroutine add_current_evaluate_spline_cl_1d (self, position_old, position_new, marker_charge, vbar, field_dofs, j_dofs, field)
 Combines add_current and evaluate_int. More...
 
subroutine update_je (self, lower, upper, index, marker_charge, sign, j_dofs, field_dofs, field)
 Helper function for add_current_evaluate. More...
 
subroutine add_current_update_v_spline_cl_1d (self, position_old, position_new, marker_charge, qoverm, bfield_dofs, vi, j_dofs)
 Add current for one particle and update v (according to H_p1 part in Hamiltonian splitting) More...
 
subroutine update_jv (self, lower, upper, index, marker_charge, qoverm, sign, vi, j_dofs, bfield_dofs)
 Helper function for add_current_update_v. More...
 
subroutine convert_x_to_xbox (self, position, xi, box)
 
subroutine sll_s_uniform_bsplines_eval_basis_clamped (spline, n_cells, degree, xi, box, spline_val)
 Helper function to evaluate uniform clamped basis splines. More...
 
subroutine init_spline_cl_1d (self, domain, n_cells, no_particles, spline_degree, smoothing_type, boundary)
 Initializer. More...
 
subroutine free_spline_cl_1d (self)
 Finalizer. More...
 
subroutine, public sll_s_new_particle_mesh_coupling_spline_cl_1d_ptr (smoother, domain, n_cells, no_particles, spline_degree, smoothing_type, boundary)
 
subroutine, public sll_s_new_particle_mesh_coupling_spline_cl_1d (smoother, domain, n_cells, no_particles, spline_degree, smoothing_type, boundary)
 

Function/Subroutine Documentation

◆ add_charge_single_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::add_charge_single_spline_cl_1d ( class( sll_t_particle_mesh_coupling_spline_cl_1d ), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), dimension(self%n_dofs), intent(inout)  rho_dofs 
)
private

Destructor.

Add charge of one particle

Parameters
[in,out]selfkernel smoother object
[in]positionPosition of the particle
[in]marker_chargeParticle weights time charge
[in,out]rho_dofsCoefficient vector of the charge distribution

Definition at line 79 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ add_current_evaluate_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::add_current_evaluate_spline_cl_1d ( class( sll_t_particle_mesh_coupling_spline_cl_1d ), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position_old,
real(kind=f64), dimension(self%dim), intent(in)  position_new,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), intent(in)  vbar,
real(kind=f64), dimension(self%n_dofs), intent(in)  field_dofs,
real(kind=f64), dimension(self%n_dofs), intent(inout)  j_dofs,
real(kind=f64), intent(out)  field 
)
private

Combines add_current and evaluate_int.

Parameters
[in,out]selfkernel smoother object
[in]position_oldPosition of the particle
[in]position_newPosition of the particle
[in]marker_chargeParticle weights time charge
[in]vbarParticle weights time charge
[in]field_dofsCoefficient vector of the current density
[in,out]j_dofsCoefficient vector of the current density
[out]fieldEfield

Definition at line 303 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ add_current_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::add_current_spline_cl_1d ( class( sll_t_particle_mesh_coupling_spline_cl_1d ), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position_old,
real(kind=f64), dimension(self%dim), intent(in)  position_new,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), dimension(self%n_dofs), intent(inout)  j_dofs 
)
private

Add current with integration over x.

Parameters
[in,out]selfkernel smoother object
[in]position_oldPosition of the particle
[in]position_newPosition of the particle
[in]marker_chargeParticle weights time charge
[in,out]j_dofsCoefficient vector of the current density

Definition at line 224 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ add_current_update_v_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::add_current_update_v_spline_cl_1d ( class(sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position_old,
real(kind=f64), dimension(self%dim), intent(in)  position_new,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), intent(in)  qoverm,
real(kind=f64), dimension(:), intent(in)  bfield_dofs,
real(kind=f64), dimension(:), intent(inout)  vi,
real(kind=f64), dimension(:), intent(inout)  j_dofs 
)
private

Add current for one particle and update v (according to H_p1 part in Hamiltonian splitting)

Parameters
[in,out]selfkernel smoother object
[in]position_oldPosition at time t
[in]position_newPosition at time t + \Delta t
[in]marker_chargeParticle weight time charge
[in]qovermcharge to mass ration
[in]bfield_dofsCoefficient of B-field expansion
[in,out]viVelocity of the particles
[in,out]j_dofsCoefficients of current expansion

Definition at line 392 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ add_particle_mass_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::add_particle_mass_spline_cl_1d ( class( sll_t_particle_mesh_coupling_spline_cl_1d ), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), dimension(:,:), intent(inout)  particle_mass 
)
private

Add the contribution of one particle to the approximate mass matrix.

Parameters
[in,out]selfkernel smoother object
[in]positionPosition of the particle
[in]marker_chargeParticle weights time charge
[in,out]particle_massCoefficient vector of the charge distribution

Definition at line 103 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ add_particle_mass_spline_cl_1d_full()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::add_particle_mass_spline_cl_1d_full ( class( sll_t_particle_mesh_coupling_spline_cl_1d ), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), dimension(:,:), intent(inout)  particle_mass 
)
private

Add the contribution of one particle to the approximate mass matrix.

Parameters
[in,out]selfkernel smoother object
[in]positionPosition of the particle
[in]marker_chargeParticle weights time charge
[in,out]particle_massCoefficient vector of the charge distribution

Definition at line 134 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ convert_x_to_xbox()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::convert_x_to_xbox ( class( sll_t_particle_mesh_coupling_spline_cl_1d ), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position,
real(kind=f64), intent(out)  xi,
integer(kind=i32), intent(out)  box 
)
private
Parameters
[in,out]selfkernel smoother object
[in]positionPosition of the particle
[out]xiPosition of the particle
[out]boxPosition of the particle

Definition at line 478 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the caller graph for this function:

◆ evaluate_field_single_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::evaluate_field_single_spline_cl_1d ( class (sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position,
real(kind=f64), dimension(self%n_dofs), intent(in)  field_dofs,
real(kind=f64), intent(out)  field_value 
)
private

Evaluate field at at position position.

Parameters
[in,out]selfKernel smoother object
[in]positionPosition of the particle
[in]field_dofsCoefficient vector for the field DoFs
[out]field_valueValue(s) of the electric fields at given position

Definition at line 165 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ evaluate_multiple_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::evaluate_multiple_spline_cl_1d ( class (sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self,
real(kind=f64), dimension(self%dim), intent(in)  position,
integer(kind=i32), dimension(:), intent(in)  components,
real(kind=f64), dimension(:,:), intent(in)  field_dofs,
real(kind=f64), dimension(:), intent(out)  field_value 
)
private

Evaluate several fields at position position.

Parameters
[in,out]selfKernel smoother object
[in]positionPosition of the particle
[in]componentsComponents of field_dofs that shall be updated
[in]field_dofsCoefficient vector for the field DoFs
[out]field_valueValue(s) of the electric fields at given position

Definition at line 192 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ free_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::free_spline_cl_1d ( class (sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self)
private

Finalizer.

Parameters
[in,out]selfKernel smoother object

Definition at line 601 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ init_spline_cl_1d()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::init_spline_cl_1d ( class( sll_t_particle_mesh_coupling_spline_cl_1d), intent(out)  self,
real(kind=f64), dimension(2), intent(in)  domain,
integer(kind=i32), intent(in)  n_cells,
integer(kind=i32), intent(in)  no_particles,
integer(kind=i32), intent(in)  spline_degree,
integer(kind=i32), intent(in)  smoothing_type,
integer(kind=i32), intent(in)  boundary 
)
private

Initializer.

Parameters
[out]selfkernel smoother object
[in]domainx_min and x_max of the domain
[in]n_cellsnumber of DoFs (spline coefficients)
[in]no_particlesnumber of particles
[in]spline_degreeDegree of smoothing kernel spline
[in]smoothing_typeDefine if Galerkin or collocation smoothing for right scaling in accumulation routines
[in]boundaryspline boundary conditions

Definition at line 536 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:

◆ sll_s_new_particle_mesh_coupling_spline_cl_1d()

subroutine, public sll_m_particle_mesh_coupling_spline_cl_1d::sll_s_new_particle_mesh_coupling_spline_cl_1d ( class( sll_c_particle_mesh_coupling_1d), intent(out), allocatable  smoother,
real(kind=f64), dimension(2), intent(in)  domain,
integer(kind=i32), intent(in)  n_cells,
integer(kind=i32), intent(in)  no_particles,
integer(kind=i32), intent(in)  spline_degree,
integer(kind=i32), intent(in)  smoothing_type,
integer(kind=i32), intent(in)  boundary 
)
Parameters
[out]smootherkernel smoother object
[in]n_cellsnumber of DoFs (spline coefficients)
[in]domainx_min and x_max of the domain
[in]no_particlesnumber of particles
[in]spline_degreeDegree of smoothing kernel spline
[in]smoothing_typeDefine if Galerkin or collocation smoothing for right scaling in accumulation routines
[in]boundaryspline boundary conditions

Definition at line 637 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

◆ sll_s_new_particle_mesh_coupling_spline_cl_1d_ptr()

subroutine, public sll_m_particle_mesh_coupling_spline_cl_1d::sll_s_new_particle_mesh_coupling_spline_cl_1d_ptr ( class( sll_c_particle_mesh_coupling_1d), intent(out), pointer  smoother,
real(kind=f64), dimension(2), intent(in)  domain,
integer(kind=i32), intent(in)  n_cells,
integer(kind=i32), intent(in)  no_particles,
integer(kind=i32), intent(in)  spline_degree,
integer(kind=i32), intent(in)  smoothing_type,
integer(kind=i32), intent(in)  boundary 
)
Parameters
[out]smootherkernel smoother object
[in]n_cellsnumber of DoFs (spline coefficients)
[in]domainx_min and x_max of the domain
[in]no_particlesnumber of particles
[in]spline_degreeDegree of smoothing kernel spline
[in]smoothing_typeDefine if Galerkin or collocation smoothing for right scaling in accumulation routines
[in]boundaryspline boundary conditions

Definition at line 612 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

◆ sll_s_uniform_bsplines_eval_basis_clamped()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::sll_s_uniform_bsplines_eval_basis_clamped ( type(sll_t_spline_pp_1d), intent(in)  spline,
integer(kind=i32), intent(in)  n_cells,
integer(kind=i32), intent(in)  degree,
real(kind=f64), intent(in)  xi,
integer(kind=i32), intent(in)  box,
real(kind=f64), dimension(:), intent(out)  spline_val 
)
private

Helper function to evaluate uniform clamped basis splines.

Definition at line 507 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

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

◆ update_j()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::update_j ( class(sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self,
real(kind=f64), intent(in)  lower,
real(kind=f64), intent(in)  upper,
integer(kind=i32), intent(in)  index,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), intent(in)  sign,
real(kind=f64), dimension(:), intent(inout)  j_dofs 
)
private

Helper function for add_current.

Parameters
[in,out]selftime splitting object

Definition at line 269 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

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

◆ update_je()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::update_je ( class(sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self,
real(kind=f64), intent(in)  lower,
real(kind=f64), intent(in)  upper,
integer(kind=i32), intent(in)  index,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), intent(in)  sign,
real(kind=f64), dimension(:), intent(inout)  j_dofs,
real(kind=f64), dimension(:), intent(in)  field_dofs,
real(kind=f64), intent(inout)  field 
)
private

Helper function for add_current_evaluate.

Parameters
[in,out]selftime splitting object

Definition at line 354 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

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

◆ update_jv()

subroutine sll_m_particle_mesh_coupling_spline_cl_1d::update_jv ( class(sll_t_particle_mesh_coupling_spline_cl_1d), intent(inout)  self,
real(kind=f64), intent(in)  lower,
real(kind=f64), intent(in)  upper,
integer(kind=i32), intent(in)  index,
real(kind=f64), intent(in)  marker_charge,
real(kind=f64), intent(in)  qoverm,
real(kind=f64), intent(in)  sign,
real(kind=f64), intent(inout)  vi,
real(kind=f64), dimension(:), intent(inout)  j_dofs,
real(kind=f64), dimension(:), intent(in)  bfield_dofs 
)
private

Helper function for add_current_update_v.

Parameters
[in,out]selftime splitting object

Definition at line 440 of file sll_m_particle_mesh_coupling_spline_cl_1d.F90.

Here is the call graph for this function:
Here is the caller graph for this function:
    Report Typos and Errors