advance_GCinterp_psi2x1t_vars Subroutine

public subroutine advance_GCinterp_psi2x1t_vars(vars, pp, tt, params, Y_R, Y_PHI, Y_Z, V_PLL, V_MU, q_cache, m_cache, flag_cache, F, P, B_R, B_PHI, B_Z, E_PHI, PSIp, curlb_R, curlb_PHI, curlb_Z, gradB_R, gradB_PHI, gradB_Z, ne)

Comment this section further with evolution equations, numerical methods, and descriptions of both. write(6,*) 'Z1',Y_Z(1)

Arguments

Type IntentOptional AttributesName
type(PARTICLES), intent(inout) :: vars
integer, intent(in) :: pp
integer(kind=ip), intent(in) :: tt

time iterator.

type(KORC_PARAMS), intent(inout) :: params

Core KORC simulation parameters.

real(kind=rp), intent(inout), DIMENSION(8):: Y_R
real(kind=rp), intent(inout), DIMENSION(8):: Y_PHI
real(kind=rp), intent(inout), DIMENSION(8):: Y_Z
real(kind=rp), intent(inout), DIMENSION(8):: V_PLL
real(kind=rp), intent(inout), DIMENSION(8):: V_MU
real(kind=rp), intent(in) :: q_cache
real(kind=rp), intent(in) :: m_cache
integer(kind=is), intent(inout), DIMENSION(8):: flag_cache
type(FIELDS), intent(in) :: F
type(PROFILES), intent(in) :: P
real(kind=rp), intent(out), DIMENSION(8):: B_R
real(kind=rp), intent(out), DIMENSION(8):: B_PHI
real(kind=rp), intent(out), DIMENSION(8):: B_Z
real(kind=rp), intent(out), DIMENSION(8):: E_PHI
real(kind=rp), intent(out), DIMENSION(8):: PSIp
real(kind=rp), intent(out), DIMENSION(8):: curlb_R
real(kind=rp), intent(out), DIMENSION(8):: curlb_PHI
real(kind=rp), intent(out), DIMENSION(8):: curlb_Z
real(kind=rp), intent(out), DIMENSION(8):: gradB_R
real(kind=rp), intent(out), DIMENSION(8):: gradB_PHI
real(kind=rp), intent(out), DIMENSION(8):: gradB_Z
real(kind=rp), intent(out), DIMENSION(8):: ne

Contents


