subroutine analytical_fields_GC_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,PSIp)
TYPE(FIELDS), INTENT(IN) :: F
REAL(rp),DIMENSION(8),INTENT(IN) :: Y_R,Y_PHI,Y_Z
REAL(rp),DIMENSION(8),INTENT(OUT) :: B_R,B_PHI,B_Z
REAL(rp),DIMENSION(8),INTENT(OUT) :: gradB_R,gradB_PHI,gradB_Z
REAL(rp),DIMENSION(8),INTENT(OUT) :: curlB_R,curlB_PHI,curlB_Z
REAL(rp),DIMENSION(8),INTENT(OUT) :: E_R,E_PHI,E_Z
REAL(rp),DIMENSION(8),INTENT(OUT) :: PSIp
REAL(rp),DIMENSION(8) :: dRBR,dRBPHI,dRBZ,dZBR,dZBPHI,dZBZ,Bmag,dRbhatPHI
REAL(rp),DIMENSION(8) :: dRbhatZ,dZbhatR,dZbhatPHI,qprof,rm,theta
REAL(rp) :: B0,E0,lam,R0,q0
integer(ip) :: cc
B0=F%Bo
E0=F%Eo
lam=F%AB%lambda
R0=F%AB%Ro
q0=F%AB%qo
!$OMP SIMD
! !$OMP& aligned(Y_R,Y_PHI,Y_Z,B_R,B_PHI,B_Z,gradB_R,gradB_PHI,gradB_Z, &
! !$OMP& curlB_R,curlB_PHI,curlB_Z,E_R,E_PHI,E_Z,PSIp)
do cc=1_idef,8_idef
rm(cc)=sqrt((Y_R(cc)-R0)*(Y_R(cc)-R0)+Y_Z(cc)*Y_Z(cc))
theta(cc)=atan2(Y_Z(cc),(Y_R(cc)-R0))
qprof(cc) = 1.0_rp + (rm(cc)*rm(cc)/(lam*lam))
PSIp(cc)=Y_R(cc)*lam**2*B0/ &
(2*q0*(R0+rm(cc)*cos(theta(cc))))* &
log(1+(rm(cc)/lam)**2)
B_R(cc)=B0*Y_Z(cc)/(q0*qprof(cc)*Y_R(cc))
B_PHI(cc)=-B0*R0/Y_R(cc)
B_Z(cc)=-B0*(Y_R(cc)-R0)/(q0*qprof(cc)*Y_R(cc))
dRBR(cc)=-B0*Y_Z(cc)/(q0*qprof(cc)*Y_R(cc))*(1./Y_R(cc)+ &
2*(Y_R(cc)-R0)/(lam*lam*qprof(cc)))
dRBPHI(cc)=B0*R0/(Y_R(cc)*Y_R(cc))
dRBZ(cc)=B0/(q0*qprof(cc)*Y_R(cc))*(-R0/Y_R(cc)+2*(Y_R(cc)- &
R0)*(Y_R(cc)-R0)/(lam*lam*qprof(cc)))
dZBR(cc)=B0/(q0*qprof(cc)*Y_R(cc))*(1-2*Y_Z(cc)*Y_Z(cc)/ &
(lam*lam*qprof(cc)))
dZBPHI(cc)=0._rp
dZBZ(cc)=B0*(Y_R(cc)-R0)/(q0*Y_R(cc))*2*Y_Z(cc)/ &
(lam*lam*qprof(cc)*qprof(cc))
Bmag(cc)=sqrt(B_R(cc)*B_R(cc)+B_PHI(cc)*B_PHI(cc)+B_Z(cc)*B_Z(cc))
gradB_R(cc)=(B_R(cc)*dRBR(cc)+B_PHI(cc)*dRBPHI(cc)+B_Z(cc)*dRBZ(cc))/ &
Bmag(cc)
gradB_PHI(cc)=0._rp
gradB_Z(cc)=(B_R(cc)*dZBR(cc)+B_PHI(cc)*dZBPHI(cc)+B_Z(cc)*dZBZ(cc))/ &
Bmag(cc)
dRbhatPHI(cc)=(Bmag(cc)*dRBPHI(cc)-B_PHI(cc)*gradB_R(cc))/ &
(Bmag(cc)*Bmag(cc))
dRbhatZ(cc)=(Bmag(cc)*dRBZ(cc)-B_Z(cc)*gradB_R(cc))/(Bmag(cc)*Bmag(cc))
dZbhatR(cc)=(Bmag(cc)*dZBR(cc)-B_R(cc)*gradB_Z(cc))/(Bmag(cc)*Bmag(cc))
dZbhatPHI(cc)=(Bmag(cc)*dZBPHI(cc)-B_PHI(cc)*gradB_Z(cc))/ &
(Bmag(cc)*Bmag(cc))
curlb_R(cc)=-dZbhatPHI(cc)
curlb_PHI(cc)=dZbhatR(cc)-dRbhatZ(cc)
curlb_Z(cc)=B_PHI(cc)/(Bmag(cc)*Y_R(cc))+dRbhatPHI(cc)
E_R(cc) = 0.0_rp
E_PHI(cc) = E0*R0/Y_R(cc)
E_Z(cc) = 0.0_rp
end do
!$OMP END SIMD
end subroutine analytical_fields_GC_p