!=================== Original source is in jsrc/ =================

MODULE field_mod
  USE prec, ONLY : rstd
  USE layout_mod
  IMPLICIT NONE

  INTEGER,PARAMETER :: type_integer=1
  INTEGER,PARAMETER :: type_real=2
  INTEGER,PARAMETER :: type_logical=3

  TYPE t_field
     CHARACTER(30)      :: name
     LOGICAL :: ondevice !< flag if field is allocated on device as well
     INTEGER :: ndim
     INTEGER :: field_type
     INTEGER :: data_type
     INTEGER :: dim3
     INTEGER :: dim4
     REAL(rstd), POINTER, CONTIGUOUS :: rval2d(:) => NULL()
     REAL(rstd), POINTER, CONTIGUOUS :: rval3d(:,:) => NULL()
     REAL(rstd), POINTER, CONTIGUOUS :: rval4d(:,:,:) => NULL()
     INTEGER, POINTER, CONTIGUOUS :: ival2d(:) => NULL()
     INTEGER, POINTER, CONTIGUOUS :: ival3d(:,:) => NULL()
     INTEGER, POINTER, CONTIGUOUS :: ival4d(:,:,:) => NULL()
     LOGICAL, POINTER, CONTIGUOUS :: lval2d(:) => NULL()
     LOGICAL, POINTER, CONTIGUOUS :: lval3d(:,:) => NULL()
     LOGICAL, POINTER, CONTIGUOUS :: lval4d(:,:,:) => NULL()
  END TYPE t_field

  INTERFACE get_val
     MODULE PROCEDURE getval_r2d
     MODULE PROCEDURE getval_r3d
     MODULE PROCEDURE getval_r4d
     MODULE PROCEDURE getval_i2d
     MODULE PROCEDURE getval_i3d
     MODULE PROCEDURE getval_i4d
     MODULE PROCEDURE getval_l2d
     MODULE PROCEDURE getval_l3d
     MODULE PROCEDURE getval_l4d
  END INTERFACE get_val

  INTERFACE ASSIGNMENT(=)
     MODULE PROCEDURE getval_r2d
     MODULE PROCEDURE getval_r3d
     MODULE PROCEDURE getval_r4d
     MODULE PROCEDURE getval_i2d
     MODULE PROCEDURE getval_i3d
     MODULE PROCEDURE getval_i4d
     MODULE PROCEDURE getval_l2d
     MODULE PROCEDURE getval_l3d
     MODULE PROCEDURE getval_l4d
  END INTERFACE ASSIGNMENT(=)

  ! This module is PUBLIC, we do not want to propagate symbols from modules usedd at module level
  PRIVATE :: allocate_field_, deallocate_field_

