check_if_in_fields_domain_p Subroutine

  subroutine check_if_in_fields_domain_p(F,Y_R,Y_PHI,Y_Z,flag)
    !! @note Subrotuine that checks if particles in the simulation are within
    !! the spatial domain where interpolants and fields are known. @endnote
    !! External fields and interpolants can have different spatial domains where
    !! they are defined. Therefore, it is necessary to
    !! check if a given particle has left these spatial domains to
    !! stop following it, otherwise this will cause an error in the simulation.
    TYPE(FIELDS), INTENT(IN)                                   :: F
    REAL(rp), DIMENSION(8),  INTENT(IN)      :: Y_R,Y_PHI,Y_Z    
    INTEGER(is), DIMENSION(8), INTENT(INOUT)  :: flag
    !! Flag that determines whether particles are followed in the
    !! simulation (flag=1), or not (flag=0).
    INTEGER                                                :: IR
    !! Variable used to localize the grid cell in the \((R,\phi,Z)\)
    !! or \((R,Z)\) grid containing the fields data that corresponds
    !! to the radial position of the particles.
    INTEGER                                                :: IPHI
    !! Variable used to localize the grid cell in the \((R,\phi,Z)\)
    !! or \((R,Z)\) grid containing the fields data that corresponds
    !! to the azimuthal position of the particles.
    INTEGER                                                :: IZ
    !! Variable used to localize the grid cell in the \((R,\phi,Z)\)
    !! or \((R,Z)\) grid containing the fields data that corresponds
    !! to the vertical position of the particles.
    INTEGER(ip)                                            :: pp
    !! Particle iterator.
    INTEGER(ip)                                            :: ss
    !! Species iterator.

    
!    write(6,'("YR:",E17.10)') Y_R
!    write(6,'("YPHI:",E17.10)') Y_PHI
!    write(6,'("YZ:",E17.10)') Y_Z

!    write(6,'("Ro:",E17.10)') fields_domain%Ro
!    write(6,'("Zo:",E17.10)') fields_domain%Zo
!    write(6,'("DR:",E17.10)') fields_domain%DR
    !    write(6,'("DZ:",E17.10)') fields_domain%DZ
!    write(6,'("DT:",E17.10)') fields_domain%DT


    
    if (ALLOCATED(fields_domain%FLAG3D)) then
       if (F%Dim2x1t) then
          !$OMP SIMD
          !       !$OMP&  aligned(IR,IPHI,IZ)
          do pp=1_idef,8_idef

             IR = INT(FLOOR((Y_R(pp)  - fields_domain%Ro + &
                  0.5_rp*fields_domain%DR)/fields_domain%DR) + 1.0_rp,idef)
             IPHI = INT(FLOOR((Y_PHI(pp)  - fields_domain%To &
                  + 0.5_rp*fields_domain%DT)/fields_domain%DT) + 1.0_rp,idef)
             IZ = INT(FLOOR((Y_Z(pp)  + ABS(fields_domain%Zo) + &
                  0.5_rp*fields_domain%DZ)/fields_domain%DZ) + 1.0_rp,idef)

!             write(6,'("IR: ",I16)') IR
!             write(6,'("IPHI: ",I16)') IPHI
!             write(6,'("IZ: ",I16)') IZ
             
             if ((fields_domain%FLAG3D(IR,IPHI,IZ).NE.1_is).OR. &
                  ((IR.GT.bfield_2X1T%NR).OR.(IZ.GT.bfield_2X1T%NZ))) then
                flag(pp) = 0_is

                !write(6,'("YR:",E17.10)') Y_R(pp)
                !write(6,'("YPHI:",E17.10)') Y_PHI(pp)
                !write(6,'("YZ:",E17.10)') Y_Z(pp)

                !write(6,'("IR: ",I16)') IR
                !write(6,'("IPHI: ",I16)') IPHI
                !write(6,'("IZ: ",I16)') IZ

                !call KORC_ABORT()

             end if

             !write(6,'("IPHI: ",I16)') IPHI
             !write(6,'("flag: ",I16)') flag(pp)


          end do
          !$OMP END SIMD
       else
          !$OMP SIMD
          !       !$OMP&  aligned(IR,IPHI,IZ)
          do pp=1_idef,8_idef

             IR = INT(FLOOR((Y_R(pp)  - fields_domain%Ro + &
                  0.5_rp*fields_domain%DR)/fields_domain%DR) + 1.0_rp,idef)
             IPHI = INT(FLOOR((Y_PHI(pp)  + 0.5_rp*fields_domain%DPHI)/ &
                  fields_domain%DPHI) + 1.0_rp,idef)
             IZ = INT(FLOOR((Y_Z(pp)  + ABS(fields_domain%Zo) + &
                  0.5_rp*fields_domain%DZ)/fields_domain%DZ) + 1.0_rp,idef)

             if ((fields_domain%FLAG3D(IR,IPHI,IZ).NE.1_is).OR. &
                  ((IR.GT.bfield_3d%NR).OR.(IZ.GT.bfield_3d%NZ))) then
                flag(pp) = 0_is

                !write(6,'("YR:",E17.10)') Y_R
                !write(6,'("YPHI:",E17.10)') Y_PHI
                !write(6,'("YZ:",E17.10)') Y_Z

                !write(6,'("IR: ",I16)') IR
                !write(6,'("IPHI: ",I16)') IPHI
                !write(6,'("IZ: ",I16)') IZ

                !call KORC_ABORT()

             end if

             !write(6,'("IPHI: ",I16)') IPHI
             !write(6,'("flag: ",I16)') flag(pp)


          end do
          !$OMP END SIMD
       end if
    else
       !$OMP SIMD
!       !$OMP& aligned(IR,IZ)
       do pp=1_idef,8_idef
          IR = INT(FLOOR((Y_R(pp)  - fields_domain%Ro + &
               0.5_rp*fields_domain%DR)/fields_domain%DR) + 1.0_rp,idef)
          IZ = INT(FLOOR((Y_Z(pp)  + ABS(fields_domain%Zo) + &
               0.5_rp*fields_domain%DZ)/fields_domain%DZ) + 1.0_rp,idef)

!          write(6,*) pp

!          write(6,'("Size of fields_domain R: ",I16)') &
!               size(fields_domain%FLAG2D,1)
!          write(6,'("Size of fields_domain Z: ",I16)') &
!               size(fields_domain%FLAG2D,2)   
          
!          if ((IR.lt.0).or.(IZ.lt.0)) then
!             write(6,'("YR:",E17.10)') Y_R(pp)
!             write(6,'("YZ:",E17.10)') Y_Z(pp)
!             write(6,'("IR: ",I16)') IR
!             write(6,'("IZ: ",I16)') IZ
!          end if
          
          if ((fields_domain%FLAG2D(IR,IZ).NE.1_is).OR. &
               ((IR.GT.bfield_2d%NR).OR.(IZ.GT.bfield_2d%NZ))) then
             flag(pp) = 0_is

!             write(6,'("Shit''s fucked.")')
          end if
       end do      
       !$OMP END SIMD
!       write(6,'("Shit''s not fucked.")')
    end if
  end subroutine check_if_in_fields_domain_p