Source Code

  subroutine advance_GCinterp_psi2x1t_vars(vars,pp,tt,params,Y_R,Y_PHI,Y_Z, &
       V_PLL,V_MU,q_cache,m_cache,flag_cache,F,P,B_R,B_PHI,B_Z,E_PHI,PSIp, &
       curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z,ne)
    !! @note Subroutine to advance GC variables \(({\bf X},p_\parallel)\)
    !! @endnote
    !! Comment this section further with evolution equations, numerical
    !! methods, and descriptions of both.
    TYPE(KORC_PARAMS), INTENT(INOUT)                              :: params
    !! Core KORC simulation parameters.
    TYPE(PARTICLES), INTENT(INOUT)     :: vars
    TYPE(PROFILES), INTENT(IN)                                 :: P
    TYPE(FIELDS), INTENT(IN)                                   :: F
    REAL(rp)                                      :: dt,time
    !! Time step used in the leapfrog step (\(\Delta t\)).

    INTEGER                                                    :: cc
    !! Chunk iterator.
    INTEGER(ip),intent(in)                                      :: tt
    !! time iterator.
    INTEGER,intent(in)                                  :: pp
    

    REAL(rp),DIMENSION(8)               :: Bmag
    REAL(rp)              :: a1 = 1./5._rp
    REAL(rp) :: a21 = 3./40._rp,a22=9./40._rp
    REAL(rp) :: a31 = 3./10._rp,a32=-9./10._rp,a33=6./5._rp
    REAL(rp) :: a41 = -11./54._rp,a42=5./2._rp,a43=-70./27._rp,a44=35./27._rp
    REAL(rp) :: a51 = 1631./55296._rp,a52=175./512._rp,a53=575./13824._rp,a54=44275./110592._rp,a55=253./4096._rp
    REAL(rp) :: b1=37./378._rp,b2=0._rp,b3=250./621._rp,b4=125./594._rp,b5=0._rp,b6=512./1771._rp

    REAL(rp),DIMENSION(8) :: k1_R,k1_PHI,k1_Z,k1_PLL,k1_MU
    REAL(rp),DIMENSION(8) :: k2_R,k2_PHI,k2_Z,k2_PLL,k2_MU
    REAL(rp),DIMENSION(8) :: k3_R,k3_PHI,k3_Z,k3_PLL,k3_MU
    REAL(rp),DIMENSION(8) :: k4_R,k4_PHI,k4_Z,k4_PLL,k4_MU
    REAL(rp),DIMENSION(8) :: k5_R,k5_PHI,k5_Z,k5_PLL,k5_MU
    REAL(rp),DIMENSION(8) :: k6_R,k6_PHI,k6_Z,k6_PLL,k6_MU
    REAL(rp),DIMENSION(8) :: Y0_R,Y0_PHI,Y0_Z
    REAL(rp),DIMENSION(8),INTENT(INOUT) :: Y_R,Y_PHI,Y_Z
    REAL(rp),DIMENSION(8),INTENT(OUT) :: B_R,B_PHI,B_Z
    REAL(rp),DIMENSION(8) :: E_R,E_Z
    REAL(rp),DIMENSION(8),INTENT(OUT) :: E_PHI
    REAL(rp),DIMENSION(8),INTENT(OUT) :: PSIp
    REAL(rp),DIMENSION(8),INTENT(OUT) :: ne
    REAL(rp),DIMENSION(8),INTENT(OUT) :: curlb_R,curlb_PHI,curlb_Z
    REAL(rp),DIMENSION(8),INTENT(OUT) :: gradB_R,gradB_PHI,gradB_Z
    REAL(rp),DIMENSION(8),INTENT(INOUT) :: V_PLL,V_MU
    REAL(rp),DIMENSION(8) :: RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU
    REAL(rp),DIMENSION(8) :: V0_PLL,V0_MU
    REAL(rp),DIMENSION(8) :: Te,Zeff

    INTEGER(is),DIMENSION(8),intent(INOUT) :: flag_cache
    REAL(rp),intent(IN)  :: q_cache,m_cache

    dt=params%dt
    time=params%init_time+(params%it-1+tt)*params%dt
    
    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU)
    do cc=1_idef,8_idef

       Y0_R(cc)=Y_R(cc)
       Y0_PHI(cc)=Y_PHI(cc)
       Y0_Z(cc)=Y_Z(cc)
       V0_PLL(cc)=V_PLL(cc)
       V0_MU(cc)=V_MU(cc)



    end do
    !$OMP END SIMD

!    write(6,*) 'R0',Y_R(1)
!    write(6,*) 'PHI0',Y_PHI(1)
!    write(6,*) 'Z0',Y_Z(1)
!    write(6,*) 'PPLL0',V_PLL(1)
!    write(6,*) 'MU0',V_MU(1)
    
    
!    call interp_fields_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)


    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne) 

!    write(6,*) 'R0',Y_R(1)
!    write(6,*) 'PHI0',Y_PHI(1)
!    write(6,*) 'Z0',Y_Z(1)
!    write(6,*) 'PPLL0',V_PLL(1)
!    write(6,*) 'MU0',V_MU(1)
    
!    write(6,*) 'BR',B_R(1)
!    write(6,*) 'BPHI',B_PHI(1)
!    write(6,*) 'BZ',B_Z(1)

!    write(6,*) 'gradBR',gradB_R(1)
!    write(6,*) 'gradBPHI',gradB_PHI(1)
!    write(6,*) 'gradBZ',gradB_Z(1)

