interp_2D_efields Subroutine

private subroutine interp_2D_efields(params, Y, E, flag)

@brief Subroutine for interpolating the pre-computed, axisymmetric electric field to the particles' position.

@param[in] Y Particles' position in cylindrical coordinates, Y(1,:) = , Y(2,:) = , and Y(3,:) = . @param[in,out] E Cartesian components of interpolated electric field components. E(1,:)=, E(2,:)=, and E(3,:)=. @param F Cylindrical components of interpolated magnetic field components. F(1,:)=, F(2,:)=, and F(3,:)=. @param flag Flag that indicates whether particles are followed in the simulation (flag=1), or not (flag=0). @param pp Particle iterator. @param ss Species iterator.

Arguments

Type IntentOptional AttributesName
type(KORC_PARAMS), intent(in) :: params

Core KORC simulation parameters.

real(kind=rp), intent(in), DIMENSION(:,:), ALLOCATABLE:: Y
real(kind=rp), intent(inout), DIMENSION(:,:), ALLOCATABLE:: E
integer(kind=is), intent(inout), DIMENSION(:), ALLOCATABLE:: flag

Contents

Source Code


Source Code

subroutine interp_2D_efields(params,Y,E,flag)
  TYPE(KORC_PARAMS), INTENT(IN)      :: params
  !! Core KORC simulation parameters.
  REAL(rp), DIMENSION(:,:), ALLOCATABLE, INTENT(IN)      :: Y
  REAL(rp), DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT)   :: E
  REAL(rp), DIMENSION(:,:), ALLOCATABLE                  :: F
  INTEGER(is), DIMENSION(:), ALLOCATABLE, INTENT(INOUT)  :: flag
  INTEGER                                                :: pp
  INTEGER                                                :: ss

!  write(6,*) 'interp E fields'
  
  if (Y(2,1).eq.0) then
     ss=1_idef
  else
     ss = size(Y,1)
  end if

  ALLOCATE(F(ss,3))
  !$OMP PARALLEL DO FIRSTPRIVATE(ss) PRIVATE(pp,ezerr) &
  !$OMP& SHARED(F,Y,E,flag,efield_2d)
  do pp=1_idef,ss
     if ( flag(pp) .EQ. 1_is ) then
        call EZspline_interp(efield_2d%R, Y(pp,1), Y(pp,3), F(pp,1), ezerr)
        call EZspline_error(ezerr)

        if (ezerr .NE. 0) then ! We flag the particle as lost
           flag(pp) = 0_is
        end if

        call EZspline_interp(efield_2d%PHI, Y(pp,1), Y(pp,3), F(pp,2), ezerr)
        call EZspline_error(ezerr)

        call EZspline_interp(efield_2d%Z, Y(pp,1), Y(pp,3), F(pp,3), ezerr)
        call EZspline_error(ezerr)

        if (.not.params%GC_coords) then
           E(pp,1) = F(pp,1)*COS(Y(pp,2)) - F(pp,2)*SIN(Y(pp,2))
           E(pp,2) = F(pp,1)*SIN(Y(pp,2)) + F(pp,2)*COS(Y(pp,2))
           E(pp,3) = F(pp,3)
        else
           E(pp,1) = F(pp,1)
           E(pp,2) = F(pp,2)
           E(pp,3) = F(pp,3)
        end if

        !write(6,*) 'EPHI',E(pp,2)
     end if     
     
  end do
  !$OMP END PARALLEL DO
  DEALLOCATE(F)
end subroutine interp_2D_efields