sll_m_particle_sampling_interface.F90

Go to the documentation of this file.

!**************************************************************
!  Copyright INRIA
!  Authors :
!     CALVI project team
!  
!  This code SeLaLib (for Semi-Lagrangian-Library) 
!  is a parallel library for simulating the plasma turbulence 
!  in a tokamak.
!  
!  This software is governed by the CeCILL-B license 
!  under French law and abiding by the rules of distribution 
!  of free software.  You can  use, modify and redistribute 
!  the software under the terms of the CeCILL-B license as 
!  circulated by CEA, CNRS and INRIA at the following URL
!  "http://www.cecill.info". 
!**************************************************************
 
 
module sll_m_particle_sampling_interface
 
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#include "sll_assert.h"
#include "sll_errors.h"
#include "sll_memory.h"
#include "sll_working_precision.h"
 
  use sll_m_particle_group_base, only: &
    sll_c_particle_group_base
 
  use sll_m_particle_sampling, only: &
       sll_t_particle_sampling
 
  use sll_m_particle_group_2d2v_lbf, only: &
    sll_t_particle_group_2d2v_lbf
 
  use sll_m_initial_distribution, only : &
       sll_c_distribution_params
 
  use sll_m_control_variate, only : &
       sll_t_control_variate
 
  implicit none
 
  public :: &
    sll_t_conservative_sampling_params, &
    sll_s_sample_particle_group, &
    sll_s_resample_particle_group
 
  private
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
  type sll_t_conservative_sampling_params
      sll_real64 :: total_charge      
    contains
      procedure :: set_total_charge
  end type sll_t_conservative_sampling_params
 
 
contains
 
  subroutine set_total_charge(self, target_total_charge)
    class(sll_t_conservative_sampling_params),  intent( inout )         :: self
    sll_real64,                                     intent( in )  :: target_total_charge
    self%total_charge = target_total_charge
  end subroutine
 
  !------------------------------------------------------------------------------------------
 
  subroutine sll_s_sample_particle_group( &
          particle_group,   &
          distribution_params, &
          random_sampler, &
          nb_weights, &
          control_variate, &
          conservative_sampling_params, &
          xmin, &
          Lx )
 
    class(sll_c_particle_group_base), pointer,  intent( inout )        :: particle_group
    class(sll_c_distribution_params),           intent( inout )           :: distribution_params
    type(sll_t_particle_sampling),              intent( inout ), optional :: random_sampler
    sll_int32,                                  intent( in ), optional :: nb_weights
    class(sll_t_control_variate),               intent( in ), optional :: control_variate
    class(sll_t_conservative_sampling_params),  intent( in ), optional :: conservative_sampling_params
    sll_real64,                                 intent( in ), optional :: xmin(:)  
    sll_real64,                                 intent( in ), optional :: lx(:)    
 
    sll_real64 :: target_total_charge
    logical    :: enforce_total_charge
 
    select type ( particle_group )
 
      type is ( sll_t_particle_group_2d2v_lbf )
 
        if( present(conservative_sampling_params) )then
          enforce_total_charge = .true.
          target_total_charge = conservative_sampling_params%total_charge
        else
          enforce_total_charge = .false.    ! no charge conservation
          target_total_charge = 0._f64      ! value does not matter then
        end if
        call particle_group%sample( target_total_charge, enforce_total_charge, distribution_params )
 
      class default
 
        sll_assert( present( random_sampler ) )
        sll_assert( present( xmin ) )
        sll_assert( present( lx ) )
        call random_sampler%sample( particle_group, distribution_params, xmin, lx )
 
 
        if (nb_weights == 1 ) then
           call random_sampler%sample ( particle_group, distribution_params, xmin, lx )
        elseif ( nb_weights == 3 ) then
           sll_assert( present( control_variate ) )
           call random_sampler%sample_cv ( particle_group, distribution_params, xmin, lx, control_variate )
        else
           sll_error("sll_s_sample_particle_group", "random sampling interface not implemented for this nb of weights")
        end if
 
      end select
 
  end subroutine sll_s_sample_particle_group
 
 
 
  subroutine sll_s_resample_particle_group( &
          particle_group,   &
          conservative_sampling_params )      !< whether charge must be conserved
 
    class(sll_c_particle_group_base),    pointer, intent( inout )         :: particle_group
    class(sll_t_conservative_sampling_params), intent( in ), optional :: conservative_sampling_params
    sll_real64 :: target_total_charge
    logical    :: enforce_total_charge
 
    select type ( particle_group )
 
    type is ( sll_t_particle_group_2d2v_lbf )
      if( present(conservative_sampling_params) )then
        enforce_total_charge = .true.
        target_total_charge = conservative_sampling_params%total_charge
      else
        enforce_total_charge = .false.    ! no charge conservation
        target_total_charge = 0._f64      ! value does not matter then
      end if
      call particle_group%resample( target_total_charge, enforce_total_charge )
 
    class default
      sll_error("sll_s_resample_particle_group", "resampling interface not implemented for this type of particle group")
 
    end select
 
  end subroutine sll_s_resample_particle_group
 
 
 
end module  sll_m_particle_sampling_interface