!    write(6,*) 'curlBR',curlB_R(1)
!    write(6,*) 'curlBPHI',curlB_PHI(1)
!    write(6,*) 'curlBZ',curlB_Z(1)

!    write(6,*) 'dt',params%dt
!    write(6,*) 'RHS_R',RHS_R(1)
!    write(6,*) 'RHS_PHI',RHS_PHI(1)
!    write(6,*) 'RHS_Z',RHS_Z(1)
!    write(6,*) 'RHS_PLL',RHS_PLL(1)
!    write(6,*) 'RHS_MU',RHS_MU(1)
    
    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU, &
!    !$OMP& RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU, &
!    !$OMP& k1_R,k1_PHI,k1_Z,k1_PLL,k1_MU)
    do cc=1_idef,8
       k1_R(cc)=dt*RHS_R(cc)              
       k1_PHI(cc)=dt*RHS_PHI(cc)    
       k1_Z(cc)=dt*RHS_Z(cc)    
       k1_PLL(cc)=dt*RHS_PLL(cc)
       k1_MU(cc)=dt*RHS_MU(cc)    
       
       Y_R(cc)=Y0_R(cc)+a1*k1_R(cc)
       Y_PHI(cc)=Y0_PHI(cc)+a1*k1_PHI(cc)
       Y_Z(cc)=Y0_Z(cc)+a1*k1_Z(cc)
       V_PLL(cc)=V0_PLL(cc)   +a1*k1_PLL(cc)
       V_MU(cc)=V0_MU(cc)   +a1*k1_MU(cc)
       
    end do
    !$OMP END SIMD

!    write(6,*) 'R1',Y_R(1)
!    write(6,*) 'PHI1',Y_PHI(1)
!!    write(6,*) 'Z1',Y_Z(1)
 !   write(6,*) 'PPLL1',V_PLL(1)
 !   write(6,*) 'MU1',V_MU(1)
    
!    call interp_fields_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)




    
    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne) 

    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU, &
!    !$OMP& RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU, &
!    !$OMP& k2_R,k2_PHI,k2_Z,k2_PLL,k2_MU)
    do cc=1_idef,8
       k2_R(cc)=dt*RHS_R(cc)    
       k2_PHI(cc)=dt*RHS_PHI (cc)   
       k2_Z(cc)=dt*RHS_Z(cc)   
       k2_PLL(cc)=dt*RHS_PLL(cc)
       k2_MU(cc)=dt*RHS_MU(cc)

       Y_R(cc)=Y0_R(cc)+a21*k1_R(cc)+a22*k2_R(cc)
       Y_PHI(cc)=Y0_PHI(cc)+a21*k1_PHI(cc)+a22*k2_PHI(cc)
       Y_Z(cc)=Y0_Z(cc)+a21*k1_Z(cc)+a22*k2_Z(cc)
       V_PLL(cc)=V0_PLL(cc)   +a21*k1_PLL(cc)+a22*k2_PLL(cc)
       V_MU(cc)=V0_MU(cc)   +a21*k1_MU(cc)+a22*k2_MU(cc)


    end do
    !$OMP END SIMD

!    call interp_fields_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)

 


    
    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne)

    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU, &
!    !$OMP& RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU, &
!    !$OMP& k3_R,k3_PHI,k3_Z,k3_PLL,k3_MU)
    do cc=1_idef,8
       k3_R(cc)=dt*RHS_R(cc)   
       k3_PHI(cc)=dt*RHS_PHI(cc)    
       k3_Z(cc)=dt*RHS_Z(cc)    
       k3_PLL(cc)=dt*RHS_PLL(cc)
       k3_MU(cc)=dt*RHS_MU(cc)

       Y_R(cc)=Y0_R(cc)+a31*k1_R(cc)+a32*k2_R(cc)+a33*k3_R(cc)
       Y_PHI(cc)=Y0_PHI(cc)+a31*k1_PHI(cc)+a32*k2_PHI(cc)+ &
            a33*k3_PHI(cc)
       Y_Z(cc)=Y0_Z(cc)+a31*k1_Z(cc)+a32*k2_Z(cc)+a33*k3_Z(cc)

       V_PLL(cc)=V0_PLL(cc)   +a31*k1_PLL(cc)+a32*k2_PLL(cc)+a33*k3_PLL(cc)
       V_MU(cc)=V0_MU(cc)   +a31*k1_MU(cc)+a32*k2_MU(cc)+a33*k3_MU(cc)

 
    end do
    !$OMP END SIMD