CONTAINS

  !====================================== PRIVATE : allocate_field_ ===================================
  SUBROUTINE allocate_field_(layout, assigned_dom, field, field_type, data_type, dim3, dim4, name, ondevice)
    USE layout_mod, ONLY : t_layout, order_ij_l_p, order_l_ij_p
    USE domain_mod
    USE omp_para
    TYPE(t_layout), INTENT(IN) :: layout(:)
    LOGICAL, INTENT(IN)        :: assigned_dom(:)
    TYPE(t_field)              :: field(:)
    INTEGER, INTENT(IN)        :: field_type
    INTEGER, INTENT(IN)        :: data_type
    INTEGER, OPTIONAL          :: dim3, dim4
    CHARACTER(*), OPTIONAL     :: name
    LOGICAL, INTENT(IN), OPTIONAL :: ondevice

    INTEGER :: ind, ij_size

    DO ind=1,SIZE(field)
       IF (.NOT. assigned_dom(ind)  .OR. .NOT. is_omp_level_master) CYCLE

       IF(PRESENT(name)) THEN
          field(ind)%name = name
       ELSE
          field(ind)%name = '(undefined)'
       END IF

       IF (PRESENT(dim4)) THEN
          field(ind)%ndim=4
          field(ind)%dim4=dim4
          field(ind)%dim3=dim3
       ELSE IF (PRESENT(dim3)) THEN
          field(ind)%ndim=3
          field(ind)%dim3=dim3
          field(ind)%dim4=1
       ELSE
          field(ind)%ndim=2
          field(ind)%dim3=1
          field(ind)%dim4=1
       ENDIF

       field(ind)%data_type=data_type
       field(ind)%field_type=field_type
       ij_size = layout(ind)%ij_size

       SELECT CASE( layout(ind)%index_order )
       CASE( order_ij_l_p )
          IF (data_type==type_integer) ALLOCATE(field(ind)%ival4d(ij_size, field(ind)%dim3, field(ind)%dim4))
          IF (data_type==type_real)    ALLOCATE(field(ind)%rval4d(ij_size, field(ind)%dim3, field(ind)%dim4))
          IF (data_type==type_logical) ALLOCATE(field(ind)%lval4d(ij_size, field(ind)%dim3, field(ind)%dim4))

          IF (field(ind)%ndim==2) THEN
             IF (data_type==type_integer) field(ind)%ival2d => field(ind)%ival4d(:,1,1)
             IF (data_type==type_real)    field(ind)%rval2d => field(ind)%rval4d(:,1,1)
             IF (data_type==type_logical) field(ind)%lval2d => field(ind)%lval4d(:,1,1)
          ENDIF

       CASE( order_l_ij_p )
          IF (data_type==type_integer) ALLOCATE(field(ind)%ival4d(field(ind)%dim3, ij_size, field(ind)%dim4))
          IF (data_type==type_real)    ALLOCATE(field(ind)%rval4d(field(ind)%dim3, ij_size, field(ind)%dim4))
          IF (data_type==type_logical) ALLOCATE(field(ind)%lval4d(field(ind)%dim3, ij_size, field(ind)%dim4))

          IF (field(ind)%ndim==2) THEN
             IF (data_type==type_integer) field(ind)%ival2d => field(ind)%ival4d(1,:,1)
             IF (data_type==type_real)    field(ind)%rval2d => field(ind)%rval4d(1,:,1)
             IF (data_type==type_logical) field(ind)%lval2d => field(ind)%lval4d(1,:,1)
          ENDIF
       END SELECT

       IF (field(ind)%ndim==3) THEN
          IF (data_type==type_integer) field(ind)%ival3d => field(ind)%ival4d(:,:,1)
          IF (data_type==type_real)    field(ind)%rval3d => field(ind)%rval4d(:,:,1)
          IF (data_type==type_logical) field(ind)%lval3d => field(ind)%lval4d(:,:,1)
       END IF


       field(ind)%ondevice = .FALSE.
       IF (PRESENT(ondevice)) THEN
          IF (ondevice) CALL create_device_field(field(ind))
       END IF

    END DO

  END SUBROUTINE allocate_field_

  !====================================== PUBLIC : allocate_field ===================================

  SUBROUTINE allocate_field_glo(field, field_type, data_type, dim3, dim4, name)
    USE layout_mod, ONLY : layout => layout_glo
    IMPLICIT NONE
    TYPE(t_field), POINTER :: field(:)
    INTEGER, INTENT(IN) :: field_type
    INTEGER, INTENT(IN) :: data_type
    INTEGER, OPTIONAL   :: dim3, dim4
    CHARACTER(*), OPTIONAL :: name
    INTEGER :: ndom
    LOGICAL, ALLOCATABLE :: assigned_dom(:)

    ndom = SIZE(layout,1)
    ALLOCATE(assigned_dom(ndom))
    ! ONLY the master thread is allowed to call this routine

    ALLOCATE(field(ndom))
    assigned_dom(:) = .TRUE.
    CALL allocate_field_(layout(:,field_type), assigned_dom, field, field_type, data_type, dim3, dim4, name)
    DEALLOCATE(assigned_dom)

  END SUBROUTINE allocate_field_glo

  SUBROUTINE allocate_field(field, field_type, data_type, dim3, dim4, name, ondevice)
    USE layout_mod, ONLY : layout
    USE domain_mod, ONLY : assigned_domain
    USE omp_para, ONLY : is_omp_level_master
    TYPE(t_field), POINTER :: field(:)
    INTEGER, INTENT(IN)    :: field_type
    INTEGER, INTENT(IN)    :: data_type
    INTEGER, OPTIONAL      :: dim3, dim4
    CHARACTER(*), OPTIONAL :: name
    LOGICAL, INTENT(IN), OPTIONAL :: ondevice
    INTEGER :: ndom

    ndom = SIZE(layout,1)
    !$OMP BARRIER
    !$OMP MASTER
    ALLOCATE(field(ndom))
    !$OMP END MASTER
    !$OMP BARRIER
    CALL allocate_field_(layout(:,field_type), assigned_domain, field, field_type, data_type, dim3, dim4, name, ondevice)
    !$OMP BARRIER

  END SUBROUTINE allocate_field

  SUBROUTINE allocate_fields(nfield, field, field_type, data_type, dim3, dim4, name, ondevice)
    USE layout_mod, ONLY : layout
    USE domain_mod, ONLY : assigned_domain
    USE omp_para, ONLY : is_omp_level_master
    INTEGER, INTENT(IN)     :: nfield
    TYPE(t_field), POINTER  :: field(:,:)
    INTEGER, INTENT(IN)     :: field_type
    INTEGER, INTENT(IN)     :: data_type
    INTEGER, OPTIONAL       :: dim3, dim4
    CHARACTER(*), OPTIONAL  :: name
    LOGICAL, INTENT(IN), OPTIONAL :: ondevice
    INTEGER :: i, ndom

    ndom = SIZE(layout,1)

    !$OMP BARRIER
    !$OMP MASTER
    ALLOCATE(field(ndom,nfield))
    !$OMP END MASTER
    !$OMP BARRIER
    DO i=1,nfield
       CALL allocate_field_(layout(:,field_type), assigned_domain, field(:,i), field_type, data_type, dim3, dim4, name, ondevice)
    END DO
    !$OMP BARRIER
  END SUBROUTINE allocate_fields

  !==================================== PRIVATE : deallocate_field_ ===================================

  SUBROUTINE deallocate_field_(assigned_dom, field)
    USE omp_para
    IMPLICIT NONE
    LOGICAL, INTENT(IN) :: assigned_dom(:)
    TYPE(t_field) :: field(:)
    INTEGER :: data_type
    INTEGER :: ind
    DO ind=1,SIZE(field)
       IF (.NOT. assigned_dom(ind)  .OR. .NOT. is_omp_level_master) CYCLE

       data_type=field(ind)%data_type
       IF (data_type==type_real) THEN
          IF (field(ind)%ondevice) THEN
             !$acc exit data delete(field(ind)%rval4d(:,:,:))
             CONTINUE
          END IF
          DEALLOCATE(field(ind)%rval4d)
       END IF
       IF (data_type==type_integer) THEN
          IF (field(ind)%ondevice) THEN
             !$acc exit data delete(field(ind)%ival4d(:,:,:))
             CONTINUE
          END IF
          DEALLOCATE(field(ind)%ival4d)
       END IF
       IF (data_type==type_logical) THEN
          IF (field(ind)%ondevice) THEN
             !$acc exit data delete(field(ind)%lval4d(:,:,:))
             CONTINUE
          END IF
          DEALLOCATE(field(ind)%lval4d)
       END IF

    ENDDO
  END SUBROUTINE deallocate_field_

  !==================================== PUBLIC : deallocate_field ===================================

  SUBROUTINE deallocate_field_glo(field)
    USE domain_mod
    IMPLICIT NONE
    TYPE(t_field),POINTER :: field(:)
    INTEGER :: ind
    LOGICAL :: assigned_dom(ndomain_glo)
    ! ONLY the master thread is allowed to call this routine
    CALL deallocate_field_(assigned_dom, field)
    DEALLOCATE(field)
  END SUBROUTINE deallocate_field_glo

  SUBROUTINE deallocate_field(field)
    USE domain_mod
    USE omp_para
    IMPLICIT NONE
    TYPE(t_field),POINTER :: field(:)
    !$OMP BARRIER
    CALL deallocate_field_(assigned_domain, field)
    !$OMP BARRIER
    !$OMP MASTER
    DEALLOCATE(field)
    !$OMP END MASTER
    !$OMP BARRIER
  END SUBROUTINE deallocate_field

  SUBROUTINE deallocate_fields(field)
    USE domain_mod
    USE omp_para
    IMPLICIT NONE
    TYPE(t_field), POINTER :: field(:,:)
    INTEGER :: i
    !$OMP BARRIER
    DO i=1,SIZE(field,2)
       CALL deallocate_field_(assigned_domain, field(:,i))
    END DO
    !$OMP BARRIER
    !$OMP MASTER
    DEALLOCATE(field)
    !$OMP END MASTER
    !$OMP BARRIER
  END SUBROUTINE deallocate_fields

  !====================================== getval ===================================

  SUBROUTINE getval_r2d(field_pt,field)
    REAL(rstd), POINTER, INTENT(INOUT) :: field_pt(:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=2 .OR. field%data_type/=type_real) THEN
       PRINT *, 'getval_r2d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%rval2d
  END SUBROUTINE getval_r2d

  SUBROUTINE getval_r3d(field_pt,field)
    REAL(rstd), POINTER, INTENT(INOUT) :: field_pt(:,:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=3 .OR. field%data_type/=type_real) THEN
       PRINT *, 'getval_r3d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%rval3d
  END SUBROUTINE getval_r3d

  SUBROUTINE getval_r4d(field_pt,field)
    REAL(rstd), POINTER, INTENT(INOUT) :: field_pt(:,:,:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=4 .OR. field%data_type/=type_real) THEN
       PRINT *, 'getval_r4d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%rval4d
  END SUBROUTINE getval_r4d

  SUBROUTINE getval_i2d(field_pt,field)
    INTEGER, POINTER, INTENT(INOUT) :: field_pt(:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=2 .OR. field%data_type/=type_integer) THEN
       PRINT *, 'getval_i2d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%ival2d
  END SUBROUTINE getval_i2d

  SUBROUTINE getval_i3d(field_pt,field)
    INTEGER, POINTER, INTENT(INOUT) :: field_pt(:,:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=3 .OR. field%data_type/=type_integer) THEN
       PRINT *, 'getval_i3d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%ival3d
  END SUBROUTINE getval_i3d

  SUBROUTINE getval_i4d(field_pt,field)
    INTEGER, POINTER, INTENT(INOUT) :: field_pt(:,:,:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=4 .OR. field%data_type/=type_integer) THEN
       PRINT *, 'getval_i4d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%ival4d
  END SUBROUTINE getval_i4d

  SUBROUTINE getval_l2d(field_pt,field)
    LOGICAL, POINTER, INTENT(INOUT) :: field_pt(:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=2 .OR. field%data_type/=type_logical) THEN
       PRINT *, 'getval_l2d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%lval2d
  END SUBROUTINE getval_l2d

  SUBROUTINE getval_l3d(field_pt,field)
    LOGICAL, POINTER, INTENT(INOUT) :: field_pt(:,:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=3 .OR. field%data_type/=type_logical) THEN
       PRINT *, 'getval_l3d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%lval3d
  END SUBROUTINE getval_l3d

  SUBROUTINE getval_l4d(field_pt,field)
    LOGICAL, POINTER, INTENT(INOUT) :: field_pt(:,:,:)
    TYPE(t_field),INTENT(IN) :: field

    IF (field%ndim/=4 .OR. field%data_type/=type_logical) THEN
       PRINT *, 'getval_l4d : bad pointer assignment with ' // TRIM(field%name)
       STOP
    END IF
    field_pt=>field%lval4d
  END SUBROUTINE getval_l4d

  !===================== Data transfer between host (CPU) and device (GPU) =========================

  SUBROUTINE update_device_field(field)
    USE domain_mod
    USE omp_para
    TYPE(t_field) :: field(:)
    INTEGER :: ind

    DO ind=1,ndomain
       IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE

       IF (.NOT. field(ind)%ondevice) CALL create_device_field(field(ind))
       IF (field(ind)%data_type==type_real) THEN
          !$acc update device(field(ind)%rval4d(:,:,:)) async
          CONTINUE
       END IF
       IF (field(ind)%data_type==type_integer) THEN
          !$acc update device(field(ind)%ival4d(:,:,:)) async
          CONTINUE
       END IF
       IF (field(ind)%data_type==type_logical) THEN
          !$acc update device(field(ind)%lval4d(:,:,:)) async
          CONTINUE
       END IF

    ENDDO
    !$OMP BARRIER
  END SUBROUTINE update_device_field

  SUBROUTINE update_host_field(field)
    USE domain_mod
    USE omp_para
    TYPE(t_field) :: field(:)
    INTEGER :: ind

    DO ind=1,ndomain
       IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE

       IF (field(ind)%ondevice) THEN
          IF (field(ind)%data_type==type_real) THEN
             !$acc update host(field(ind)%rval4d(:,:,:)) async
             CONTINUE
          END IF
          IF (field(ind)%data_type==type_integer) THEN
             !$acc update host(field(ind)%ival4d(:,:,:)) async
             CONTINUE
          END IF
          IF (field(ind)%data_type==type_logical) THEN
             !$acc update host(field(ind)%lval4d(:,:,:)) async
             CONTINUE
          END IF

       END IF
    ENDDO
    !$acc wait
    !$OMP BARRIER
  END SUBROUTINE update_host_field

  SUBROUTINE create_device_field(field)
    TYPE(t_field) :: field

    IF (field%ondevice) THEN
       PRINT *, "Field is already on device !"
       STOP 1
    END IF
    IF (field%data_type==type_real) THEN
       !$acc enter data create(field%rval4d(:,:,:)) async
    END IF
    IF (field%data_type==type_integer) THEN
       !$acc enter data create(field%ival4d(:,:,:)) async
    END IF
    IF (field%data_type==type_logical) THEN
       !$acc enter data create(field%lval4d(:,:,:)) async
    END IF

    field%ondevice = .TRUE.
  END SUBROUTINE create_device_field

  !==================================== slice tools  ===================================

  SUBROUTINE extract_slice(field_in, field_out, l)
    USE domain_mod
    USE omp_para
    IMPLICIT NONE
    TYPE(t_field) :: field_in(:)
    TYPE(t_field) :: field_out(:)
    INTEGER, INTENT(IN) :: l

    INTEGER :: ind
    INTEGER :: data_type

    !$OMP BARRIER
    DO ind = 1, ndomain
       data_type = field_in(ind)%data_type
       IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE

       IF (field_in(ind)%ndim==3 .AND. field_out(ind)%ndim==2) THEN
          IF (data_type==type_integer)  field_out(ind)%ival2d=field_in(ind)%ival3d(:,l)
          IF (data_type==type_real)     field_out(ind)%rval2d=field_in(ind)%rval3d(:,l)
          IF (data_type==type_logical)  field_out(ind)%lval2d=field_in(ind)%lval3d(:,l)
       ELSE IF  (field_in(ind)%ndim==4 .AND. field_out(ind)%ndim==3) THEN
          IF (data_type==type_integer)  field_out(ind)%ival3d=field_in(ind)%ival4d(:,:,l)
          IF (data_type==type_real)     field_out(ind)%rval3d=field_in(ind)%rval4d(:,:,l)
          IF (data_type==type_logical)  field_out(ind)%lval3d=field_in(ind)%lval4d(:,:,l)
       ELSE
          PRINT *, 'extract_slice : cannot extract slice, dimension incompatible'
          STOP
       ENDIF
    ENDDO
    !$OMP BARRIER
  END  SUBROUTINE extract_slice


  SUBROUTINE insert_slice(field_in, field_out, l)
    USE domain_mod
    USE omp_para
    IMPLICIT NONE
    TYPE(t_field) :: field_in(:)
    TYPE(t_field) :: field_out(:)
    INTEGER,INTENT(IN) :: l

    INTEGER :: ind
    INTEGER :: data_type

    !$OMP BARRIER
    DO ind=1,ndomain
       data_type=field_in(ind)%data_type
       IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE

       IF (field_in(ind)%ndim==2 .AND. field_out(ind)%ndim==3) THEN
          IF (data_type==type_integer)  field_out(ind)%ival3d(:,l)=field_in(ind)%ival2d(:)
          IF (data_type==type_real)     field_out(ind)%rval3d(:,l)=field_in(ind)%rval2d(:)
          IF (data_type==type_logical)  field_out(ind)%lval3d(:,l)=field_in(ind)%lval2d(:)
       ELSE IF  (field_in(ind)%ndim==4 .AND. field_out(ind)%ndim==3) THEN
          IF (data_type==type_integer)  field_out(ind)%ival4d(:,:,l)=field_out(ind)%ival3d(:,:)
          IF (data_type==type_real)     field_out(ind)%rval4d(:,:,l)=field_out(ind)%rval3d(:,:)
          IF (data_type==type_logical)  field_out(ind)%lval4d(:,:,l)=field_out(ind)%lval3d(:,:)
       ELSE
          PRINT *, 'extract_slice : cannot insert slice, dimension incompatible'
          STOP
       ENDIF
    ENDDO
    !$OMP BARRIER

  END SUBROUTINE insert_slice

END MODULE field_mod
