Report Typos and Errors    
Semi-Lagrangian Library
Modular library for kinetic and gyrokinetic simulations of plasmas in fusion energy devices.
sll_m_particle_group_2d2v_lbf.F90
Go to the documentation of this file.
1 !**************************************************************
2 ! Copyright INRIA
3 ! Authors :
4 ! MCP
5 !
6 ! This code SeLaLib (for Semi-Lagrangian-Library)
7 ! is a parallel library for simulating the plasma turbulence
8 ! in a tokamak.
9 !
10 ! This software is governed by the CeCILL-B license
11 ! under French law and abiding by the rules of distribution
12 ! of free software. You can use, modify and redistribute
13 ! the software under the terms of the CeCILL-B license as
14 ! circulated by CEA, CNRS and INRIA at the following URL
15 ! "http://www.cecill.info".
16 !**************************************************************
17 
19 
21 
23 
24 module sll_m_particle_group_2d2v_lbf
25 
26 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
27 #include "sll_assert.h"
28 #include "sll_errors.h"
29 #include "sll_memory.h"
30 #include "sll_working_precision.h"
31 
32  use sll_m_cartesian_meshes, only: &
33  sll_f_new_cartesian_mesh_4d, &
34  sll_t_cartesian_mesh_4d
35 
36  use sll_m_particle_group_base, only: &
37  sll_c_particle_group_base, &
38  sll_t_species
39 
40  use sll_m_constants, only: &
41  sll_p_pi
42 
43  use sll_m_initial_distribution, only: &
44  sll_c_distribution_params
45 
46  use sll_m_sparse_grid_4d, only: &
47  sll_t_sparse_grid_4d
48 
49  implicit none
50 
51  public :: &
52  sll_t_particle_group_2d2v_lbf, &
53  sll_s_new_particle_group_2d2v_lbf_ptr
54 
55  private
56 
57 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
58 
60  sll_int32, parameter :: sll_p_lbf_remapping_grid = 1
61  sll_int32, parameter :: sll_p_lbf_given_grid = 2
62 
64  type sll_t_int_list_element
65  sll_int32 :: value
66  type(sll_t_int_list_element), pointer :: next
67  end type sll_t_int_list_element
68 
70  type sll_t_int_list_element_ptr
71  type(sll_t_int_list_element), pointer :: pointed_element
72  end type sll_t_int_list_element_ptr
73 
74 
76  type, extends(sll_c_particle_group_base) :: sll_t_particle_group_2d2v_lbf
77 
80  sll_int32 :: n_particles_x
81  sll_int32 :: n_particles_y
82  sll_int32 :: n_particles_vx
83  sll_int32 :: n_particles_vy
84  sll_real64, allocatable :: particle_array(:,:)
85  sll_real64 :: common_weight
87 
94  type(sll_t_cartesian_mesh_4d), pointer :: lbf_grid
95  logical :: domain_is_periodic(2)
97 
98 
101  sll_int32 :: remapped_f_interpolation_degree
102  type(sll_t_sparse_grid_4d) :: sparse_grid_interpolator
103  sll_int32, dimension(4) :: sparse_grid_max_levels
104  sll_real64, dimension(:), allocatable :: remapped_f_sparse_grid_coefficients
105  sll_real64, dimension(:), allocatable :: tmp_f_values_on_remapping_sparse_grid
107 
108  contains
109 
112  procedure :: get_x => get_x_2d2v_lbf
113  procedure :: get_v => get_v_2d2v_lbf
114  procedure :: get_charge => get_charge_2d2v_lbf
115  procedure :: get_mass => get_mass_2d2v_lbf
116  procedure :: get_weights => get_weights_2d2v_lbf
117  procedure :: get_common_weight => get_common_weight_2d2v_lbf
119 
122  procedure :: set_x => set_x_2d2v_lbf
123  procedure :: set_v => set_v_2d2v_lbf
124  procedure :: set_weights => set_weights_2d2v_lbf
125  procedure :: set_common_weight => set_common_weight_2d2v_lbf
127 
129  procedure :: init => initialize_particle_group_2d2v_lbf
130  procedure :: free => delete_particle_group_2d2v_lbf
131 
134  ! procedure :: write_hat_density_on_remapping_grid !> this evaluates an analytic f0
135  ! procedure :: write_landau_density_on_remapping_grid !> this evaluates an analytic f0
136  !
137  procedure, private :: write_known_density_on_remapping_grid
138  procedure, private :: reset_particles_positions
139  procedure, private :: reset_particles_weights_with_direct_interpolation
140  procedure, private :: reconstruct_f_lbf
141  procedure, private :: reconstruct_f_lbf_on_remapping_grid
142  procedure :: reconstruct_f_lbf_on_given_grid
143  procedure :: sample
144  procedure :: resample
146 
147  procedure, private :: interpolate_value_of_remapped_f
148  procedure, private :: get_ltp_deformation_matrix
149  procedure, private :: get_neighbor_index
150 
151  end type sll_t_particle_group_2d2v_lbf
152 
153 
154 contains
155 
157 
158  !----------------------------------------------------------------------------------------------------------------------------
160  pure function get_x_2d2v_lbf( self, i ) result( r )
161  class( sll_t_particle_group_2d2v_lbf ), intent( in ) :: self
162  sll_int32 , intent( in ) :: i
163  sll_real64 :: r(3)
164 
165  r = 1.0_f64
166  r(1:2) = self%particle_array(1:2, i)
167 
168  end function get_x_2d2v_lbf
169 
170  !----------------------------------------------------------------------------------------------------------------------------
172  pure function get_v_2d2v_lbf( self, i ) result( r )
173  class( sll_t_particle_group_2d2v_lbf ), intent( in ) :: self
174  sll_int32 , intent( in ) :: i
175  sll_real64 :: r(3)
176 
177  r = 1.0_f64
178  r(1:2) = self%particle_array(3:4, i)
179 
180  end function get_v_2d2v_lbf
181 
182 
183 
184  !----------------------------------------------------------------------!
186  pure function get_charge_2d2v_lbf( self, i , i_weight) result (r)
187  class( sll_t_particle_group_2d2v_lbf ), intent( in ) :: self
188  sll_int32 , intent( in ) :: i
189  sll_int32, optional , intent( in ) :: i_weight
190  sll_real64 :: r
191 
192  sll_int32 :: i_wi
193 
194  i_wi = 1
195  if(present(i_weight)) i_wi = i_weight
196  r = self%species%q * self%particle_array(4+i_wi, i) * self%common_weight
197 
198  end function get_charge_2d2v_lbf
199 
200 
201  !----------------------------------------------------------------------!
203  pure function get_mass_2d2v_lbf( self, i, i_weight) result (r)
204  class( sll_t_particle_group_2d2v_lbf ), intent( in ) :: self
205  sll_int32 , intent( in ) :: i
206  sll_int32, optional , intent( in ) :: i_weight
207  sll_real64 :: r
208 
209  sll_int32 :: i_wi
210 
211  i_wi = 1
212  if(present(i_weight)) i_wi = i_weight
213  r = self%species%m * self%particle_array(4+i_wi, i) * self%common_weight
214 
215  end function get_mass_2d2v_lbf
216 
217 
218  !----------------------------------------------------------------------!
220  pure function get_weights_2d2v_lbf( self, i) result (r)
221  class( sll_t_particle_group_2d2v_lbf ), intent( in ) :: self
222  sll_int32 , intent( in ) :: i
223  sll_real64 :: r(self%n_weights)
224 
225  r = self%particle_array(5:4+self%n_weights, i)
226 
227  end function get_weights_2d2v_lbf
228 
229  !----------------------------------------------------------------------!
231  pure function get_common_weight_2d2v_lbf( self )result(r)
232  class( sll_t_particle_group_2d2v_lbf ), intent( in ) :: self
233  sll_real64 :: r
234 
235  r = self%common_weight
236 
237  end function get_common_weight_2d2v_lbf
238 
240 
241  !----------------------------------------------------------------------!
243  subroutine set_x_2d2v_lbf( self, i, x )
244  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
245  sll_int32 , intent( in ) :: i
246  sll_real64 , intent( in ) :: x(3)
247 
248  self%particle_array(1:2, i) = x(1:2)
249 
250  end subroutine set_x_2d2v_lbf
251 
252  !----------------------------------------------------------------------!
254  subroutine set_v_2d2v_lbf( self, i, x )
255  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
256  sll_int32 , intent( in ) :: i
257  sll_real64 , intent( in ) :: x(3)
258 
259  self%particle_array(3:4, i) = x(1:2)
260 
261  end subroutine set_v_2d2v_lbf
262 
263  !----------------------------------------------------------------------!
265  subroutine set_weights_2d2v_lbf( self, i, x )
266  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
267  sll_int32 , intent( in ) :: i
268  sll_real64 , intent( in ) :: x(self%n_weights)
269 
270  self%particle_array(5:4+self%n_weights, i) = x
271 
272  end subroutine set_weights_2d2v_lbf
273 
274  !----------------------------------------------------------------------!
276  subroutine set_common_weight_2d2v_lbf( self, x )
277  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
278  sll_real64 , intent( in ) :: x
279 
280  self%common_weight = x
281 
282  end subroutine set_common_weight_2d2v_lbf
283 
284 
287  subroutine initialize_particle_group_2d2v_lbf ( &
288  self, &
289  species_charge, &
290  species_mass, &
291  domain_is_x_periodic, &
292  domain_is_y_periodic, &
293  remap_degree, &
294  remapping_grid_eta_min, &
295  remapping_grid_eta_max, &
296  remapping_sparse_grid_max_levels, &
297  n_particles_x, &
298  n_particles_y, &
299  n_particles_vx, &
300  n_particles_vy, &
301  n_weights &
302  )
303 
304  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
305 
306  sll_real64, intent(in) :: species_charge
307  sll_real64, intent(in) :: species_mass
308  logical, intent(in) :: domain_is_x_periodic
309  logical, intent(in) :: domain_is_y_periodic
310  sll_int32, intent(in) :: remap_degree
311  sll_real64, dimension(4), intent(in) :: remapping_grid_eta_min
312  sll_real64, dimension(4), intent(in) :: remapping_grid_eta_max
313  sll_int32, dimension(4), intent(in) :: remapping_sparse_grid_max_levels
314  sll_int32, intent(in) :: n_particles_x
315  sll_int32, intent(in) :: n_particles_y
316  sll_int32, intent(in) :: n_particles_vx
317  sll_int32, intent(in) :: n_particles_vy
318  sll_int32, intent(in) :: n_weights
319 
320  character(len=*), parameter :: this_fun_name = "initialize_particle_group_2d2v_lbf"
321  sll_int32 :: ierr
322 
323  print*, "[", this_fun_name, "] - initializing the particle group... "
324 
325  self%domain_is_periodic(1) = domain_is_x_periodic
326  self%domain_is_periodic(2) = domain_is_y_periodic
327 
328  self%n_particles_x = n_particles_x
329  self%n_particles_y = n_particles_y
330  self%n_particles_vx = n_particles_vx
331  self%n_particles_vy = n_particles_vy
332  self%n_particles = n_particles_x * n_particles_y * n_particles_vx * n_particles_vy
333  self%n_total_particles = self%n_particles
334  self%n_weights = n_weights
335 
337  allocate(self%species, stat=ierr)
338  sll_assert( ierr == 0)
339  call self%species%init( species_charge, species_mass )
340 
341  sll_allocate( self%particle_array(4+n_weights, self%n_particles), ierr)
342  sll_assert( ierr == 0)
343 
347  self%lbf_grid => sll_f_new_cartesian_mesh_4d( &
348  n_particles_x, &
349  n_particles_y, &
350  n_particles_vx, &
351  n_particles_vy, &
352  remapping_grid_eta_min(1), &
353  remapping_grid_eta_max(1), &
354  remapping_grid_eta_min(2), &
355  remapping_grid_eta_max(2), &
356  remapping_grid_eta_min(3), &
357  remapping_grid_eta_max(3), &
358  remapping_grid_eta_min(4), &
359  remapping_grid_eta_max(4) )
360 
362  self%remapped_f_interpolation_degree = remap_degree
363 
364  print*, "[", this_fun_name, "] - sparse grid levels for the remapping tool:", remapping_sparse_grid_max_levels
365  self%sparse_grid_max_levels = remapping_sparse_grid_max_levels
366  call self%sparse_grid_interpolator%init( &
367  self%sparse_grid_max_levels, &
368  self%remapped_f_interpolation_degree, &
369  self%remapped_f_interpolation_degree+1, &
370  0, &
371  remapping_grid_eta_min, &
372  remapping_grid_eta_max &
373  )
374 
376  sll_allocate( self%remapped_f_sparse_grid_coefficients(self%sparse_grid_interpolator%size_basis), ierr )
377  self%remapped_f_sparse_grid_coefficients = 0.0_f64
378 
380  sll_allocate( self%tmp_f_values_on_remapping_sparse_grid(self%sparse_grid_interpolator%size_basis), ierr)
381 
382  end subroutine initialize_particle_group_2d2v_lbf
383 
384  !-----------------------------------------------------------------------------------------------------------------------------
386  subroutine sll_s_new_particle_group_2d2v_lbf_ptr(&
387  particle_group, &
388  species_charge, &
389  species_mass, &
390  domain_is_x_periodic, &
391  domain_is_y_periodic, &
392  remap_degree, &
393  remapping_grid_eta_min, &
394  remapping_grid_eta_max, &
395  remapping_sparse_grid_max_levels, &
396  n_particles_x, &
397  n_particles_y, &
398  n_particles_vx, &
399  n_particles_vy &
400  )
401  class( sll_c_particle_group_base ), pointer, intent( out ) :: particle_group
402  sll_real64, intent(in) :: species_charge
403  sll_real64, intent(in) :: species_mass
404  logical, intent(in) :: domain_is_x_periodic
405  logical, intent(in) :: domain_is_y_periodic
406  sll_int32, intent(in) :: remap_degree
407  sll_real64, dimension(4), intent(in) :: remapping_grid_eta_min
408  sll_real64, dimension(4), intent(in) :: remapping_grid_eta_max
409  sll_int32, dimension(4), intent(in) :: remapping_sparse_grid_max_levels
410  sll_int32, intent(in) :: n_particles_x
411  sll_int32, intent(in) :: n_particles_y
412  sll_int32, intent(in) :: n_particles_vx
413  sll_int32, intent(in) :: n_particles_vy
414 
415  sll_int32 :: ierr
416  sll_int32 :: n_weights
417 
418  sll_allocate( sll_t_particle_group_2d2v_lbf :: particle_group, ierr)
419  n_weights = 1
420 
421  select type( particle_group )
422  type is ( sll_t_particle_group_2d2v_lbf )
423  call particle_group%init( &
424  species_charge, &
425  species_mass, &
426  domain_is_x_periodic, &
427  domain_is_y_periodic, &
428  remap_degree, &
429  remapping_grid_eta_min, &
430  remapping_grid_eta_max, &
431  remapping_sparse_grid_max_levels, &
432  n_particles_x, &
433  n_particles_y, &
434  n_particles_vx, &
435  n_particles_vy, &
436  n_weights &
437  )
438  end select
439 
440  end subroutine sll_s_new_particle_group_2d2v_lbf_ptr
441 
443  subroutine delete_particle_group_2d2v_lbf( self )
444  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
445 
446  deallocate(self%particle_array)
447  deallocate(self%remapped_f_sparse_grid_coefficients)
448  deallocate(self%tmp_f_values_on_remapping_sparse_grid)
449 
450  end subroutine
451 
454 
455 
456 
461  subroutine sample( self, target_total_charge, enforce_total_charge, init_f_params, rand_seed, rank, world_size )
462  class(sll_t_particle_group_2d2v_lbf), intent( inout ) :: self
463  sll_real64, intent( in ) :: target_total_charge
464  logical, intent( in ) :: enforce_total_charge
465  class(sll_c_distribution_params), intent( in ) :: init_f_params
466  sll_int32, dimension(:) , intent( in ), optional :: rand_seed
467  sll_int32 , intent( in ), optional :: rank, world_size
468 
469  if( present(rand_seed) )then
470  sll_assert( present( rank ) )
471  sll_assert( present( world_size ) )
472  call self%resample( target_total_charge, enforce_total_charge, init_f_params, rand_seed, rank, world_size )
473  else
474  call self%resample( target_total_charge, enforce_total_charge, init_f_params )
475  end if
476  end subroutine
477 
483  subroutine resample( self, target_total_charge, enforce_total_charge, init_f_params, rand_seed, rank, world_size )
484  class(sll_t_particle_group_2d2v_lbf), intent( inout ) :: self
485  sll_real64, intent( in ) :: target_total_charge
486  logical, intent( in ) :: enforce_total_charge
487  class(sll_c_distribution_params), intent( in ), optional :: init_f_params ! < use for initial sampling
488  sll_int32, dimension(:) , intent( in ), optional :: rand_seed
489  sll_int32 , intent( in ), optional :: rank, world_size
490 
491  logical :: initial_step
492 
493  initial_step = present(init_f_params)
494 
496  print *, "particle_group_2d2v_lbf%resample -- step A (LBF approximation of transported density on remapping grid)"
497 
499  if( initial_step )then
501  call self%write_known_density_on_remapping_grid(init_f_params)
502 
503  else
505  call self%reconstruct_f_lbf_on_remapping_grid()
506  end if
507 
509  call self%sparse_grid_interpolator%compute_hierarchical_surplus( self%tmp_f_values_on_remapping_sparse_grid )
510  self%remapped_f_sparse_grid_coefficients = self%tmp_f_values_on_remapping_sparse_grid
511 
513  print *, "particle_group_2d2v_lbf%resample -- step B (deterministic resampling of particles using the remapped density)"
514  call self%reset_particles_positions
516  call self%reset_particles_weights_with_direct_interpolation( target_total_charge, enforce_total_charge )
517 
519  if(present(rank))then
520  sll_assert(present(world_size))
521  sll_assert(present(rand_seed))
522  end if
523 
524  end subroutine resample
525 
528  subroutine write_known_density_on_remapping_grid( &
529  self, &
530  distribution_params &
531  )
532 
533  class(sll_t_particle_group_2d2v_lbf), intent(inout) :: self
534  class(sll_c_distribution_params), intent(in) :: distribution_params
535 
536  sll_int32 :: j
537  sll_real64 :: x_j(1:2)
538  sll_real64 :: v_j(1:2)
539 
540  sll_assert( size(self%tmp_f_values_on_remapping_sparse_grid,1) == self%sparse_grid_interpolator%size_basis )
541  self%tmp_f_values_on_remapping_sparse_grid = 0.0_f64
542 
543  do j=1, self%sparse_grid_interpolator%size_basis
544  x_j = self%sparse_grid_interpolator%hierarchy(j)%coordinate(1:2)
545  v_j = self%sparse_grid_interpolator%hierarchy(j)%coordinate(3:4)
546  self%tmp_f_values_on_remapping_sparse_grid(j) = distribution_params%eval_xv_density( x_j, v_j )
547  end do
548 
549  end subroutine write_known_density_on_remapping_grid
550 
551 
556  subroutine reset_particles_positions( self )
557  class(sll_t_particle_group_2d2v_lbf),intent(inout) :: self
558 
559  sll_int32 :: k
560  sll_int32 :: k_check
561  sll_real64, dimension(3) :: coords
562 
563  sll_int32 :: j_x
564  sll_int32 :: j_y
565  sll_int32 :: j_vx
566  sll_int32 :: j_vy
567  sll_real64 :: h_x
568  sll_real64 :: h_y
569  sll_real64 :: h_vx
570  sll_real64 :: h_vy
571  sll_real64 :: x_min
572  sll_real64 :: y_min
573  sll_real64 :: vx_min
574  sll_real64 :: vy_min
575  sll_real64 :: x_j
576  sll_real64 :: y_j
577  sll_real64 :: vx_j
578  sll_real64 :: vy_j
579 
580  sll_assert( self%n_particles_x == self%lbf_grid%num_cells1 )
581  sll_assert( self%n_particles_y == self%lbf_grid%num_cells2 )
582  sll_assert( self%n_particles_vx == self%lbf_grid%num_cells3 )
583  sll_assert( self%n_particles_vy == self%lbf_grid%num_cells4 )
584 
585  h_x = self%lbf_grid%delta_eta1
586  h_y = self%lbf_grid%delta_eta2
587  h_vx = self%lbf_grid%delta_eta3
588  h_vy = self%lbf_grid%delta_eta4
589 
590  x_min = self%lbf_grid%eta1_min
591  y_min = self%lbf_grid%eta2_min
592  vx_min = self%lbf_grid%eta3_min
593  vy_min = self%lbf_grid%eta4_min
594 
595  k_check = 0
596  do j_x = 1, self%n_particles_x
597  x_j = x_min + (j_x-0.5) * h_x
598 
599  do j_y = 1, self%n_particles_y
600  y_j = y_min + (j_y-0.5) * h_y
601 
602  do j_vx = 1, self%n_particles_vx
603  vx_j = vx_min + (j_vx-0.5) * h_vx
604 
605  do j_vy = 1, self%n_particles_vy
606  vy_j = vy_min + (j_vy-0.5) * h_vy
607 
608  k_check = k_check + 1
609  call convert_4d_index_to_1d( &
610  k, &
611  j_x, j_y, j_vx, j_vy, &
612  self%n_particles_x, self%n_particles_y, self%n_particles_vx, self%n_particles_vy &
613  )
614  sll_assert(k == k_check)
615 
616  coords(1) = x_j
617  coords(2) = y_j
618  call self%set_x( k, coords )
619 
620  coords(1) = vx_j
621  coords(2) = vy_j
622  call self%set_v( k, coords )
623 
624  end do
625  end do
626  end do
627  end do
628 
629  end subroutine reset_particles_positions
630 
633  function get_neighbor_index( self, k, dim, dir ) result(k_ngb)
634  class(sll_t_particle_group_2d2v_lbf),intent(inout) :: self
635 
636  sll_int32, intent(in) :: k
637  sll_int32, intent(in) :: dim
638  sll_int32, intent(in) :: dir
639  sll_int32 :: k_ngb
640 
641  sll_int32 :: j_x
642  sll_int32 :: j_y
643  sll_int32 :: j_vx
644  sll_int32 :: j_vy
645  sll_int32 :: j_ngb
646 
647  call convert_1d_index_to_4d( &
648  j_x, j_y, j_vx, j_vy, &
649  k, &
650  self%n_particles_x, self%n_particles_y, &
651  self%n_particles_vx, self%n_particles_vy &
652  )
653 
654  select case( dim )
655  case( 1 )
656  j_ngb = j_x + dir
657  ! correct neighbor index if out of bounds
658  if(j_ngb < 1 .or. j_ngb > self%n_particles_x )then
659  if( self%domain_is_periodic(1) )then
660  j_ngb = modulo(j_ngb-1, self%n_particles_x)+1
661  else
662  sll_assert( (j_x == 1 .and. dir == -1) .or. (j_x == self%n_particles_x .and. dir == 1) )
663  k_ngb = k ! no neighbor, return the index itself
664  return
665  sll_error( "get_neighbor_index", "CHECK 1 -- should not be here")
666  end if
667  end if
668  call convert_4d_index_to_1d( &
669  k_ngb, &
670  j_ngb, j_y, j_vx, j_vy, &
671  self%n_particles_x, self%n_particles_y, &
672  self%n_particles_vx, self%n_particles_vy &
673  )
674  return
675 
676  case( 2 )
677  j_ngb = j_y + dir
678  ! correct neighbor index if out of bounds
679  if(j_ngb < 1 .or. j_ngb > self%n_particles_y )then
680  if( self%domain_is_periodic(2) )then
681  j_ngb = modulo(j_ngb-1, self%n_particles_y)+1
682  else
683  sll_assert( (j_y == 1 .and. dir == -1) .or. (j_y == self%n_particles_y .and. dir == 1) )
684  k_ngb = k ! no neighbor, return the index itself
685  return
686  sll_error( "get_neighbor_index", "CHECK 2 -- should not be here")
687  end if
688  end if
689  call convert_4d_index_to_1d( &
690  k_ngb, &
691  j_x, j_ngb, j_vx, j_vy, &
692  self%n_particles_x, self%n_particles_y, &
693  self%n_particles_vx, self%n_particles_vy &
694  )
695  return
696 
697  case( 3 )
698  j_ngb = j_vx + dir
699  ! correct neighbor index if out of bounds
700  if(j_ngb < 1 .or. j_ngb > self%n_particles_vx )then
701  ! no periodicity along vx
702  sll_assert( (j_vx == 1 .and. dir == -1) .or. (j_vx == self%n_particles_vx .and. dir == 1) )
703  k_ngb = k ! no neighbor, return the index itself
704  return
705  sll_error( "get_neighbor_index", "CHECK 3 -- should not be here")
706  end if
707  call convert_4d_index_to_1d( &
708  k_ngb, &
709  j_x, j_y, j_ngb, j_vy, &
710  self%n_particles_x, self%n_particles_y, &
711  self%n_particles_vx, self%n_particles_vy &
712  )
713  return
714 
715  case( 4 )
716  j_ngb = j_vy + dir
717  ! correct neighbor index if out of bounds
718  if(j_ngb < 1 .or. j_ngb > self%n_particles_vy )then
719  ! no periodicity along vy
720  sll_assert( (j_vy == 1 .and. dir == -1) .or. (j_vy == self%n_particles_vy .and. dir == 1) )
721  k_ngb = k ! no neighbor, return the index itself
722  return
723  sll_error( "get_neighbor_index", "CHECK 4 -- should not be here")
724  end if
725  call convert_4d_index_to_1d( &
726  k_ngb, &
727  j_x, j_y, j_vx, j_ngb, &
728  self%n_particles_x, self%n_particles_y, &
729  self%n_particles_vx, self%n_particles_vy &
730  )
731  return
732 
733  case default
734  sll_error("get_neighbor_index", "wrong value for dim")
735  end select
736 
737  end function
738 
742  subroutine reset_particles_weights_with_direct_interpolation( &
743  self, &
744  target_total_charge, &
745  enforce_total_charge &
746  )
747  class( sll_t_particle_group_2d2v_lbf ), intent( inout ) :: self
748  sll_real64, intent( in ) :: target_total_charge
749  logical, intent( in ) :: enforce_total_charge
750  sll_real64 :: eta(4)
751  sll_int32 :: i_part
752  sll_real64 :: phase_space_volume
753  sll_real64 :: particle_density_factor
754  sll_real64 :: point_density, total_computed_density
755  sll_real64 :: total_computed_charge
756  sll_real64 :: charge_correction_factor
757 
759  call self%set_common_weight(1.0_f64)
760 
761  phase_space_volume = (self%lbf_grid%eta4_max - self%lbf_grid%eta4_min) &
762  * (self%lbf_grid%eta3_max - self%lbf_grid%eta3_min) &
763  * (self%lbf_grid%eta2_max - self%lbf_grid%eta2_min) &
764  * (self%lbf_grid%eta1_max - self%lbf_grid%eta1_min)
765  particle_density_factor = phase_space_volume / self%n_particles
766  total_computed_density = 0.0d0
767  do i_part = 1, self%n_particles
768  eta = self%particle_array(1:4, i_part)
769  point_density = particle_density_factor * self%interpolate_value_of_remapped_f(eta)
770  self%particle_array(5, i_part) = point_density
771  total_computed_density = total_computed_density + point_density
772  end do
773 
774  if( enforce_total_charge )then
775  total_computed_charge = self%species%q * total_computed_density
776  if( total_computed_density == 0 )then
777  print *, "WARNING (876786587689) -- total computed_density is zero, which is strange..."
778  print *, " -- (no charge correction in this case) "
779  else
780  charge_correction_factor = target_total_charge / total_computed_charge
781 
782  print *, "[Enforcing charge in reset_particles_weights_with_direct_interpolation] ... "
783  print *, " ... target_total_charge, total_computed_charge, charge_correction_factor = ", &
784  target_total_charge, total_computed_charge, charge_correction_factor
785 
786  ! todo: try with self%particle_array(5, :) = self%particle_array(5, :) * charge_correction_factor
787  do i_part = 1, self%n_particles
788  self%particle_array(5, i_part) = self%particle_array(5, i_part) * charge_correction_factor
789  end do
790  end if
791  end if
792 
793  end subroutine reset_particles_weights_with_direct_interpolation
794 
798  function interpolate_value_of_remapped_f ( self, eta ) result(val)
799  class(sll_t_particle_group_2d2v_lbf), intent(inout) :: self
800  sll_real64, dimension(4), intent(in) :: eta
801  sll_real64 :: val
802 
803  val = self%sparse_grid_interpolator%interpolate_from_interpolant_value(self%remapped_f_sparse_grid_coefficients, eta)
804 
805  end function
806 
808 #define UPDATE_CLOSEST_PARTICLE_ARRAYS_USING_NEIGHBOR_CELLS(djx,djy,djvx,djvy) \
809  do; \
810  k_neighbor = closest_particle(j_x+(djx), j_y+(djy), j_vx+(djvx), j_vy+(djvy)); \
811 ; \
812  if(k_neighbor /= 0) then; \
813  do; \
814  coords = self%get_x(k_neighbor) ; \
815  x = coords(1) ; \
816  y = coords(2) ; \
817  coords = self%get_v(k_neighbor) ; \
818  vx = coords(1) ; \
819  vy = coords(2) ; \
820  call periodic_correction(self,x,y) ; \
821  call update_closest_particle_arrays(k_neighbor, \
822  x - self%lbf_grid%eta1_min, y - self%lbf_grid%eta2_min, vx - self%lbf_grid%eta3_min, vy - self%lbf_grid%eta4_min, \
823  j_x, j_y, j_vx, j_vy, \
824  h_flow_grid_x, \
825  h_flow_grid_y, \
826  h_flow_grid_vx, \
827  h_flow_grid_vy, \
828  closest_particle, \
829  closest_particle_distance) ; \
830  exit; \
831  end do; \
832  end if; \
833  exit; \
834  end do
835 
836 
837 
843 
844  subroutine reconstruct_f_lbf( &
845  self, &
846  reconstruction_set_type, &
847  given_grid_4d, &
848  given_array_4d, &
849  reconstruct_f_on_last_node, &
850  target_total_charge, &
851  enforce_total_charge &
852  )
853  class(sll_t_particle_group_2d2v_lbf), intent(inout) :: self
854  sll_int32, intent(in) :: reconstruction_set_type
855  type(sll_t_cartesian_mesh_4d), pointer, intent(in) :: given_grid_4d
856  sll_real64, dimension(:,:,:,:), pointer, intent(inout) :: given_array_4d
857  logical, intent(in) :: reconstruct_f_on_last_node(4)
858  sll_real64, intent(in) :: target_total_charge
859  logical, intent(in) :: enforce_total_charge
860 
861  sll_real64 :: deposited_charge
862  sll_real64 :: charge_correction_factor
863 
865  sll_int32, dimension(:,:,:,:), allocatable :: closest_particle
866  sll_real64, dimension(:,:,:,:), allocatable :: closest_particle_distance
867 
870  type(sll_t_int_list_element_ptr), dimension(:,:,:,:), allocatable :: nodes_in_flow_cell
871  type(sll_t_int_list_element), pointer :: new_int_list_element, head
872 
873  sll_int32 :: k ! particle index
874  sll_int32 :: k_neighbor
875  sll_int32 :: i_x, i_y, i_vx, i_vy
876  sll_int32 :: j_x, j_y, j_vx, j_vy
877  sll_int32 :: m_x, m_y, m_vx, m_vy
878  sll_int32 :: j_tmp
879 
880  sll_int32 :: i_min_x
881  sll_int32 :: i_min_y
882  sll_int32 :: i_min_vx
883  sll_int32 :: i_min_vy
884 
885  sll_int32 :: i_max_x
886  sll_int32 :: i_max_y
887  sll_int32 :: i_max_vx
888  sll_int32 :: i_max_vy
889 
890  sll_int32 :: node_index
891  sll_int32 :: k_particle_closest_to_first_corner
892 
894  type(sll_t_cartesian_mesh_4d),pointer :: g
895 
896  sll_int32 :: g_num_points_x
897  sll_int32 :: g_num_points_y
898  sll_int32 :: g_num_points_vx
899  sll_int32 :: g_num_points_vy
900 
901  sll_real64 :: g_grid_x_min
902  sll_real64 :: g_grid_y_min
903  sll_real64 :: g_grid_vx_min
904  sll_real64 :: g_grid_vy_min
905 
906  sll_real64 :: h_g_grid_x
907  sll_real64 :: h_g_grid_y
908  sll_real64 :: h_g_grid_vx
909  sll_real64 :: h_g_grid_vy
910 
912  sll_real64 :: h_flow_grid_x
913  sll_real64 :: h_flow_grid_y
914  sll_real64 :: h_flow_grid_vx
915  sll_real64 :: h_flow_grid_vy
916 
917  sll_int32 :: flow_grid_num_cells_x
918  sll_int32 :: flow_grid_num_cells_y
919  sll_int32 :: flow_grid_num_cells_vx
920  sll_int32 :: flow_grid_num_cells_vy
921 
922  sll_real64 :: flow_grid_x_min
923  sll_real64 :: flow_grid_y_min
924  sll_real64 :: flow_grid_vx_min
925  sll_real64 :: flow_grid_vy_min
926 
927  sll_int32 :: flow_grid_j_x_min
928  sll_int32 :: flow_grid_j_x_max
929  sll_int32 :: flow_grid_j_y_min
930  sll_int32 :: flow_grid_j_y_max
931  sll_int32 :: flow_grid_j_vx_min
932  sll_int32 :: flow_grid_j_vx_max
933  sll_int32 :: flow_grid_j_vy_min
934  sll_int32 :: flow_grid_j_vy_max
935 
936  sll_real64 :: x
937  sll_real64 :: y
938  sll_real64 :: vx
939  sll_real64 :: vy
940 
941  sll_real64 :: closest_particle_distance_to_first_corner
942  sll_real64 :: particle_distance_to_first_corner
943 
944  sll_real64 :: mesh_period_x
945  sll_real64 :: mesh_period_y
946 
948  sll_real64 :: d11,d12,d13,d14
949  sll_real64 :: d21,d22,d23,d24
950  sll_real64 :: d31,d32,d33,d34
951  sll_real64 :: d41,d42,d43,d44
952 
953  sll_real64, dimension(3) :: coords
954  sll_real64, dimension(4) :: eta
955 
957  sll_real64 :: x_k,y_k,vx_k,vy_k
958  sll_real64 :: x_k_t0,y_k_t0,vx_k_t0,vy_k_t0
959 
960  sll_real64 :: x_to_xk, y_to_yk, vx_to_vxk, vy_to_vyk
961 
962  sll_real64 :: d1_x, d1_y, d1_vx, d1_vy
963  sll_real64 :: d2_x, d2_y, d2_vx, d2_vy
964  sll_real64 :: d3_x, d3_y, d3_vx, d3_vy
965  sll_real64 :: d4_x, d4_y, d4_vx, d4_vy
966 
967  sll_int32 :: nodes_number
968  sll_real64 :: reconstructed_f_value
969 
971  sll_real64 :: x_t0
972  sll_real64 :: y_t0
973  sll_real64 :: vx_t0
974  sll_real64 :: vy_t0
975 
977  sll_real64 :: x_t0_to_xk_t0
978  sll_real64 :: y_t0_to_yk_t0
979  sll_real64 :: vx_t0_to_vxk_t0
980  sll_real64 :: vy_t0_to_vyk_t0
981 
982  sll_int32 :: ierr
983 
985  flow_grid_x_min = self%lbf_grid%eta1_min
986  flow_grid_y_min = self%lbf_grid%eta2_min
987  flow_grid_vx_min = self%lbf_grid%eta3_min
988  flow_grid_vy_min = self%lbf_grid%eta4_min
989 
990  h_flow_grid_x = self%lbf_grid%delta_eta1
991  h_flow_grid_y = self%lbf_grid%delta_eta2
992  h_flow_grid_vx = self%lbf_grid%delta_eta3
993  h_flow_grid_vy = self%lbf_grid%delta_eta4
994 
995  flow_grid_num_cells_x = self%lbf_grid%num_cells1
996  flow_grid_num_cells_y = self%lbf_grid%num_cells2
997  flow_grid_num_cells_vx = self%lbf_grid%num_cells3
998  flow_grid_num_cells_vy = self%lbf_grid%num_cells4
999 
1002  flow_grid_j_x_min = 1
1003  flow_grid_j_x_max = flow_grid_num_cells_x
1004  flow_grid_j_y_min = 1
1005  flow_grid_j_y_max = flow_grid_num_cells_y
1006  flow_grid_j_vx_min = 1
1007  flow_grid_j_vx_max = flow_grid_num_cells_vx
1008  flow_grid_j_vy_min = 1
1009  flow_grid_j_vy_max = flow_grid_num_cells_vy
1010 
1011  deposited_charge = 0.0_f64
1012 
1014  if( reconstruction_set_type == sll_p_lbf_remapping_grid )then
1015 
1017  allocate(nodes_in_flow_cell(flow_grid_num_cells_x, &
1018  flow_grid_num_cells_y, &
1019  flow_grid_num_cells_vx, &
1020  flow_grid_num_cells_vy) &
1021  , stat=ierr)
1022  call sll_s_test_error_code(ierr, 'Memory allocation Failure.', __file__, __line__)
1023 
1024  do j_x = 1, flow_grid_num_cells_x
1025  do j_y = 1, flow_grid_num_cells_y
1026  do j_vx = 1, flow_grid_num_cells_vx
1027  do j_vy = 1, flow_grid_num_cells_vy
1028  nullify(nodes_in_flow_cell(j_x,j_y,j_vx,j_vy)%pointed_element)
1029  end do
1030  end do
1031  end do
1032  end do
1033 
1034  nodes_number = self%sparse_grid_interpolator%size_basis
1035  sll_assert( size(self%tmp_f_values_on_remapping_sparse_grid,1) == nodes_number )
1036  self%tmp_f_values_on_remapping_sparse_grid = 0.0_f64
1037 
1039  do node_index = 1, nodes_number
1040 
1042  x = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(1)
1043  y = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(2)
1044  vx = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(3)
1045  vy = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(4)
1046 
1048  call get_1d_cell_containing_point(j_x, x, flow_grid_x_min, h_flow_grid_x )
1049  call get_1d_cell_containing_point(j_y, y, flow_grid_y_min, h_flow_grid_y )
1050  call get_1d_cell_containing_point(j_vx, vx, flow_grid_vx_min, h_flow_grid_vx)
1051  call get_1d_cell_containing_point(j_vy, vy, flow_grid_vy_min, h_flow_grid_vy)
1052 
1054  if( j_x >= flow_grid_j_x_min .and. j_x <= flow_grid_j_x_max .and. &
1055  j_y >= flow_grid_j_y_min .and. j_y <= flow_grid_j_y_max .and. &
1056  j_vx >= flow_grid_j_vx_min .and. j_vx <= flow_grid_j_vx_max .and. &
1057  j_vy >= flow_grid_j_vy_min .and. j_vy <= flow_grid_j_vy_max )then
1058 
1060  sll_allocate( new_int_list_element, ierr )
1061  new_int_list_element%value = node_index
1062  head => nodes_in_flow_cell(j_x,j_y,j_vx,j_vy)%pointed_element
1063  nodes_in_flow_cell(j_x,j_y,j_vx,j_vy)%pointed_element => sll_f_add_element_in_int_list(head, new_int_list_element)
1064  end if
1065 
1066  end do
1067 
1068  nullify(g)
1069  g_num_points_x = 0
1070  g_num_points_y = 0
1071  g_num_points_vx = 0
1072  g_num_points_vy = 0
1073 
1074  else
1075  sll_assert( reconstruction_set_type == sll_p_lbf_given_grid )
1076 
1078  g => given_grid_4d
1079 
1081  call get_1d_cell_containing_point( j_tmp, g%eta1_min, flow_grid_x_min, h_flow_grid_x )
1082  flow_grid_j_x_min = max(flow_grid_j_x_min, j_tmp)
1083  call get_1d_cell_containing_point( j_tmp, g%eta1_max, flow_grid_x_min, h_flow_grid_x )
1084  flow_grid_j_x_max = min(flow_grid_j_x_max, j_tmp)
1085 
1086  call get_1d_cell_containing_point( j_tmp, g%eta2_min, flow_grid_y_min, h_flow_grid_y )
1087  flow_grid_j_y_min = max(flow_grid_j_y_min, j_tmp)
1088  call get_1d_cell_containing_point( j_tmp, g%eta2_max, flow_grid_y_min, h_flow_grid_y )
1089  flow_grid_j_y_max = min(flow_grid_j_y_max, j_tmp)
1090 
1091  call get_1d_cell_containing_point( j_tmp, g%eta3_min, flow_grid_vx_min, h_flow_grid_vx )
1092  flow_grid_j_vx_min = max(flow_grid_j_vx_min, j_tmp)
1093  call get_1d_cell_containing_point( j_tmp, g%eta3_max, flow_grid_vx_min, h_flow_grid_vx )
1094  flow_grid_j_vx_max = min(flow_grid_j_vx_max, j_tmp)
1095 
1096  call get_1d_cell_containing_point( j_tmp, g%eta4_min, flow_grid_vy_min, h_flow_grid_vy )
1097  flow_grid_j_vy_min = max(flow_grid_j_vy_min, j_tmp)
1098  call get_1d_cell_containing_point( j_tmp, g%eta4_max, flow_grid_vy_min, h_flow_grid_vy )
1099  flow_grid_j_vy_max = min(flow_grid_j_vy_max, j_tmp)
1100 
1101  g_num_points_x = g%num_cells1
1102  if( reconstruct_f_on_last_node(1) ) g_num_points_x = g_num_points_x + 1
1103  sll_assert( size(given_array_4d, 1) == g_num_points_x )
1104 
1105  g_num_points_y = g%num_cells2
1106  if( reconstruct_f_on_last_node(2) ) g_num_points_y = g_num_points_y + 1
1107  sll_assert( size(given_array_4d, 2) == g_num_points_y )
1108 
1109  g_num_points_vx = g%num_cells3
1110  if( reconstruct_f_on_last_node(3) ) g_num_points_vx = g_num_points_vx + 1
1111  sll_assert( size(given_array_4d, 3) == g_num_points_vx )
1112 
1113  g_num_points_vy = g%num_cells4
1114  if( reconstruct_f_on_last_node(4) ) g_num_points_vy = g_num_points_vy + 1
1115  sll_assert( size(given_array_4d, 4) == g_num_points_vy )
1116  given_array_4d(:,:,:,:) = 0.0_f64
1117 
1118  h_g_grid_x = g%delta_eta1
1119  h_g_grid_y = g%delta_eta2
1120  h_g_grid_vx = g%delta_eta3
1121  h_g_grid_vy = g%delta_eta4
1122 
1123  g_grid_x_min = g%eta1_min
1124  g_grid_y_min = g%eta2_min
1125  g_grid_vx_min = g%eta3_min
1126  g_grid_vy_min = g%eta4_min
1127 
1128  end if
1129 
1134  sll_allocate(closest_particle(flow_grid_num_cells_x,flow_grid_num_cells_y,flow_grid_num_cells_vx,flow_grid_num_cells_vy),ierr)
1135  closest_particle(:,:,:,:) = 0
1136 
1137  allocate(closest_particle_distance(flow_grid_num_cells_x, &
1138  flow_grid_num_cells_y, &
1139  flow_grid_num_cells_vx, &
1140  flow_grid_num_cells_vy) , stat=ierr)
1141  call sll_s_test_error_code(ierr, 'Memory allocation Failure.', __file__, __line__)
1142  closest_particle_distance(:,:,:,:) = 0.0_f64
1143 
1144  closest_particle_distance_to_first_corner = 1d30
1145  k_particle_closest_to_first_corner = 0
1146 
1147  do k=1, self%n_particles
1148 
1150  coords = self%get_x(k)
1151  x_k = coords(1)
1152  y_k = coords(2)
1153  coords = self%get_v(k)
1154  vx_k = coords(1)
1155  vy_k = coords(2)
1156 
1158  call get_1d_cell_containing_point(j_x, x_k, flow_grid_x_min, h_flow_grid_x )
1159  call get_1d_cell_containing_point(j_y, y_k, flow_grid_y_min, h_flow_grid_y )
1160  call get_1d_cell_containing_point(j_vx, vx_k, flow_grid_vx_min, h_flow_grid_vx)
1161  call get_1d_cell_containing_point(j_vy, vy_k, flow_grid_vy_min, h_flow_grid_vy)
1162 
1164  if( j_x >= flow_grid_j_x_min .and. j_x <= flow_grid_j_x_max .and. &
1165  j_y >= flow_grid_j_y_min .and. j_y <= flow_grid_j_y_max .and. &
1166  j_vx >= flow_grid_j_vx_min .and. j_vx <= flow_grid_j_vx_max .and. &
1167  j_vy >= flow_grid_j_vy_min .and. j_vy <= flow_grid_j_vy_max )then
1168 
1169  call update_closest_particle_arrays( &
1170  k, &
1171  x_k - flow_grid_x_min, &
1172  y_k - flow_grid_y_min, &
1173  vx_k - flow_grid_vx_min, &
1174  vy_k - flow_grid_vy_min, &
1175  j_x, j_y, j_vx, j_vy, &
1176  h_flow_grid_x, &
1177  h_flow_grid_y, &
1178  h_flow_grid_vx, &
1179  h_flow_grid_vy, &
1180  closest_particle, &
1181  closest_particle_distance &
1182  )
1183  end if
1184 
1185  particle_distance_to_first_corner = abs(x_k - flow_grid_x_min ) &
1186  + abs(y_k - flow_grid_y_min ) &
1187  + abs(vx_k - flow_grid_vx_min) &
1188  + abs(vy_k - flow_grid_vy_min)
1189  if( particle_distance_to_first_corner < closest_particle_distance_to_first_corner )then
1190  closest_particle_distance_to_first_corner = particle_distance_to_first_corner
1191  k_particle_closest_to_first_corner = k
1192  end if
1193  end do
1194 
1195  closest_particle(1,1,1,1) = k_particle_closest_to_first_corner
1196 
1197  if( .not. ( self%domain_is_periodic(1) .and. self%domain_is_periodic(2) ) )then
1198  print*, "WARNING -- STOP -- verify that the non-periodic case is well implemented"
1199  stop
1200  end if
1201 
1208 
1209  if(self%domain_is_periodic(1)) then
1211  mesh_period_x = self%lbf_grid%eta1_max - self%lbf_grid%eta1_min
1212  else
1213  mesh_period_x = 0.0_f64
1214  end if
1215 
1216  if(self%domain_is_periodic(2)) then
1218  mesh_period_y = self%lbf_grid%eta2_max - self%lbf_grid%eta2_min
1219  else
1220  mesh_period_y = 0.0_f64
1221  end if
1222 
1223  do j_x = flow_grid_j_x_min, flow_grid_j_x_max
1224  do j_y = flow_grid_j_y_min, flow_grid_j_y_max
1225  do j_vx = flow_grid_j_vx_min, flow_grid_j_vx_max
1226  do j_vy = flow_grid_j_vy_min, flow_grid_j_vy_max
1227 
1229  k = closest_particle(j_x,j_y,j_vx,j_vy)
1230 
1231  if(k == 0) then
1232  if( j_x > 1 )then
1233  update_closest_particle_arrays_using_neighbor_cells(-1,0,0,0)
1234  end if
1235  if( j_x < flow_grid_num_cells_x )then
1236  update_closest_particle_arrays_using_neighbor_cells( 1,0,0,0)
1237  end if
1238 
1239  if( j_y > 1 )then
1240  update_closest_particle_arrays_using_neighbor_cells(0,-1,0,0)
1241  end if
1242  if( j_y < flow_grid_num_cells_y )then
1243  update_closest_particle_arrays_using_neighbor_cells(0, 1,0,0)
1244  end if
1245 
1246  if( j_vx > 1 )then
1247  update_closest_particle_arrays_using_neighbor_cells(0,0,-1,0)
1248  end if
1249  if( j_vx < flow_grid_num_cells_vx )then
1250  update_closest_particle_arrays_using_neighbor_cells(0,0, 1,0)
1251  end if
1252 
1253  if( j_vy > 1 )then
1254  update_closest_particle_arrays_using_neighbor_cells(0,0,0,-1)
1255  end if
1256  if( j_vy < flow_grid_num_cells_vy )then
1257  update_closest_particle_arrays_using_neighbor_cells(0,0,0, 1)
1258  end if
1259  end if
1260 
1261  k = closest_particle(j_x,j_y,j_vx,j_vy)
1262  sll_assert(k /= 0)
1263 
1265 
1267 
1268  call self%get_ltp_deformation_matrix ( &
1269  k, &
1270  mesh_period_x, &
1271  mesh_period_y, &
1272  h_flow_grid_x, &
1273  h_flow_grid_y, &
1274  h_flow_grid_vx, &
1275  h_flow_grid_vy, &
1276  x_k,y_k,vx_k,vy_k, &
1277  d11,d12,d13,d14, &
1278  d21,d22,d23,d24, &
1279  d31,d32,d33,d34, &
1280  d41,d42,d43,d44 &
1281  )
1282 
1284  call convert_1d_index_to_4d( &
1285  m_x,m_y,m_vx,m_vy, &
1286  k, &
1287  self%n_particles_x, self%n_particles_y, &
1288  self%n_particles_vx, self%n_particles_vy &
1289  )
1291  x_k_t0 = flow_grid_x_min + (m_x-0.5) * h_flow_grid_x
1292  y_k_t0 = flow_grid_y_min + (m_y-0.5) * h_flow_grid_y
1293  vx_k_t0 = flow_grid_vx_min + (m_vx-0.5) * h_flow_grid_vx
1294  vy_k_t0 = flow_grid_vy_min + (m_vy-0.5) * h_flow_grid_vy
1295 
1296 
1300 
1301  if( reconstruction_set_type == sll_p_lbf_remapping_grid )then
1302 
1303  new_int_list_element => nodes_in_flow_cell(j_x,j_y,j_vx,j_vy)%pointed_element
1304 
1305  do while( associated(new_int_list_element) )
1306  node_index = new_int_list_element%value
1307 
1308  x = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(1)
1309  y = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(2)
1310  vx = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(3)
1311  vy = self%sparse_grid_interpolator%hierarchy(node_index)%coordinate(4)
1312 
1314  x_to_xk = x - x_k
1315  y_to_yk = y - y_k
1316  vx_to_vxk = vx - vx_k
1317  vy_to_vyk = vy - vy_k
1318 
1320  x_t0_to_xk_t0 = d11 * x_to_xk + d12 * y_to_yk + d13 * vx_to_vxk + d14 * vy_to_vyk
1321  y_t0_to_yk_t0 = d21 * x_to_xk + d22 * y_to_yk + d23 * vx_to_vxk + d24 * vy_to_vyk
1322  vx_t0_to_vxk_t0 = d31 * x_to_xk + d32 * y_to_yk + d33 * vx_to_vxk + d34 * vy_to_vyk
1323  vy_t0_to_vyk_t0 = d41 * x_to_xk + d42 * y_to_yk + d43 * vx_to_vxk + d44 * vy_to_vyk
1324 
1325  x_t0 = x_t0_to_xk_t0 + x_k_t0
1326  y_t0 = y_t0_to_yk_t0 + y_k_t0
1327  vx_t0 = vx_t0_to_vxk_t0 + vx_k_t0
1328  vy_t0 = vy_t0_to_vyk_t0 + vy_k_t0
1329 
1331  call periodic_correction(self,x_t0,y_t0)
1332 
1334  eta(1) = x_t0
1335  eta(2) = y_t0
1336  eta(3) = vx_t0
1337  eta(4) = vy_t0
1338 
1340  reconstructed_f_value = self%interpolate_value_of_remapped_f(eta)
1341  self%tmp_f_values_on_remapping_sparse_grid(node_index) = reconstructed_f_value
1342 
1343  new_int_list_element => new_int_list_element%next
1344 
1345  end do
1346 
1349 
1350  else
1351  sll_assert( reconstruction_set_type == sll_p_lbf_given_grid )
1352 
1361 
1362  i_min_x = int(ceiling( (flow_grid_x_min - g_grid_x_min + (j_x-1) * h_flow_grid_x) / h_g_grid_x + 1 ))
1363  i_min_y = int(ceiling( (flow_grid_y_min - g_grid_y_min + (j_y-1) * h_flow_grid_y) / h_g_grid_y + 1 ))
1364  i_min_vx = int(ceiling( (flow_grid_vx_min - g_grid_vx_min + (j_vx-1) * h_flow_grid_vx)/ h_g_grid_vx + 1 ))
1365  i_min_vy = int(ceiling( (flow_grid_vy_min - g_grid_vy_min + (j_vy-1) * h_flow_grid_vy)/ h_g_grid_vy + 1 ))
1366 
1367  i_max_x = int(ceiling( (flow_grid_x_min - g_grid_x_min + (j_x) * h_flow_grid_x) / h_g_grid_x + 1 )) - 1
1368  i_max_y = int(ceiling( (flow_grid_y_min - g_grid_y_min + (j_y) * h_flow_grid_y) / h_g_grid_y + 1 )) - 1
1369  i_max_vx = int(ceiling( (flow_grid_vx_min - g_grid_vx_min + (j_vx) * h_flow_grid_vx)/ h_g_grid_vx + 1 )) - 1
1370  i_max_vy = int(ceiling( (flow_grid_vy_min - g_grid_vy_min + (j_vy) * h_flow_grid_vy)/ h_g_grid_vy + 1 )) - 1
1371 
1372  i_min_x = max(i_min_x, 1)
1373  i_min_y = max(i_min_y, 1)
1374  i_min_vx = max(i_min_vx, 1)
1375  i_min_vy = max(i_min_vy, 1)
1376 
1377  i_max_x = min(i_max_x, g_num_points_x)
1378  i_max_y = min(i_max_y, g_num_points_y)
1379  i_max_vx = min(i_max_vx, g_num_points_vx)
1380  i_max_vy = min(i_max_vy, g_num_points_vy)
1381 
1382  do i_x = i_min_x, i_max_x
1383  x = g_grid_x_min + (i_x-1)*h_g_grid_x
1384  x_to_xk = x - x_k
1385  d1_x = d11 * x_to_xk
1386  d2_x = d21 * x_to_xk
1387  d3_x = d31 * x_to_xk
1388  d4_x = d41 * x_to_xk
1389 
1390  do i_y = i_min_y, i_max_y
1391  y = g_grid_y_min + (i_y-1)*h_g_grid_y
1392  y_to_yk = y - y_k
1393  d1_y = d12 * y_to_yk
1394  d2_y = d22 * y_to_yk
1395  d3_y = d32 * y_to_yk
1396  d4_y = d42 * y_to_yk
1397 
1398  do i_vx = i_min_vx, i_max_vx
1399  vx = g_grid_vx_min + (i_vx-1)*h_g_grid_vx
1400  vx_to_vxk = vx - vx_k
1401  d1_vx = d13 * vx_to_vxk
1402  d2_vx = d23 * vx_to_vxk
1403  d3_vx = d33 * vx_to_vxk
1404  d4_vx = d43 * vx_to_vxk
1405 
1406 
1407  do i_vy = i_min_vy, i_max_vy
1408  vy = g_grid_vy_min + (i_vy-1)*h_g_grid_vy
1409  vy_to_vyk = vy - vy_k
1410  d1_vy = d14 * vy_to_vyk
1411  d2_vy = d24 * vy_to_vyk
1412  d3_vy = d34 * vy_to_vyk
1413  d4_vy = d44 * vy_to_vyk
1414 
1417 
1418  x_t0_to_xk_t0 = d1_x + d1_y + d1_vx + d1_vy
1419  y_t0_to_yk_t0 = d2_x + d2_y + d2_vx + d2_vy
1420  vx_t0_to_vxk_t0 = d3_x + d3_y + d3_vx + d3_vy
1421  vy_t0_to_vyk_t0 = d4_x + d4_y + d4_vx + d4_vy
1422 
1423  x_t0 = x_t0_to_xk_t0 + x_k_t0
1424  y_t0 = y_t0_to_yk_t0 + y_k_t0
1425  vx_t0 = vx_t0_to_vxk_t0 + vx_k_t0
1426  vy_t0 = vy_t0_to_vyk_t0 + vy_k_t0
1427 
1429  call periodic_correction(self,x_t0,y_t0)
1430 
1432  eta(1) = x_t0
1433  eta(2) = y_t0
1434  eta(3) = vx_t0
1435  eta(4) = vy_t0
1436 
1438  reconstructed_f_value = self%interpolate_value_of_remapped_f(eta)
1439 
1440  ! [DEBUG]
1441  if( .false. )then
1442  if( reconstruction_set_type == sll_p_lbf_given_grid )then
1443  print*, "[WRITE ON GIVEN GRID] reconstructing "
1444  print*, "on: eta = ", eta
1445  print*, "value: reconstructed_f_value = ", reconstructed_f_value
1446  end if
1447  end if
1448 
1450 
1451  if( reconstructed_f_value /= 0 )then
1452 
1453  sll_assert( 1 <= i_x .and. i_x <= g_num_points_x )
1454  sll_assert( 1 <= i_y .and. i_y <= g_num_points_y )
1455  sll_assert( 1 <= i_vx .and. i_vx <= g_num_points_vx )
1456  sll_assert( 1 <= i_vy .and. i_vy <= g_num_points_vy )
1457 
1458  if( abs(given_array_4d(i_x, i_y, i_vx, i_vy)) > 0.00001 )then
1459  print *, "Warning -- 987698666999979979 -- stored value should be zero"
1460  print *, "given_array_4d(i_x, i_y, i_vx, i_vy) = ", given_array_4d(i_x, i_y, i_vx, i_vy)
1461  print *, "with (i_x, i_y, i_vx, i_vy) = ", i_x, i_y, i_vx, i_vy
1462  end if
1463  given_array_4d(i_x, i_y, i_vx, i_vy) = reconstructed_f_value
1464 
1465  else
1466  print *, "Warning -- 654654535466545434564 -- just reconstructed a Zero value"
1467  end if
1469  end do
1470  end do
1471  end do
1472  end do
1473  end if
1475  end do
1476  end do
1477  end do
1478  end do
1479 
1480  if( enforce_total_charge )then
1481  sll_assert( reconstruction_set_type == sll_p_lbf_given_grid )
1482  if( deposited_charge == 0 )then
1483  print *, "WARNING (76576537475) -- total deposited charge is zero, which is strange..."
1484  print *, " -- (no charge correction in this case) "
1485  else
1486  charge_correction_factor = target_total_charge / deposited_charge
1487  print *, "[Enforcing charge]: target_total_charge, deposited_charge, charge_correction_factor = ", &
1488  target_total_charge, deposited_charge, charge_correction_factor
1489  do i_x = 1, g_num_points_x
1490  do i_y = 1, g_num_points_y
1491  do i_vx = 1, g_num_points_vx
1492  do i_vy = 1, g_num_points_vy
1493  given_array_4d(i_x, i_y, i_vx, i_vy) = given_array_4d(i_x, i_y, i_vx, i_vy) &
1494  * charge_correction_factor
1495  end do
1496  end do
1497  end do
1498  end do
1499  end if
1500  end if
1501 
1502  end subroutine reconstruct_f_lbf
1503 
1506  subroutine reconstruct_f_lbf_on_remapping_grid( self )
1507  class(sll_t_particle_group_2d2v_lbf), intent(inout) :: self
1508  type(sll_t_cartesian_mesh_4d), pointer :: void_grid_4d
1509  sll_real64, dimension(:,:,:,:), pointer :: void_array_4d
1510  logical :: dummy_reconstruct_f_on_last_node(4)
1511  sll_real64 :: dummy_target_total_charge
1512  logical :: dummy_enforce_total_charge
1513 
1514  nullify(void_grid_4d)
1515  nullify(void_array_4d)
1516  dummy_reconstruct_f_on_last_node = .false.
1517  dummy_target_total_charge = 0.0_f64
1518  dummy_enforce_total_charge = .false.
1519 
1520  call self%reconstruct_f_lbf( &
1521  sll_p_lbf_remapping_grid, &
1522  void_grid_4d, &
1523  void_array_4d, &
1524  dummy_reconstruct_f_on_last_node, &
1525  dummy_target_total_charge, &
1526  dummy_enforce_total_charge )
1527 
1528  end subroutine reconstruct_f_lbf_on_remapping_grid
1529 
1530 
1533  subroutine reconstruct_f_lbf_on_given_grid( &
1534  self, &
1535  given_grid_4d, &
1536  given_array_4d, &
1537  reconstruct_f_on_last_node, &
1538  target_total_charge, &
1539  enforce_total_charge &
1540  )
1541  class(sll_t_particle_group_2d2v_lbf), intent(inout) :: self
1542  type(sll_t_cartesian_mesh_4d), pointer, intent(in) :: given_grid_4d
1543  sll_real64, dimension(:,:,:,:), pointer, intent(inout) :: given_array_4d
1544  logical, intent(in) :: reconstruct_f_on_last_node(4)
1545  sll_real64, intent(in) :: target_total_charge
1546  logical, intent(in) :: enforce_total_charge
1547 
1548  call self%reconstruct_f_lbf( &
1549  sll_p_lbf_given_grid, &
1550  given_grid_4d, &
1551  given_array_4d, &
1552  reconstruct_f_on_last_node, &
1553  target_total_charge, &
1554  enforce_total_charge )
1555 
1556  end subroutine reconstruct_f_lbf_on_given_grid
1557 
1558 
1585  ! <<get_ltp_deformation_matrix>>
1586  subroutine get_ltp_deformation_matrix ( &
1587  self, &
1588  k, &
1589  mesh_period_x, &
1590  mesh_period_y, &
1591  h_particles_x, &
1592  h_particles_y, &
1593  h_particles_vx, &
1594  h_particles_vy, &
1595  x_k, y_k, &
1596  vx_k, vy_k, &
1597  d11,d12,d13,d14, &
1598  d21,d22,d23,d24, &
1599  d31,d32,d33,d34, &
1600  d41,d42,d43,d44 &
1601  )
1602 
1603  class(sll_t_particle_group_2d2v_lbf),intent(inout) :: self
1604  sll_int32, intent(in) :: k
1605 
1606  sll_real64, intent(in) :: mesh_period_x
1607  sll_real64, intent(in) :: mesh_period_y
1608 
1609  sll_real64, intent(in) :: h_particles_x
1610  sll_real64, intent(in) :: h_particles_y
1611  sll_real64, intent(in) :: h_particles_vx
1612  sll_real64, intent(in) :: h_particles_vy
1613 
1614  sll_real64, intent(out) :: x_k, y_k
1615  sll_real64, intent(out) :: vx_k, vy_k
1616  sll_real64, intent(out) :: d11,d12,d13,d14
1617  sll_real64, intent(out) :: d21,d22,d23,d24
1618  sll_real64, intent(out) :: d31,d32,d33,d34
1619  sll_real64, intent(out) :: d41,d42,d43,d44
1620 
1621  sll_int32 :: k_ngb
1622  sll_real64 :: x_k_left, x_k_right
1623  sll_real64 :: y_k_left, y_k_right
1624  sll_real64 :: vx_k_left, vx_k_right
1625  sll_real64 :: vy_k_left, vy_k_right
1626  sll_real64, dimension(3) :: coords
1627 
1628  sll_real64 :: j11,j12,j13,j14
1629  sll_real64 :: j21,j22,j23,j24
1630  sll_real64 :: j31,j32,j33,j34
1631  sll_real64 :: j41,j42,j43,j44
1632  sll_real64 :: factor, det_j, inv_det_j
1633 
1634  ! for clarity:
1635  sll_int32, parameter :: dim_x = 1
1636  sll_int32, parameter :: dim_y = 2
1637  sll_int32, parameter :: dim_vx = 3
1638  sll_int32, parameter :: dim_vy = 4
1639  sll_int32, parameter :: left_ngb = -1
1640  sll_int32, parameter :: right_ngb = 1
1641 
1642  logical domain_is_x_periodic
1643  logical domain_is_y_periodic
1644 
1645  domain_is_x_periodic = self%domain_is_periodic(1)
1646  domain_is_y_periodic = self%domain_is_periodic(2)
1647 
1648  coords(:) = 0.0_f64
1649 
1650  coords = self%get_x(k)
1651  x_k = coords(1)
1652  y_k = coords(2)
1653 
1654  coords = self%get_v(k)
1655  vx_k = coords(1)
1656  vy_k = coords(2)
1657 
1659 
1661  factor = 1./(2*h_particles_x)
1662 
1663  k_ngb = self%get_neighbor_index(k, dim_x, right_ngb)
1664  if( k_ngb == k )then
1666  factor = 2*factor
1667  x_k_right = x_k
1668  y_k_right = y_k
1669  vx_k_right = vx_k
1670  vy_k_right = vy_k
1671  else
1672  coords = self%get_x(k_ngb)
1673  x_k_right = coords(1)
1674  y_k_right = coords(2)
1675  coords = self%get_v(k_ngb)
1676  vx_k_right = coords(1)
1677  vy_k_right = coords(2)
1678 
1679  if( domain_is_x_periodic .and. x_k_right < x_k - 0.5*mesh_period_x ) x_k_right = x_k_right + mesh_period_x
1680  if( domain_is_x_periodic .and. x_k_right > x_k + 0.5*mesh_period_x ) x_k_right = x_k_right - mesh_period_x
1681  if( domain_is_y_periodic .and. y_k_right < y_k - 0.5*mesh_period_y ) y_k_right = y_k_right + mesh_period_y
1682  if( domain_is_y_periodic .and. y_k_right > y_k + 0.5*mesh_period_y ) y_k_right = y_k_right - mesh_period_y
1683  end if
1684 
1685 
1686  k_ngb = self%get_neighbor_index(k, dim_x, left_ngb)
1687  if( k_ngb == k )then
1689  factor = 2*factor
1690  x_k_left = x_k
1691  y_k_left = y_k
1692  vx_k_left = vx_k
1693  vy_k_left = vy_k
1694  else
1695  coords = self%get_x(k_ngb)
1696  x_k_left = coords(1)
1697  y_k_left = coords(2)
1698  coords = self%get_v(k_ngb)
1699  vx_k_left = coords(1)
1700  vy_k_left = coords(2)
1701 
1702  if( domain_is_x_periodic .and. x_k_left < x_k - 0.5*mesh_period_x ) x_k_left = x_k_left + mesh_period_x
1703  if( domain_is_x_periodic .and. x_k_left > x_k + 0.5*mesh_period_x ) x_k_left = x_k_left - mesh_period_x
1704  if( domain_is_y_periodic .and. y_k_left < y_k - 0.5*mesh_period_y ) y_k_left = y_k_left + mesh_period_y
1705  if( domain_is_y_periodic .and. y_k_left > y_k + 0.5*mesh_period_y ) y_k_left = y_k_left - mesh_period_y
1706  end if
1707 
1708  j11 = factor * ( x_k_right - x_k_left )
1709  j21 = factor * ( y_k_right - y_k_left )
1710  j31 = factor * ( vx_k_right - vx_k_left )
1711  j41 = factor * ( vy_k_right - vy_k_left )
1712 
1714  factor = 1./(2*h_particles_y)
1715 
1716  k_ngb = self%get_neighbor_index(k, dim_y, right_ngb)
1717  if( k_ngb == k )then
1719  factor = 2*factor
1720  x_k_right = x_k
1721  y_k_right = y_k
1722  vx_k_right = vx_k
1723  vy_k_right = vy_k
1724  else
1725  coords = self%get_x(k_ngb)
1726  x_k_right = coords(1)
1727  y_k_right = coords(2)
1728  coords = self%get_v(k_ngb)
1729  vx_k_right = coords(1)
1730  vy_k_right = coords(2)
1731 
1732  if( domain_is_x_periodic .and. x_k_right < x_k - 0.5*mesh_period_x ) x_k_right = x_k_right + mesh_period_x
1733  if( domain_is_x_periodic .and. x_k_right > x_k + 0.5*mesh_period_x ) x_k_right = x_k_right - mesh_period_x
1734  if( domain_is_y_periodic .and. y_k_right < y_k - 0.5*mesh_period_y ) y_k_right = y_k_right + mesh_period_y
1735  if( domain_is_y_periodic .and. y_k_right > y_k + 0.5*mesh_period_y ) y_k_right = y_k_right - mesh_period_y
1736  end if
1737 
1738  k_ngb = self%get_neighbor_index(k, dim_y, left_ngb)
1739  if( k_ngb == k )then
1741  factor = 2*factor
1742  x_k_left = x_k
1743  y_k_left = y_k
1744  vx_k_left = vx_k
1745  vy_k_left = vy_k
1746  else
1747  coords = self%get_x(k_ngb)
1748  x_k_left = coords(1)
1749  y_k_left = coords(2)
1750  coords = self%get_v(k_ngb)
1751  vx_k_left = coords(1)
1752  vy_k_left = coords(2)
1753 
1754  if( domain_is_x_periodic .and. x_k_left < x_k - 0.5*mesh_period_x ) x_k_left = x_k_left + mesh_period_x
1755  if( domain_is_x_periodic .and. x_k_left > x_k + 0.5*mesh_period_x ) x_k_left = x_k_left - mesh_period_x
1756  if( domain_is_y_periodic .and. y_k_left < y_k - 0.5*mesh_period_y ) y_k_left = y_k_left + mesh_period_y
1757  if( domain_is_y_periodic .and. y_k_left > y_k + 0.5*mesh_period_y ) y_k_left = y_k_left - mesh_period_y
1758  end if
1759 
1760  j12 = factor * ( x_k_right - x_k_left )
1761  j22 = factor * ( y_k_right - y_k_left )
1762  j32 = factor * ( vx_k_right - vx_k_left )
1763  j42 = factor * ( vy_k_right - vy_k_left )
1764 
1765 
1767  factor = 1./(2*h_particles_vx)
1768 
1769  k_ngb = self%get_neighbor_index(k, dim_vx, right_ngb)
1770  if( k_ngb == k )then
1772  factor = 2*factor
1773  x_k_right = x_k
1774  y_k_right = y_k
1775  vx_k_right = vx_k
1776  vy_k_right = vy_k
1777  else
1778  coords = self%get_x(k_ngb)
1779  x_k_right = coords(1)
1780  y_k_right = coords(2)
1781  coords = self%get_v(k_ngb)
1782  vx_k_right = coords(1)
1783  vy_k_right = coords(2)
1784 
1785  if( domain_is_x_periodic .and. x_k_right < x_k - 0.5*mesh_period_x ) x_k_right = x_k_right + mesh_period_x
1786  if( domain_is_x_periodic .and. x_k_right > x_k + 0.5*mesh_period_x ) x_k_right = x_k_right - mesh_period_x
1787  if( domain_is_y_periodic .and. y_k_right < y_k - 0.5*mesh_period_y ) y_k_right = y_k_right + mesh_period_y
1788  if( domain_is_y_periodic .and. y_k_right > y_k + 0.5*mesh_period_y ) y_k_right = y_k_right - mesh_period_y
1789  end if
1790 
1791  k_ngb = self%get_neighbor_index(k, dim_vx, left_ngb)
1792  if( k_ngb == k )then
1794  factor = 2*factor
1795  x_k_left = x_k
1796  y_k_left = y_k
1797  vx_k_left = vx_k
1798  vy_k_left = vy_k
1799  else
1800  coords = self%get_x(k_ngb)
1801  x_k_left = coords(1)
1802  y_k_left = coords(2)
1803  coords = self%get_v(k_ngb)
1804  vx_k_left = coords(1)
1805  vy_k_left = coords(2)
1806 
1807  if( domain_is_x_periodic .and. x_k_left < x_k - 0.5*mesh_period_x ) x_k_left = x_k_left + mesh_period_x
1808  if( domain_is_x_periodic .and. x_k_left > x_k + 0.5*mesh_period_x ) x_k_left = x_k_left - mesh_period_x
1809  if( domain_is_y_periodic .and. y_k_left < y_k - 0.5*mesh_period_y ) y_k_left = y_k_left + mesh_period_y
1810  if( domain_is_y_periodic .and. y_k_left > y_k + 0.5*mesh_period_y ) y_k_left = y_k_left - mesh_period_y
1811  end if
1812 
1813  j13 = factor * ( x_k_right - x_k_left )
1814  j23 = factor * ( y_k_right - y_k_left )
1815  j33 = factor * ( vx_k_right - vx_k_left )
1816  j43 = factor * ( vy_k_right - vy_k_left )
1817 
1818 
1820  factor = 1./(2*h_particles_vy)
1821 
1822  k_ngb = self%get_neighbor_index(k, dim_vy, right_ngb)
1823  if( k_ngb == k )then
1825  factor = 2*factor
1826  x_k_right = x_k
1827  y_k_right = y_k
1828  vx_k_right = vx_k
1829  vy_k_right = vy_k
1830  else
1831  coords = self%get_x(k_ngb)
1832  x_k_right = coords(1)
1833  y_k_right = coords(2)
1834  coords = self%get_v(k_ngb)
1835  vx_k_right = coords(1)
1836  vy_k_right = coords(2)
1837 
1838  if( domain_is_x_periodic .and. x_k_right < x_k - 0.5*mesh_period_x ) x_k_right = x_k_right + mesh_period_x
1839  if( domain_is_x_periodic .and. x_k_right > x_k + 0.5*mesh_period_x ) x_k_right = x_k_right - mesh_period_x
1840  if( domain_is_y_periodic .and. y_k_right < y_k - 0.5*mesh_period_y ) y_k_right = y_k_right + mesh_period_y
1841  if( domain_is_y_periodic .and. y_k_right > y_k + 0.5*mesh_period_y ) y_k_right = y_k_right - mesh_period_y
1842  end if
1843 
1844  k_ngb = self%get_neighbor_index(k, dim_vy, left_ngb)
1845  if( k_ngb == k )then
1847  factor = 2*factor
1848  x_k_left = x_k
1849  y_k_left = y_k
1850  vx_k_left = vx_k
1851  vy_k_left = vy_k
1852  else
1853  coords = self%get_x(k_ngb)
1854  x_k_left = coords(1)
1855  y_k_left = coords(2)
1856  coords = self%get_v(k_ngb)
1857  vx_k_left = coords(1)
1858  vy_k_left = coords(2)
1859 
1860  if( domain_is_x_periodic .and. x_k_left < x_k - 0.5*mesh_period_x ) x_k_left = x_k_left + mesh_period_x
1861  if( domain_is_x_periodic .and. x_k_left > x_k + 0.5*mesh_period_x ) x_k_left = x_k_left - mesh_period_x
1862  if( domain_is_y_periodic .and. y_k_left < y_k - 0.5*mesh_period_y ) y_k_left = y_k_left + mesh_period_y
1863  if( domain_is_y_periodic .and. y_k_left > y_k + 0.5*mesh_period_y ) y_k_left = y_k_left - mesh_period_y
1864  end if
1865 
1866  j14 = factor * ( x_k_right - x_k_left )
1867  j24 = factor * ( y_k_right - y_k_left )
1868  j34 = factor * ( vx_k_right - vx_k_left )
1869  j44 = factor * ( vy_k_right - vy_k_left )
1870 
1871 
1873 
1874  det_j = j11 * j22 * j33 * j44 &
1875  + j11 * j23 * j34 * j42 &
1876  + j11 * j24 * j32 * j43 &
1877  + j12 * j21 * j34 * j43 &
1878  + j12 * j23 * j31 * j44 &
1879  + j12 * j24 * j33 * j41 &
1880  + j13 * j21 * j32 * j44 &
1881  + j13 * j22 * j34 * j41 &
1882  + j13 * j24 * j31 * j42 &
1883  + j14 * j21 * j33 * j42 &
1884  + j14 * j22 * j31 * j43 &
1885  + j14 * j23 * j32 * j41 &
1886  - j11 * j22 * j34 * j43 &
1887  - j11 * j23 * j32 * j44 &
1888  - j11 * j24 * j33 * j42 &
1889  - j12 * j21 * j33 * j44 &
1890  - j12 * j23 * j34 * j41 &
1891  - j12 * j24 * j31 * j43 &
1892  - j13 * j21 * j34 * j42 &
1893  - j13 * j22 * j31 * j44 &
1894  - j13 * j24 * j32 * j41 &
1895  - j14 * j21 * j32 * j43 &
1896  - j14 * j22 * j33 * j41 &
1897  - j14 * j23 * j31 * j42
1898 
1899  inv_det_j = 1./det_j
1900 
1901  d11 = j22 * j33 * j44 &
1902  + j23 * j34 * j42 &
1903  + j24 * j32 * j43 &
1904  - j22 * j34 * j43 &
1905  - j23 * j32 * j44 &
1906  - j24 * j33 * j42
1907  d11 = inv_det_j * d11
1908 
1909  d12 = j12 * j34 * j43 &
1910  + j13 * j32 * j44 &
1911  + j14 * j33 * j42 &
1912  - j12 * j33 * j44 &
1913  - j13 * j34 * j42 &
1914  - j14 * j32 * j43
1915  d12 = inv_det_j * d12
1916 
1917  d13 = j12 * j23 * j44 &
1918  + j13 * j24 * j42 &
1919  + j14 * j22 * j43 &
1920  - j12 * j24 * j43 &
1921  - j13 * j22 * j44 &
1922  - j14 * j23 * j42
1923  d13 = inv_det_j * d13
1924 
1925  d14 = j12 * j24 * j33 &
1926  + j13 * j22 * j34 &
1927  + j14 * j23 * j32 &
1928  - j12 * j23 * j34 &
1929  - j13 * j24 * j32 &
1930  - j14 * j22 * j33
1931  d14 = inv_det_j * d14
1932 
1933  d21 = j21 * j34 * j43 &
1934  + j23 * j31 * j44 &
1935  + j24 * j33 * j41 &
1936  - j21 * j33 * j44 &
1937  - j23 * j34 * j41 &
1938  - j24 * j31 * j43
1939  d21 = inv_det_j * d21
1940 
1941  d22 = j11 * j33 * j44 &
1942  + j13 * j34 * j41 &
1943  + j14 * j31 * j43 &
1944  - j11 * j34 * j43 &
1945  - j13 * j31 * j44 &
1946  - j14 * j33 * j41
1947  d22 = inv_det_j * d22
1948 
1949  d23 = j11 * j24 * j43 &
1950  + j13 * j21 * j44 &
1951  + j14 * j23 * j41 &
1952  - j11 * j23 * j44 &
1953  - j13 * j24 * j41 &
1954  - j14 * j21 * j43
1955  d23 = inv_det_j * d23
1956 
1957  d24 = j11 * j23 * j34 &
1958  + j13 * j24 * j31 &
1959  + j14 * j21 * j33 &
1960  - j11 * j24 * j33 &
1961  - j13 * j21 * j34 &
1962  - j14 * j23 * j31
1963  d24 = inv_det_j * d24
1964 
1965  d31 = j21 * j32 * j44 &
1966  + j22 * j34 * j41 &
1967  + j24 * j31 * j42 &
1968  - j21 * j34 * j42 &
1969  - j22 * j31 * j44 &
1970  - j24 * j32 * j41
1971  d31 = inv_det_j * d31
1972 
1973  d32 = j11 * j34 * j42 &
1974  + j12 * j31 * j44 &
1975  + j14 * j32 * j41 &
1976  - j11 * j32 * j44 &
1977  - j12 * j34 * j41 &
1978  - j14 * j31 * j42
1979  d32 = inv_det_j * d32
1980 
1981  d33 = j11 * j22 * j44 &
1982  + j12 * j24 * j41 &
1983  + j14 * j21 * j42 &
1984  - j11 * j24 * j42 &
1985  - j12 * j21 * j44 &
1986  - j14 * j22 * j41
1987  d33 = inv_det_j * d33
1988 
1989  d34 = j11 * j24 * j32 &
1990  + j12 * j21 * j34 &
1991  + j14 * j22 * j31 &
1992  - j11 * j22 * j34 &
1993  - j12 * j24 * j31 &
1994  - j14 * j21 * j32
1995  d34 = inv_det_j * d34
1996 
1997  d41 = j21 * j33 * j42 &
1998  + j22 * j31 * j43 &
1999  + j23 * j32 * j41 &
2000  - j21 * j32 * j43 &
2001  - j22 * j33 * j41 &
2002  - j23 * j31 * j42
2003  d41 = inv_det_j * d41
2004 
2005  d42 = j11 * j32 * j43 &
2006  + j12 * j33 * j41 &
2007  + j13 * j31 * j42 &
2008  - j11 * j33 * j42 &
2009  - j12 * j31 * j43 &
2010  - j13 * j32 * j41
2011  d42 = inv_det_j * d42
2012 
2013  d43 = j11 * j23 * j42 &
2014  + j12 * j21 * j43 &
2015  + j13 * j22 * j41 &
2016  - j11 * j22 * j43 &
2017  - j12 * j23 * j41 &
2018  - j13 * j21 * j42
2019  d43 = inv_det_j * d43
2020 
2021  d44 = j11 * j22 * j33 &
2022  + j12 * j23 * j31 &
2023  + j13 * j21 * j32 &
2024  - j11 * j23 * j32 &
2025  - j12 * j21 * j33 &
2026  - j13 * j22 * j31
2027  d44 = inv_det_j * d44
2028 
2029  end subroutine get_ltp_deformation_matrix
2030 
2033  subroutine periodic_correction(p_group, x, y)
2034  class(sll_t_particle_group_2d2v_lbf), intent(in) :: p_group
2035  sll_real64, intent(inout) :: x
2036  sll_real64, intent(inout) :: y
2037  sll_real64 :: eta_max, eta_min
2038 
2039  if( p_group%domain_is_periodic(1) )then
2040  eta_min = p_group%lbf_grid%eta1_min
2041  eta_max = p_group%lbf_grid%eta1_max
2042  if( (x < eta_min) .or. (x >= eta_max) ) x = eta_min + modulo(x - eta_min, eta_max - eta_min)
2043  end if
2044 
2045  if( p_group%domain_is_periodic(2) )then
2046  eta_min = p_group%lbf_grid%eta2_min
2047  eta_max = p_group%lbf_grid%eta2_max
2048  if( (y < eta_min) .or. (y >= eta_max) ) y = eta_min + modulo(y - eta_min, eta_max - eta_min)
2049  end if
2050  end subroutine periodic_correction
2051 
2054  function sll_f_add_element_in_int_list(head, new_element) result( new_list )
2055  type( sll_t_int_list_element ), pointer :: head, new_element
2056  type( sll_t_int_list_element ), pointer :: new_list
2057  new_element%next => head
2058  new_list => new_element
2059  end function
2060 
2064  subroutine convert_4d_index_to_1d( &
2065  k, &
2066  j_x, j_y, j_vx, j_vy, &
2067  n_parts_x, n_parts_y, n_parts_vx, n_parts_vy &
2068  )
2069  sll_int32, intent(out):: k
2070  sll_int32, intent(in) :: j_x
2071  sll_int32, intent(in) :: j_y
2072  sll_int32, intent(in) :: j_vx
2073  sll_int32, intent(in) :: j_vy
2074  sll_int32, intent(in) :: n_parts_x
2075  sll_int32, intent(in) :: n_parts_y
2076  sll_int32, intent(in) :: n_parts_vx
2077  sll_int32, intent(in) :: n_parts_vy
2078 
2079  if( j_x <= 0 .or. j_x > n_parts_x &
2080  .or. j_y <= 0 .or. j_y > n_parts_y &
2081  .or. j_vx <= 0 .or. j_vx > n_parts_vx &
2082  .or. j_vy <= 0 .or. j_vy > n_parts_vy )then
2083  k = 0
2084  else
2085  k = 1+ (j_vy-1) + (j_vx-1) * n_parts_vy + (j_y-1) * n_parts_vy * n_parts_vx + (j_x-1) * n_parts_vy * n_parts_vx * n_parts_y
2086  end if
2087 
2088  sll_assert( k <= n_parts_x * n_parts_y * n_parts_vx * n_parts_vy )
2089  end subroutine
2090 
2094  subroutine convert_1d_index_to_4d( &
2095  j_x, j_y, j_vx, j_vy, &
2096  k, &
2097  n_parts_x, n_parts_y, n_parts_vx, n_parts_vy &
2098  )
2099  sll_int32, intent(out) :: j_x
2100  sll_int32, intent(out) :: j_y
2101  sll_int32, intent(out) :: j_vx
2102  sll_int32, intent(out) :: j_vy
2103  sll_int32, intent(in) :: k
2104  sll_int32, intent(in) :: n_parts_x
2105  sll_int32, intent(in) :: n_parts_y
2106  sll_int32, intent(in) :: n_parts_vx
2107  sll_int32, intent(in) :: n_parts_vy
2108  sll_int32 :: k_aux
2109 
2110  k_aux = k-1
2112  j_vy = mod(k_aux, n_parts_vy) + 1
2113 
2114  k_aux = (k_aux - (j_vy-1)) / n_parts_vy
2116  j_vx = mod(k_aux, n_parts_vx) + 1
2117 
2118  k_aux = (k_aux - (j_vx-1)) / n_parts_vx
2120  j_y = mod(k_aux, n_parts_y) + 1
2121 
2122  k_aux = (k_aux - (j_y-1)) / n_parts_y
2124  j_x = k_aux + 1
2125 
2126  sll_assert(j_x <= n_parts_x)
2127  end subroutine
2128 
2133  subroutine get_1d_cell_containing_point( j, x, x_min, h )
2134  sll_int32, intent(out) :: j
2135  sll_real64, intent(in) :: x
2136  sll_real64, intent(in) :: x_min
2137  sll_real64, intent(in) :: h
2138 
2139  j = int( (x - x_min) / h ) + 1
2140  end subroutine
2141 
2145  subroutine update_closest_particle_arrays(k_part, &
2146  x_aux, y_aux, vx_aux, vy_aux, &
2147  i, j, l, m, &
2148  h_cell_x, h_cell_y, h_cell_vx, h_cell_vy, &
2149  closest_particle, &
2150  closest_particle_distance)
2151 
2152  sll_int32, intent(in) :: k_part
2153  sll_int32, intent(in) :: i, j, l, m
2154  sll_real64, intent(in) :: x_aux
2155  sll_real64, intent(in) :: y_aux
2156  sll_real64, intent(in) :: vx_aux
2157  sll_real64, intent(in) :: vy_aux
2158  sll_real64, intent(in) :: h_cell_x, h_cell_y, h_cell_vx, h_cell_vy
2159  sll_int32, dimension(:,:,:,:), intent(inout) :: closest_particle
2160  sll_real64, dimension(:,:,:,:), intent(inout) :: closest_particle_distance
2161  sll_real64 :: square_dist_to_cell_center
2162 
2163  square_dist_to_cell_center = (x_aux - (i + 0.5) * h_cell_x )**2. &
2164  + (y_aux - (j + 0.5) * h_cell_y )**2. &
2165  + (vx_aux - (l + 0.5) * h_cell_vx )**2. &
2166  + (vy_aux - (m + 0.5) * h_cell_vy )**2.
2167 
2169  if(closest_particle(i,j,l,m) == 0 .or. square_dist_to_cell_center < closest_particle_distance(i,j,l,m)) then
2170  closest_particle(i,j,l,m) = k_part
2171  closest_particle_distance(i,j,l,m) = square_dist_to_cell_center
2172  end if
2173 
2174  end subroutine update_closest_particle_arrays
2175 
2176 end module sll_m_particle_group_2d2v_lbf
sll_m_cartesian_meshes
Cartesian mesh basic types.
Definition: sll_m_cartesian_meshes.F90:20
sll_m_cartesian_meshes::sll_f_new_cartesian_mesh_4d
type(sll_t_cartesian_mesh_4d) function, pointer, public sll_f_new_cartesian_mesh_4d(num_cells1, num_cells2, num_cells3, num_cells4, eta1_min, eta1_max, eta2_min, eta2_max, eta3_min, eta3_max, eta4_min, eta4_max)
allocates the memory space for a new 4D cartesian mesh on the heap,
Definition: sll_m_cartesian_meshes.F90:866
sll_m_constants
Fortran module where set some physical and mathematical constants.
Definition: sll_m_constants.F90:23
sll_m_constants::sll_p_pi
real(kind=f64), parameter, public sll_p_pi
Definition: sll_m_constants.F90:51
sll_m_initial_distribution
Parameters to define common initial distributions.
Definition: sll_m_initial_distribution.F90:5
sll_m_particle_group_2d2v_lbf
Module for a particle group with linearized-backward-flow (lbf) resamplings.
Definition: sll_m_particle_group_2d2v_lbf.F90:24
sll_m_particle_group_2d2v_lbf::set_v_2d2v_lbf
subroutine set_v_2d2v_lbf(self, i, x)
Set the velocity of a particle.
Definition: sll_m_particle_group_2d2v_lbf.F90:255
sll_m_particle_group_2d2v_lbf::reset_particles_weights_with_direct_interpolation
subroutine reset_particles_weights_with_direct_interpolation(self, target_total_charge, enforce_total_charge)
reset_particles_weights_with_direct_interpolation ---------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:747
sll_m_particle_group_2d2v_lbf::sample
subroutine sample(self, target_total_charge, enforce_total_charge, init_f_params, rand_seed, rank, world_size)
sample (layer for resample procedure) ---------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:462
sll_m_particle_group_2d2v_lbf::get_x_2d2v_lbf
pure real(kind=f64) function, dimension(3) get_x_2d2v_lbf(self, i)
Getters ---------------------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:161
sll_m_particle_group_2d2v_lbf::reset_particles_positions
subroutine reset_particles_positions(self)
reset_particles_positions ---------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:557
sll_m_particle_group_2d2v_lbf::get_neighbor_index
integer(kind=i32) function get_neighbor_index(self, k, dim, dir)
get_neighbor_index ----------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:634
sll_m_particle_group_2d2v_lbf::delete_particle_group_2d2v_lbf
subroutine delete_particle_group_2d2v_lbf(self)
Destructor ------------------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:444
sll_m_particle_group_2d2v_lbf::write_known_density_on_remapping_grid
subroutine write_known_density_on_remapping_grid(self, distribution_params)
write_known_density_on_remapping_grid ---------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:532
sll_m_particle_group_2d2v_lbf::set_x_2d2v_lbf
subroutine set_x_2d2v_lbf(self, i, x)
Setters ---------------------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:244
sll_m_particle_group_2d2v_lbf::get_1d_cell_containing_point
subroutine get_1d_cell_containing_point(j, x, x_min, h)
get_1d_cell_containing_point ------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:2134
sll_m_particle_group_2d2v_lbf::get_ltp_deformation_matrix
subroutine get_ltp_deformation_matrix(self, k, mesh_period_x, mesh_period_y, h_particles_x, h_particles_y, h_particles_vx, h_particles_vy, x_k, y_k, vx_k, vy_k, d11, d12, d13, d14, d21, d22, d23, d24, d31, d32, d33, d34, d41, d42, d43, d44)
get_ltp_deformation_matrix --------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:1602
sll_m_particle_group_2d2v_lbf::sll_p_lbf_given_grid
integer(kind=i32), parameter sll_p_lbf_given_grid
Definition: sll_m_particle_group_2d2v_lbf.F90:61
sll_m_particle_group_2d2v_lbf::convert_4d_index_to_1d
subroutine convert_4d_index_to_1d(k, j_x, j_y, j_vx, j_vy, n_parts_x, n_parts_y, n_parts_vx, n_parts_vy)
convert_4d_index_to_1d ------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:2069
sll_m_particle_group_2d2v_lbf::update_closest_particle_arrays
subroutine update_closest_particle_arrays(k_part, x_aux, y_aux, vx_aux, vy_aux, i, j, l, m, h_cell_x, h_cell_y, h_cell_vx, h_cell_vy, closest_particle, closest_particle_distance)
update_closest_particle_arrays ----------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:2151
sll_m_particle_group_2d2v_lbf::reconstruct_f_lbf
subroutine reconstruct_f_lbf(self, reconstruction_set_type, given_grid_4d, given_array_4d, reconstruct_f_on_last_node, target_total_charge, enforce_total_charge)
macro for the lbf approximation of a transported density (reconstruct_f_lbf) below ------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:853
sll_m_particle_group_2d2v_lbf::get_common_weight_2d2v_lbf
pure real(kind=f64) function get_common_weight_2d2v_lbf(self)
Set the common weight.
Definition: sll_m_particle_group_2d2v_lbf.F90:232
sll_m_particle_group_2d2v_lbf::get_weights_2d2v_lbf
pure real(kind=f64) function, dimension(self%n_weights) get_weights_2d2v_lbf(self, i)
Get weights of a particle.
Definition: sll_m_particle_group_2d2v_lbf.F90:221
sll_m_particle_group_2d2v_lbf::initialize_particle_group_2d2v_lbf
subroutine initialize_particle_group_2d2v_lbf(self, species_charge, species_mass, domain_is_x_periodic, domain_is_y_periodic, remap_degree, remapping_grid_eta_min, remapping_grid_eta_max, remapping_sparse_grid_max_levels, n_particles_x, n_particles_y, n_particles_vx, n_particles_vy, n_weights)
Initializer -----------------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:303
sll_m_particle_group_2d2v_lbf::sll_s_new_particle_group_2d2v_lbf_ptr
subroutine, public sll_s_new_particle_group_2d2v_lbf_ptr(particle_group, species_charge, species_mass, domain_is_x_periodic, domain_is_y_periodic, remap_degree, remapping_grid_eta_min, remapping_grid_eta_max, remapping_sparse_grid_max_levels, n_particles_x, n_particles_y, n_particles_vx, n_particles_vy)
Constructor for abstract type.
Definition: sll_m_particle_group_2d2v_lbf.F90:401
sll_m_particle_group_2d2v_lbf::sll_p_lbf_remapping_grid
integer(kind=i32), parameter sll_p_lbf_remapping_grid
possible values for the parameter reconstruction_set_type in the reconstruction routine
Definition: sll_m_particle_group_2d2v_lbf.F90:60
sll_m_particle_group_2d2v_lbf::convert_1d_index_to_4d
subroutine convert_1d_index_to_4d(j_x, j_y, j_vx, j_vy, k, n_parts_x, n_parts_y, n_parts_vx, n_parts_vy)
convert_1d_index_to_4d ------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:2099
sll_m_particle_group_2d2v_lbf::resample
subroutine resample(self, target_total_charge, enforce_total_charge, init_f_params, rand_seed, rank, world_size)
resample (and sample) -------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:484
sll_m_particle_group_2d2v_lbf::set_weights_2d2v_lbf
subroutine set_weights_2d2v_lbf(self, i, x)
Set the weights of a particle.
Definition: sll_m_particle_group_2d2v_lbf.F90:266
sll_m_particle_group_2d2v_lbf::get_mass_2d2v_lbf
pure real(kind=f64) function get_mass_2d2v_lbf(self, i, i_weight)
Get mass of a particle ( m * particle_weight)
Definition: sll_m_particle_group_2d2v_lbf.F90:204
sll_m_particle_group_2d2v_lbf::reconstruct_f_lbf_on_given_grid
subroutine reconstruct_f_lbf_on_given_grid(self, given_grid_4d, given_array_4d, reconstruct_f_on_last_node, target_total_charge, enforce_total_charge)
reconstruct_f_lbf_on_remapping_grid -----------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:1541
sll_m_particle_group_2d2v_lbf::get_charge_2d2v_lbf
pure real(kind=f64) function get_charge_2d2v_lbf(self, i, i_weight)
Get charge of a particle (q * particle_weight)
Definition: sll_m_particle_group_2d2v_lbf.F90:187
sll_m_particle_group_2d2v_lbf::reconstruct_f_lbf_on_remapping_grid
subroutine reconstruct_f_lbf_on_remapping_grid(self)
reconstruct_f_lbf_on_remapping_grid -----------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:1507
sll_m_particle_group_2d2v_lbf::get_v_2d2v_lbf
pure real(kind=f64) function, dimension(3) get_v_2d2v_lbf(self, i)
Get the velocity of a particle.
Definition: sll_m_particle_group_2d2v_lbf.F90:173
sll_m_particle_group_2d2v_lbf::interpolate_value_of_remapped_f
real(kind=f64) function interpolate_value_of_remapped_f(self, eta)
interpolate_value_of_remapped_f ---------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:799
sll_m_particle_group_2d2v_lbf::set_common_weight_2d2v_lbf
subroutine set_common_weight_2d2v_lbf(self, x)
Set the common weight.
Definition: sll_m_particle_group_2d2v_lbf.F90:277
sll_m_particle_group_2d2v_lbf::periodic_correction
subroutine periodic_correction(p_group, x, y)
periodic_correction ---------------------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:2034
sll_m_particle_group_2d2v_lbf::sll_f_add_element_in_int_list
type(sll_t_int_list_element) function, pointer sll_f_add_element_in_int_list(head, new_element)
sll_f_add_element_in_int_list -----------------------------------------------------------------------...
Definition: sll_m_particle_group_2d2v_lbf.F90:2055
sll_m_particle_group_base
Definition: sll_m_particle_group_base.F90:1
sll_m_sparse_grid_4d
Implementation of a 4D sparse grid with interpolation routines.
Definition: sll_m_sparse_grid_4d.F90:6
sll_m_cartesian_meshes::sll_t_cartesian_mesh_4d
4D cartesian mesh
Definition: sll_m_cartesian_meshes.F90:132
sll_m_initial_distribution::sll_c_distribution_params
Abstract data type for parameters of initial distribution.
Definition: sll_m_initial_distribution.F90:48
sll_m_particle_group_2d2v_lbf::sll_t_int_list_element_ptr
a pointer type is needed for arrays of lists
Definition: sll_m_particle_group_2d2v_lbf.F90:70
sll_m_particle_group_2d2v_lbf::sll_t_int_list_element
linked lists of integers, used in some of the module routines
Definition: sll_m_particle_group_2d2v_lbf.F90:64
sll_m_particle_group_2d2v_lbf::sll_t_particle_group_2d2v_lbf
Group of sll_t_particle_group_2d2v_lbf.
Definition: sll_m_particle_group_2d2v_lbf.F90:76
sll_m_particle_group_base::sll_c_particle_group_base
Definition: sll_m_particle_group_base.F90:42
sll_m_particle_group_base::sll_t_species
Definition: sll_m_particle_group_base.F90:27
sll_m_sparse_grid_4d::sll_t_sparse_grid_4d
Sparse grid object for 4d with interpolation routines. Note in 4d we have only an implementation of a...
Definition: sll_m_sparse_grid_4d.F90:26
    Report Typos and Errors
Generated on Mon Oct 23 2023 19:15:40 for Semi-Lagrangian Library by   doxygen 1.9.1