!    call interp_fields_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)




    
    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne)     

    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU, &
!    !$OMP& RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU, &
!    !$OMP& k4_R,k4_PHI,k4_Z,k4_PLL,k4_MU)
    do cc=1_idef,8
       k4_R(cc)=dt*RHS_R(cc)   
       k4_PHI(cc)=dt*RHS_PHI(cc)    
       k4_Z(cc)=dt*RHS_Z(cc)    
       k4_PLL(cc)=dt*RHS_PLL(cc)
       k4_MU(cc)=dt*RHS_MU(cc)

       Y_R(cc)=Y0_R(cc)+a41*k1_R(cc)+a42*k2_R(cc)+a43*k3_R(cc)+ &
            a44*k4_R(cc)
       Y_PHI(cc)=Y0_PHI(cc)+a41*k1_PHI(cc)+a42*k2_PHI(cc)+ &
            a43*k3_PHI(cc)+a44*k4_PHI(cc)
       Y_Z(cc)=Y0_Z(cc)+a41*k1_Z(cc)+a42*k2_Z(cc)+a43*k3_Z(cc)+ &
            a44*k4_Z(cc)
       V_PLL(cc)=V0_PLL(cc)   +a41*k1_PLL(cc)+a42*k2_PLL(cc)+ &
            a43*k3_PLL(cc)+a44*k4_PLL(cc)
       V_MU(cc)=V0_MU(cc)   +a41*k1_MU(cc)+a42*k2_MU(cc)+ &
            a43*k3_MU(cc)+a44*k4_MU(cc)

     
    end do
    !$OMP END SIMD


!    call interp_fields_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)




    
    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne)   

    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU, &
!    !$OMP& RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU, &
!    !$OMP& k5_R,k5_PHI,k5_Z,k5_PLL,k5_MU)
    do cc=1_idef,8
       k5_R(cc)=dt*RHS_R(cc)    
       k5_PHI(cc)=dt*RHS_PHI(cc)    
       k5_Z(cc)=dt*RHS_Z(cc)    
       k5_PLL(cc)=dt*RHS_PLL(cc)
       k5_MU(cc)=dt*RHS_MU(cc)

       Y_R(cc)=Y0_R(cc)+a51*k1_R(cc)+a52*k2_R(cc)+a53*k3_R(cc)+ &
            a54*k4_R(cc)+a55*k5_R(cc)
       Y_PHI(cc)=Y0_PHI(cc)+a51*k1_PHI(cc)+a52*k2_PHI(cc)+ &
            a53*k3_PHI(cc)+a54*k4_PHI(cc)+a55*k5_PHI(cc)
       Y_Z(cc)=Y0_Z(cc)+a51*k1_Z(cc)+a52*k2_Z(cc)+a53*k3_Z(cc)+ &
            a54*k4_Z(cc)+a55*k5_Z(cc)
       V_PLL(cc)=V0_PLL(cc)   +a51*k1_PLL(cc)+a52*k2_PLL(cc)+ &
            a53*k3_PLL(cc)+a54*k4_PLL(cc)+a55*k5_PLL(cc)
       V_MU(cc)=V0_MU(cc)   +a51*k1_MU(cc)+a52*k2_MU(cc)+ &
            a53*k3_MU(cc)+a54*k4_MU(cc)+a55*k5_MU(cc)

     
    end do
    !$OMP END SIMD

!    call interp_fields_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)

    
    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne)         

    !$OMP SIMD
!    !$OMP& aligned(Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU,Y_R,Y_PHI,Y_Z,V_PLL,V_MU, &
!    !$OMP& RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU, &
!    !$OMP& k6_R,k6_PHI,k6_Z,k6_PLL,k6_MU)
    do cc=1_idef,8
       k6_R(cc)=dt*RHS_R(cc)    
       k6_PHI(cc)=dt*RHS_PHI(cc)    
       k6_Z(cc)=dt*RHS_Z(cc)    
       k6_PLL(cc)=dt*RHS_PLL(cc)
       k6_MU(cc)=dt*RHS_MU(cc)

       Y_R(cc)=Y0_R(cc)+b1*k1_R(cc)+b2*k2_R(cc)+ &
            b3*k3_R(cc)+b4*k4_R(cc)+b5*k5_R(cc)+b6*k6_R(cc)
       Y_PHI(cc)=Y0_PHI(cc)+b1*k1_PHI(cc)+b2*k2_PHI(cc)+ &
            b3*k3_PHI(cc)+b4*k4_PHI(cc)+b5*k5_PHI(cc)+b6*k6_PHI(cc)
       Y_Z(cc)=Y0_Z(cc)+b1*k1_Z(cc)+b2*k2_Z(cc)+ &
            b3*k3_Z(cc)+b4*k4_Z(cc)+b5*k5_Z(cc)+b6*k6_Z(cc)
       V_PLL(cc)=V0_PLL(cc)+b1*k1_PLL(cc)+b2*k2_PLL(cc)+ &
            b3*k3_PLL(cc)+b4*k4_PLL(cc)+b5*k5_PLL(cc)+b6*k6_PLL(cc)
       V_MU(cc)=V0_MU(cc)+b1*k1_MU(cc)+b2*k2_MU(cc)+ &
            b3*k3_MU(cc)+b4*k4_MU(cc)+b5*k5_MU(cc)+b6*k6_MU(cc)

     
    end do
    !$OMP END SIMD

    !$OMP SIMD
    !    !$OMP& aligned(Y_R,Y_PHI,Y_Z,V_PLL,V_MU,Y0_R,Y0_PHI,Y0_Z,V0_PLL,V0_MU)
    do cc=1_idef,8

       if (flag_cache(cc).eq.0_is) then
          Y_R(cc)=Y0_R(cc)
          Y_PHI(cc)=Y0_PHI(cc)
          Y_Z(cc)=Y0_Z(cc)
          V_PLL(cc)=V0_PLL(cc)
          V_MU(cc)=V0_MU(cc)
       end if          
 
    end do
    !$OMP END SIMD

    call calculate_GCfields_2x1t_p(F,Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,E_R,E_PHI, &
         E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R,gradB_PHI,gradB_Z, &
         flag_cache,PSIp,time)

    call add_analytical_E_p(params,tt,F,E_PHI,Y_R)
    
       
    
    call GCEoM1_p(tt,P,F,params,RHS_R,RHS_PHI,RHS_Z,RHS_PLL,RHS_MU,B_R,B_PHI, &
         B_Z,E_R,E_PHI,E_Z,curlb_R,curlb_PHI,curlb_Z,gradB_R, &
         gradB_PHI,gradB_Z,V_PLL,V_MU,Y_R,Y_Z,q_cache,m_cache,PSIp,ne) 

    !$OMP SIMD
    do cc=1_idef,8
       vars%RHS(pp-1+cc,1)=RHS_R(cc)
       vars%RHS(pp-1+cc,2)=RHS_PHI(cc)
       vars%RHS(pp-1+cc,3)=RHS_Z(cc)
       vars%RHS(pp-1+cc,4)=RHS_PLL(cc)
       vars%RHS(pp-1+cc,5)=RHS_MU(cc)
    end do
    !$OMP END SIMD       
    
    if (params%collisions) then       
       
       call include_CoulombCollisions_GC_p(tt,params,Y_R,Y_PHI,Y_Z, &
            V_PLL,V_MU,m_cache,flag_cache,F,P,E_PHI,ne,PSIp)

    end if


  end subroutine advance_GCinterp_psi2x1t_vars