Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • nemo/nemo
  • sparonuz/nemo
  • hatfield/nemo
  • extdevs/nemo
4 results
Show changes
Hide whitespace changes
Inline Side-by-side
Showing
with 685 additions and 584 deletions
src/OCE/LBC/lbcnfd.F90
View file @ 3efbcb17
...@@ -23,8 +23,11 @@ MODULE lbcnfd ...@@ -23,8 +23,11 @@ MODULE lbcnfd
PRIVATE PRIVATE
INTERFACE lbc_nfd ! called by mpp_nfd, lbc_lnk_pt2pt or lbc_lnk_neicoll INTERFACE lbc_nfd ! called by mpp_nfd, lbc_lnk_pt2pt or lbc_lnk_neicoll
MODULE PROCEDURE lbc_nfd_sp, lbc_nfd_ext_sp MODULE PROCEDURE lbc_nfd_sp, lbc_nfd_dp
MODULE PROCEDURE lbc_nfd_dp, lbc_nfd_ext_dp END INTERFACE
INTERFACE lbc_nfd_ext ! called by mpp_lnk_2d_icb
MODULE PROCEDURE lbc_nfd_ext_sp, lbc_nfd_ext_dp
END INTERFACE END INTERFACE
INTERFACE mpp_nfd ! called by lbc_lnk_pt2pt or lbc_lnk_neicoll INTERFACE mpp_nfd ! called by lbc_lnk_pt2pt or lbc_lnk_neicoll
...@@ -33,11 +36,13 @@ MODULE lbcnfd ...@@ -33,11 +36,13 @@ MODULE lbcnfd
PUBLIC mpp_nfd ! mpi north fold conditions PUBLIC mpp_nfd ! mpi north fold conditions
PUBLIC lbc_nfd ! north fold conditions PUBLIC lbc_nfd ! north fold conditions
PUBLIC lbc_nfd_ext ! north fold conditions, called by mpp_lnk_2d_icb
INTEGER, PUBLIC :: nfd_nbnei INTEGER, PUBLIC :: nfd_nbnei
INTEGER, PUBLIC, ALLOCATABLE, DIMENSION (: ) :: nfd_rknei INTEGER, PUBLIC, ALLOCATABLE, DIMENSION (: ) :: nfd_rknei
INTEGER, PUBLIC, ALLOCATABLE, DIMENSION (:,:) :: nfd_rksnd INTEGER, PUBLIC, ALLOCATABLE, DIMENSION (:,:,:) :: nfd_rksnd
INTEGER, PUBLIC, ALLOCATABLE, DIMENSION (:,:) :: nfd_jisnd INTEGER, PUBLIC, ALLOCATABLE, DIMENSION (:,:,:) :: nfd_jisnd
LOGICAL, PUBLIC, ALLOCATABLE, DIMENSION (:,: ) :: lnfd_same
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
......
src/OCE/LBC/lib_mpp.F90
View file @ 3efbcb17
...@@ -142,10 +142,10 @@ MODULE lib_mpp ...@@ -142,10 +142,10 @@ MODULE lib_mpp
INTEGER :: MPI_SUMDD INTEGER :: MPI_SUMDD
! Neighbourgs informations ! Neighbourgs informations
INTEGER, PARAMETER, PUBLIC :: n_hlsmax = 3 INTEGER, PARAMETER, PUBLIC :: n_hlsmax = 2
INTEGER, DIMENSION( 8), PUBLIC :: mpinei !: 8-neighbourg MPI indexes (starting at 0, -1 if no neighbourg) INTEGER, DIMENSION( 8), PUBLIC :: mpinei !: 8-neighbourg MPI indexes (starting at 0, -1 if no neighbourg)
INTEGER, DIMENSION(n_hlsmax,8), PUBLIC :: mpiSnei !: 8-neighbourg Send MPI indexes (starting at 0, -1 if no neighbourg) INTEGER, DIMENSION(0:n_hlsmax,8), PUBLIC :: mpiSnei !: 8-neighbourg Send MPI indexes (starting at 0, -1 if no neighbourg)
INTEGER, DIMENSION(n_hlsmax,8), PUBLIC :: mpiRnei !: 8-neighbourg Recv MPI indexes (starting at 0, -1 if no neighbourg) INTEGER, DIMENSION(0:n_hlsmax,8), PUBLIC :: mpiRnei !: 8-neighbourg Recv MPI indexes (starting at 0, -1 if no neighbourg)
INTEGER, PARAMETER, PUBLIC :: jpwe = 1 !: WEst INTEGER, PARAMETER, PUBLIC :: jpwe = 1 !: WEst
INTEGER, PARAMETER, PUBLIC :: jpea = 2 !: EAst INTEGER, PARAMETER, PUBLIC :: jpea = 2 !: EAst
INTEGER, PARAMETER, PUBLIC :: jpso = 3 !: SOuth INTEGER, PARAMETER, PUBLIC :: jpso = 3 !: SOuth
...@@ -1127,7 +1127,7 @@ CONTAINS ...@@ -1127,7 +1127,7 @@ CONTAINS
INTEGER :: ierr INTEGER :: ierr
LOGICAL, PARAMETER :: ireord = .FALSE. LOGICAL, PARAMETER :: ireord = .FALSE.
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
#if ! defined key_mpi_off && ! defined key_mpi2 #if ! defined key_mpi_off
iScnt4 = COUNT( mpiSnei(khls,1:4) >= 0 ) iScnt4 = COUNT( mpiSnei(khls,1:4) >= 0 )
iRcnt4 = COUNT( mpiRnei(khls,1:4) >= 0 ) iRcnt4 = COUNT( mpiRnei(khls,1:4) >= 0 )
...@@ -1141,10 +1141,19 @@ CONTAINS ...@@ -1141,10 +1141,19 @@ CONTAINS
iSnei8 = PACK( mpiSnei(khls,1:8), mask = mpiSnei(khls,1:8) >= 0 ) iSnei8 = PACK( mpiSnei(khls,1:8), mask = mpiSnei(khls,1:8) >= 0 )
iRnei8 = PACK( mpiRnei(khls,1:8), mask = mpiRnei(khls,1:8) >= 0 ) iRnei8 = PACK( mpiRnei(khls,1:8), mask = mpiRnei(khls,1:8) >= 0 )
! Isolated processes (i.e., processes WITH no outgoing or incoming edges, that is, processes that have specied
! indegree and outdegree as zero and thus DO not occur as source or destination rank in the graph specication)
! are allowed.
# if ! defined key_mpi2
CALL MPI_Dist_graph_create_adjacent( mpi_comm_oce, iScnt4, iSnei4, MPI_UNWEIGHTED, iRcnt4, iRnei4, MPI_UNWEIGHTED, & CALL MPI_Dist_graph_create_adjacent( mpi_comm_oce, iScnt4, iSnei4, MPI_UNWEIGHTED, iRcnt4, iRnei4, MPI_UNWEIGHTED, &
& MPI_INFO_NULL, ireord, mpi_nc_com4(khls), ierr ) & MPI_INFO_NULL, ireord, mpi_nc_com4(khls), ierr )
CALL MPI_Dist_graph_create_adjacent( mpi_comm_oce, iScnt8, iSnei8, MPI_UNWEIGHTED, iRcnt8, iRnei8, MPI_UNWEIGHTED, & CALL MPI_Dist_graph_create_adjacent( mpi_comm_oce, iScnt8, iSnei8, MPI_UNWEIGHTED, iRcnt8, iRnei8, MPI_UNWEIGHTED, &
& MPI_INFO_NULL, ireord, mpi_nc_com8(khls), ierr) & MPI_INFO_NULL, ireord, mpi_nc_com8(khls), ierr)
# else
mpi_nc_com4(khls) = -1
mpi_nc_com8(khls) = -1
# endif
DEALLOCATE( iSnei4, iRnei4, iSnei8, iRnei8 ) DEALLOCATE( iSnei4, iRnei4, iSnei8, iRnei8 )
#endif #endif
...@@ -1307,7 +1316,7 @@ CONTAINS ...@@ -1307,7 +1316,7 @@ CONTAINS
IF ( ncomm_sequence(ji,1) .GT. 1 .AND. ncomm_sequence(ji,2) .GT. 1 ) jj = jj + 1 IF ( ncomm_sequence(ji,1) .GT. 1 .AND. ncomm_sequence(ji,2) .GT. 1 ) jj = jj + 1
jh = MAX (jh, ncomm_sequence(ji,1)*ncomm_sequence(ji,2)) jh = MAX (jh, ncomm_sequence(ji,1)*ncomm_sequence(ji,2))
END DO END DO
WRITE(numcom,'(A,I3)') ' 3D Exchanged halos : ', jk WRITE(numcom,'(A,I3)') ' 3D or 4D Exchanged halos : ', jk
WRITE(numcom,'(A,I3)') ' Multi arrays exchanged halos : ', jf WRITE(numcom,'(A,I3)') ' Multi arrays exchanged halos : ', jf
WRITE(numcom,'(A,I3)') ' from which 3D : ', jj WRITE(numcom,'(A,I3)') ' from which 3D : ', jj
WRITE(numcom,'(A,I10)') ' Array max size : ', jh*jpi*jpj WRITE(numcom,'(A,I10)') ' Array max size : ', jh*jpi*jpj
......
src/OCE/LBC/mpp_lbc_north_icb_generic.h90
View file @ 3efbcb17
...@@ -92,7 +92,7 @@ ...@@ -92,7 +92,7 @@
! 2. North-Fold boundary conditions ! 2. North-Fold boundary conditions
! ---------------------------------- ! ----------------------------------
CALL lbc_nfd( ztab_e(:,1-kextj:ipj+kextj), cd_type, psgn, kextj ) CALL lbc_nfd_ext( ztab_e(:,1-kextj:ipj+kextj), cd_type, psgn, kextj )
ij = 1 - kextj ij = 1 - kextj
!! Scatter back to pt2d !! Scatter back to pt2d
......
src/OCE/LBC/mpp_lnk_icb_generic.h90
View file @ 3efbcb17
...@@ -87,7 +87,7 @@ ...@@ -87,7 +87,7 @@
IF( l_IdoNFold ) THEN IF( l_IdoNFold ) THEN
! !
SELECT CASE ( jpni ) SELECT CASE ( jpni )
CASE ( 1 ) ; CALL lbc_nfd ( pt2d(1:jpi,1:jpj+kextj), cd_type, psgn, kextj ) CASE ( 1 ) ; CALL lbc_nfd_ext ( pt2d(1:jpi,1:jpj+kextj), cd_type, psgn, kextj )
CASE DEFAULT ; CALL LBCNORTH ( pt2d(1:jpi,1:jpj+kextj), cd_type, psgn, kextj ) CASE DEFAULT ; CALL LBCNORTH ( pt2d(1:jpi,1:jpj+kextj), cd_type, psgn, kextj )
END SELECT END SELECT
! !
......
src/OCE/LBC/mpp_loc_generic.h90
View file @ 3efbcb17
...@@ -47,6 +47,7 @@ ...@@ -47,6 +47,7 @@
! !
INTEGER :: ierror, ii, idim INTEGER :: ierror, ii, idim
INTEGER :: index0 INTEGER :: index0
INTEGER :: ihls, ipiglo, ipjglo
INTEGER , DIMENSION(:), ALLOCATABLE :: ilocs INTEGER , DIMENSION(:), ALLOCATABLE :: ilocs
REAL(PRECISION) :: zmin ! local minimum REAL(PRECISION) :: zmin ! local minimum
REAL(PRECISION), DIMENSION(2,1) :: zain, zaout REAL(PRECISION), DIMENSION(2,1) :: zain, zaout
...@@ -60,6 +61,9 @@ ...@@ -60,6 +61,9 @@
ENDIF ENDIF
! !
idim = SIZE(kindex) idim = SIZE(kindex)
ihls = ( SIZE(ARRAY_IN(:,:,:), 1) - Ni_0 ) / 2
ipiglo = Ni0glo + 2*ihls
ipjglo = Nj0glo + 2*ihls
! !
IF ( ANY( MASK_IN(:,:,:) ) ) THEN ! there is at least 1 valid point... IF ( ANY( MASK_IN(:,:,:) ) ) THEN ! there is at least 1 valid point...
! !
...@@ -68,9 +72,9 @@ ...@@ -68,9 +72,9 @@
ilocs = LOC_OPERATION( ARRAY_IN(:,:,:) , mask= MASK_IN(:,:,:) ) ilocs = LOC_OPERATION( ARRAY_IN(:,:,:) , mask= MASK_IN(:,:,:) )
zmin = ARRAY_IN(ilocs(1),ilocs(2),ilocs(3)) zmin = ARRAY_IN(ilocs(1),ilocs(2),ilocs(3))
! !
kindex(1) = mig( ilocs(1) ) kindex(1) = mig( ilocs(1), ihls )
#if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */ #if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */
kindex(2) = mjg( ilocs(2) ) kindex(2) = mjg( ilocs(2), ihls )
#endif #endif
#if defined DIM_3d /* avoid warning when kindex has 2 elements */ #if defined DIM_3d /* avoid warning when kindex has 2 elements */
kindex(3) = ilocs(3) kindex(3) = ilocs(3)
...@@ -80,10 +84,10 @@ ...@@ -80,10 +84,10 @@
! !
index0 = kindex(1)-1 ! 1d index starting at 0 index0 = kindex(1)-1 ! 1d index starting at 0
#if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */ #if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */
index0 = index0 + jpiglo * (kindex(2)-1) index0 = index0 + ipiglo * (kindex(2)-1)
#endif #endif
#if defined DIM_3d /* avoid warning when kindex has 2 elements */ #if defined DIM_3d /* avoid warning when kindex has 2 elements */
index0 = index0 + jpiglo * jpjglo * (kindex(3)-1) index0 = index0 + ipiglo * ipjglo * (kindex(3)-1)
#endif #endif
ELSE ELSE
! special case for land processors ! special case for land processors
...@@ -105,20 +109,20 @@ ...@@ -105,20 +109,20 @@
pmin = zaout(1,1) pmin = zaout(1,1)
index0 = NINT( zaout(2,1) ) index0 = NINT( zaout(2,1) )
#if defined DIM_3d /* avoid warning when kindex has 2 elements */ #if defined DIM_3d /* avoid warning when kindex has 2 elements */
kindex(3) = index0 / (jpiglo*jpjglo) kindex(3) = index0 / (ipiglo*ipjglo)
index0 = index0 - kindex(3) * (jpiglo*jpjglo) index0 = index0 - kindex(3) * (ipiglo*ipjglo)
#endif #endif
#if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */ #if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */
kindex(2) = index0 / jpiglo kindex(2) = index0 / ipiglo
index0 = index0 - kindex(2) * jpiglo index0 = index0 - kindex(2) * ipiglo
#endif #endif
kindex(1) = index0 kindex(1) = index0
kindex(:) = kindex(:) + 1 ! start indices at 1 kindex(:) = kindex(:) + 1 ! start indices at 1
IF( .NOT. llhalo ) THEN IF( .NOT. llhalo ) THEN
kindex(1) = kindex(1) - nn_hls kindex(1) = kindex(1) - ihls
#if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */ #if defined DIM_2d || defined DIM_3d /* avoid warning when kindex has 1 element */
kindex(2) = kindex(2) - nn_hls kindex(2) = kindex(2) - ihls
#endif #endif
ENDIF ENDIF
......
src/OCE/LBC/mpp_nfd_generic.h90
View file @ 3efbcb17
SUBROUTINE mpp_nfd_/**/PRECISION( ptab, cd_nat, psgn, kfillmode, pfillval, khls, kfld ) SUBROUTINE mpp_nfd_/**/PRECISION( ptab, cd_nat, psgn, kfillmode, pfillval, kfld )
TYPE(PTR_4d_/**/PRECISION), DIMENSION(:), INTENT(inout) :: ptab ! pointer of arrays on which apply the b.c. TYPE(PTR_4d_/**/PRECISION), DIMENSION(:), INTENT(inout) :: ptab ! pointer of arrays on which apply the b.c.
CHARACTER(len=1), DIMENSION(:), INTENT(in ) :: cd_nat ! nature of array grid-points CHARACTER(len=1), DIMENSION(kfld), INTENT(in ) :: cd_nat ! nature of array grid-points
REAL(PRECISION), DIMENSION(:), INTENT(in ) :: psgn ! sign used across the north fold boundary REAL(PRECISION), DIMENSION(kfld), INTENT(in ) :: psgn ! sign used across the north fold boundary
INTEGER , INTENT(in ) :: kfillmode ! filling method for halo over land INTEGER , INTENT(in ) :: kfillmode ! filling method for halo over land
REAL(PRECISION) , INTENT(in ) :: pfillval ! background value (used at closed boundaries) REAL(PRECISION) , INTENT(in ) :: pfillval ! background value (used at closed boundaries)
INTEGER , INTENT(in ) :: khls ! halo size, default = nn_hls INTEGER , INTENT(in ) :: kfld ! number of pt3d arrays
INTEGER , INTENT(in ) :: kfld ! number of pt3d arrays
! !
LOGICAL :: ll_add_line LOGICAL :: ll_add_line
INTEGER :: ji, jj, jk, jl, jf, jr, jg, jn ! dummy loop indices INTEGER :: ji, jj, jk, jl, jf, jr, jg, jn ! dummy loop indices
INTEGER :: ipi, ipj, ipj2, ipk, ipl, ipf ! dimension of the input array INTEGER :: ierr, ibuffsize, impp, ipi0
INTEGER :: ierr, ibuffsize, iis0, iie0, impp INTEGER :: ii1, ii2, ij1, ij2, ij3, iig, inei
INTEGER :: ii1, ii2, ij1, ij2, iis, iie, iib, iig, iin INTEGER :: i0max, ilntot, iisht, ijsht, ihsz
INTEGER :: i0max INTEGER :: iproc, ijnr, ipjtot, iFT, iFU, i012
INTEGER :: ij, iproc, ipni, ijnr INTEGER, DIMENSION(kfld) :: ipi, ipj, ipj1, ipj2, ipk, ipl ! dimension of the input array
INTEGER, DIMENSION (:), ALLOCATABLE :: ireq_s, ireq_r ! for mpi_isend when avoiding mpi_allgather INTEGER, DIMENSION(kfld) :: ihls ! halo size
INTEGER :: ipjtot ! sum of lines for all multi fields INTEGER, DIMENSION(:) , ALLOCATABLE :: ireq_s, ireq_r ! for mpi_isend when avoiding mpi_allgather
INTEGER :: i012 ! 0, 1 or 2 INTEGER, DIMENSION(:) , ALLOCATABLE :: ipjfld ! number of sent lines for each field
INTEGER , DIMENSION(:,:) , ALLOCATABLE :: ijsnd ! j-position of sent lines for each field REAL(PRECISION) :: zhuge, zztmp
INTEGER , DIMENSION(:,:) , ALLOCATABLE :: ijbuf ! j-position of send buffer lines for each field REAL(PRECISION), DIMENSION(:,:) , ALLOCATABLE :: zbufs ! buffer, receive and work arrays
INTEGER , DIMENSION(:,:) , ALLOCATABLE :: ijrcv ! j-position of recv buffer lines for each field REAL(PRECISION), DIMENSION(:,:,:), ALLOCATABLE :: zbufr ! buffer, receive and work arrays
INTEGER , DIMENSION(:,:) , ALLOCATABLE :: ii1st, iiend REAL(PRECISION), DIMENSION(:,:) , ALLOCATABLE :: znorthloc
INTEGER , DIMENSION(:) , ALLOCATABLE :: ipjfld ! number of sent lines for each field REAL(PRECISION), DIMENSION(:,:,:), ALLOCATABLE :: znorthall
REAL(PRECISION), DIMENSION(:,:,:,:) , ALLOCATABLE :: zbufs ! buffer, receive and work arrays TYPE(PTR_4D_/**/PRECISION), DIMENSION(1) :: ztabglo ! array or pointer of arrays on which apply the b.c.
REAL(PRECISION), DIMENSION(:,:,:,:,:) , ALLOCATABLE :: zbufr ! buffer, receive and work arrays
REAL(PRECISION), DIMENSION(:,:,:,:,:) , ALLOCATABLE :: znorthloc
REAL(PRECISION), DIMENSION(:,:,:,:,:,:), ALLOCATABLE :: znorthglo
TYPE(PTR_4D_/**/PRECISION), DIMENSION(:), ALLOCATABLE :: ztabglo ! array or pointer of arrays on which apply the b.c.
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
ipk = SIZE(ptab(1)%pt4d,3) zhuge = HUGE(0._/**/PRECISION) ! avoid to call the huge function inside do loops
ipl = SIZE(ptab(1)%pt4d,4) !
ipf = kfld DO jf = 1, kfld
ipi(jf) = SIZE(ptab(jf)%pt4d,1)
ipj(jf) = SIZE(ptab(jf)%pt4d,2)
ipk(jf) = SIZE(ptab(jf)%pt4d,3)
ipl(jf) = SIZE(ptab(jf)%pt4d,4)
ihls(jf) = ( ipi(jf) - Ni_0 ) / 2
END DO
! !
IF( ln_nnogather ) THEN !== no allgather exchanges ==! IF( ln_nnogather ) THEN !== no allgather exchanges ==!
...@@ -61,74 +62,43 @@ ...@@ -61,74 +62,43 @@
! also force it if not restart during the first 2 steps (leap frog?) ! also force it if not restart during the first 2 steps (leap frog?)
ll_add_line = l_full_nf_update .OR. ( ncom_stp <= nit000+1 .AND. .NOT. ln_rstart ) ll_add_line = l_full_nf_update .OR. ( ncom_stp <= nit000+1 .AND. .NOT. ln_rstart )
ALLOCATE(ipjfld(ipf)) ! how many lines do we exchange for each field? ALLOCATE(ipjfld(kfld)) ! how many lines do we send for each field?
IF( ll_add_line ) THEN IF( ll_add_line ) THEN
DO jf = 1, ipf ! Loop over the number of arrays to be processed DO jf = 1, kfld ! Loop over the number of arrays to be processed
ipjfld(jf) = khls + COUNT( (/ c_NFtype == 'T' .OR. cd_nat(jf) == 'V' .OR. cd_nat(jf) == 'F' /) ) ipjfld(jf) = ihls(jf) + COUNT( (/ c_NFtype == 'T' .OR. cd_nat(jf) == 'V' .OR. cd_nat(jf) == 'F' /) ) &
& + COUNT( (/ c_NFtype == 'F' .AND. cd_nat(jf) == 'T' .AND. ihls(jf) == 0 /) )
END DO END DO
ELSE ELSE
ipjfld(:) = khls ipjfld(:) = ihls(:)
ENDIF ENDIF
ipj = MAXVAL(ipjfld(:)) ! Max 2nd dimension of message transfers
ipjtot = SUM( ipjfld(:)) ! Total number of lines to be exchanged
! Index of modifying lines in input
ALLOCATE( ijsnd(ipj, ipf), ijbuf(ipj, ipf), ijrcv(ipj, ipf), ii1st(ipj, ipf), iiend(ipj, ipf) )
ij1 = 0
DO jf = 1, ipf ! Loop over the number of arrays to be processed
!
DO jj = 1, khls ! first khls lines (starting from top) must be fully defined
ii1st(jj, jf) = 1
iiend(jj, jf) = jpi
END DO
!
! what do we do with line khls+1 (starting from top)
IF( c_NFtype == 'T' ) THEN ! * North fold T-point pivot
SELECT CASE ( cd_nat(jf) )
CASE ('T','W') ; i012 = 1 ; ii1st(khls+1, jf) = mi0(jpiglo/2+2) ; iiend(khls+1, jf) = mi1(jpiglo-khls)
CASE ('U' ) ; i012 = 1 ; ii1st(khls+1, jf) = mi0(jpiglo/2+1) ; iiend(khls+1, jf) = mi1(jpiglo-khls)
CASE ('V' ) ; i012 = 2 ; ii1st(khls+1, jf) = 1 ; iiend(khls+1, jf) = jpi
CASE ('F' ) ; i012 = 2 ; ii1st(khls+1, jf) = 1 ; iiend(khls+1, jf) = jpi
END SELECT
ENDIF
IF( c_NFtype == 'F' ) THEN ! * North fold F-point pivot
SELECT CASE ( cd_nat(jf) )
CASE ('T','W') ; i012 = 0 ! we don't touch line khls+1
CASE ('U' ) ; i012 = 0 ! we don't touch line khls+1
CASE ('V' ) ; i012 = 1 ; ii1st(khls+1, jf) = mi0(jpiglo/2+1) ; iiend(khls+1, jf) = mi1(jpiglo-khls )
CASE ('F' ) ; i012 = 1 ; ii1st(khls+1, jf) = mi0(jpiglo/2+1) ; iiend(khls+1, jf) = mi1(jpiglo-khls-1)
END SELECT
ENDIF
!
DO jj = 1, ipjfld(jf)
ij1 = ij1 + 1
ijsnd(jj,jf) = jpj - 2*khls + jj - i012 ! sent lines (from bottom of sent lines)
ijbuf(jj,jf) = ij1 ! gather all lines in the snd/rcv buffers
ijrcv(jj,jf) = jpj - jj + 1 ! recv lines (from the top -> reverse order for jj)
END DO
!
END DO
! !
i0max = jpimax - 2 * khls ! we are not sending the halos i0max = MAXVAL( nfni_0, mask = nfproc /= -1 ) ! largest value of Ni_0 among processors (we are not sending halos)
ALLOCATE( zbufs(i0max,ipjtot,ipk,ipl), ireq_s(nfd_nbnei) ) ! store all the data to be sent in a buffer array ilntot = SUM( ipjfld(:) * ipk(:) * ipl(:) )
ibuffsize = i0max * ipjtot * ipk * ipl ALLOCATE( zbufs(i0max,ilntot), ireq_s(nfd_nbnei) ) ! store all the data to be sent in a buffer array
ibuffsize = i0max * ilntot ! must be the same for all processors -> use i0max
! !
! fill the send buffer with all the lines ! fill the send buffer with all the lines
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk ij1 = 0
DO jj = 1, ipjfld(jf) DO jf = 1, kfld
ij1 = ijbuf(jj,jf) !
ij2 = ijsnd(jj,jf) i012 = COUNT( (/ c_NFtype == 'T' /) ) + COUNT( (/ cd_nat(jf) == 'V' .OR. cd_nat(jf) == 'F' /) ) &
DO ji = Nis0, Nie0 ! should not use any other value & + COUNT( (/ ihls(jf) == 0 .AND. cd_nat(jf) == 'T' .AND. c_NFtype == 'F' /) ) ! 0, 1 OR 2
iib = ji - Nis0 + 1 ijsht = ipj(jf) - 2*ihls(jf) - i012 ! j-position of the sent lines (from bottom of sent lines)
zbufs(iib,ij1,jk,jl) = ptab(jf)%pt4d(ji,ij2,jk,jl) !
END DO DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf)
DO ji = Ni_0+1, i0max ! avoid sending uninitialized values (make sure we don't use it) DO jj = 1, ipjfld(jf)
zbufs(ji,ij1,jk,jl) = HUGE(0._/**/PRECISION) ! make sure we don't use it... ij1 = ij1 + 1
ij2 = jj + ijsht
DO ji = 1, Ni_0 ! use only inner domain
ii2 = ji + ihls(jf)
zbufs(ji,ij1) = ptab(jf)%pt4d(ii2,ij2,jk,jl)
END DO
DO ji = Ni_0+1, i0max ! avoid sending uninitialized values and make sure we don't use it
zbufs(ji,ij1) = zhuge
END DO
END DO END DO
END DO END DO ; END DO
END DO ; END DO ; END DO END DO ! jf
! !
! start waiting time measurement ! start waiting time measurement
IF( ln_timing ) CALL tic_tac(.TRUE.) IF( ln_timing ) CALL tic_tac(.TRUE.)
...@@ -136,68 +106,62 @@ ...@@ -136,68 +106,62 @@
! send the same buffer data to all neighbourgs as soon as possible ! send the same buffer data to all neighbourgs as soon as possible
DO jn = 1, nfd_nbnei DO jn = 1, nfd_nbnei
iproc = nfd_rknei(jn) iproc = nfd_rknei(jn)
IF( iproc /= narea-1 .AND. iproc /= -1 ) THEN IF( iproc /= narea-1 .AND. iproc /= -1 ) THEN ! it is neither me nor a land-only neighbourg
#if ! defined key_mpi_off #if ! defined key_mpi_off
CALL MPI_Isend( zbufs, ibuffsize, MPI_TYPE, iproc, 5, mpi_comm_oce, ireq_s(jn), ierr ) CALL MPI_Isend( zbufs, ibuffsize, MPI_TYPE, iproc, 5, mpi_comm_oce, ireq_s(jn), ierr )
#endif #endif
ELSE ELSE
ireq_s(jn) = MPI_REQUEST_NULL ireq_s(jn) = MPI_REQUEST_NULL ! must be defined for mpi_waitall
ENDIF ENDIF
END DO END DO
! !
ALLOCATE( zbufr(i0max,ipjtot,ipk,ipl,nfd_nbnei), ireq_r(nfd_nbnei) ) ALLOCATE( zbufr(i0max,ilntot,nfd_nbnei), ireq_r(nfd_nbnei) )
! !
DO jn = 1, nfd_nbnei DO jn = 1, nfd_nbnei ! 1st loop: first get data which does not need any communication
! ! ! -> this gives more time to receive the communications
iproc = nfd_rknei(jn) iproc = nfd_rknei(jn)
! !
IF( iproc == -1 ) THEN ! No neighbour (land proc that was suppressed) IF( iproc == -1 ) THEN ! No neighbour (land-only neighbourg that was suppressed)
! !
ireq_r(jn) = MPI_REQUEST_NULL ! no message to be received ireq_r(jn) = MPI_REQUEST_NULL ! no message to be received, must be defined for mpi_waitall
zbufr(:,:,:,:,jn) = HUGE(0._/**/PRECISION) ! default: define it and make sure we don't use it...
SELECT CASE ( kfillmode ) SELECT CASE ( kfillmode )
CASE ( jpfillnothing ) ! no filling CASE ( jpfillnothing ) ! no filling
CASE ( jpfillcopy ) ! filling with inner domain values CASE ( jpfillcopy ) ! filling with my inner domain values
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk ! ! trick: we use only the 1st value, see init_nfdcom
DO jj = 1, ipjfld(jf) zbufr(1,:,jn) = zbufs(1,:) ! chose to take the 1st inner domain point
ij1 = ijbuf(jj,jf)
ij2 = ijsnd(jj,jf) ! we will use only the first value, see init_nfdcom
zbufr(1,ij1,jk,jl,jn) = ptab(jf)%pt4d(Nis0,ij2,jk,jl) ! chose to take the 1st inner domain point
END DO
END DO ; END DO ; END DO
CASE ( jpfillcst ) ! filling with constant value CASE ( jpfillcst ) ! filling with constant value
zbufr(1,:,:,:,jn) = pfillval ! we will use only the first value, see init_nfdcom zbufr(1,:,jn) = pfillval ! trick: we use only the 1st value, see init_nfdcom
END SELECT END SELECT
! !
ELSE IF( iproc == narea-1 ) THEN ! get data from myself! ELSE IF( iproc == narea-1 ) THEN ! I get data from myself!
! !
ireq_r(jn) = MPI_REQUEST_NULL ! no message to be received ireq_r(jn) = MPI_REQUEST_NULL ! no message to be received, must be defined for mpi_waitall
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk zbufr(:,:,jn) = zbufs(:,:) ! we can directly do: received buffer = sent buffer!
DO jj = 1, ipjfld(jf)
ij1 = ijbuf(jj,jf)
ij2 = ijsnd(jj,jf)
DO ji = Nis0, Nie0 ! should not use any other value
iib = ji - Nis0 + 1
zbufr(iib,ij1,jk,jl,jn) = ptab(jf)%pt4d(ji,ij2,jk,jl)
END DO
END DO
END DO ; END DO ; END DO
! !
ELSE ! get data from a neighbour trough communication ENDIF
!
END DO ! nfd_nbnei
!
DO jn = 1, nfd_nbnei ! 2nd loop: now get data from a neighbour trough communication
!
iproc = nfd_rknei(jn)
IF( iproc /= narea-1 .AND. iproc /= -1 ) THEN ! it is neither me nor a land-only neighbourg
#if ! defined key_mpi_off #if ! defined key_mpi_off
CALL MPI_Irecv( zbufr(:,:,:,:,jn), ibuffsize, MPI_TYPE, iproc, 5, mpi_comm_oce, ireq_r(jn), ierr ) CALL MPI_Irecv( zbufr(:,:,jn), ibuffsize, MPI_TYPE, iproc, 5, mpi_comm_oce, ireq_r(jn), ierr )
#endif #endif
ENDIF ENDIF
!
END DO ! nfd_nbnei END DO ! nfd_nbnei
! !
#if ! defined key_mpi_off
CALL mpi_waitall(nfd_nbnei, ireq_r, MPI_STATUSES_IGNORE, ierr) ! wait for all Irecv CALL mpi_waitall(nfd_nbnei, ireq_r, MPI_STATUSES_IGNORE, ierr) ! wait for all Irecv
#endif
! !
IF( ln_timing ) CALL tic_tac(.FALSE.) IF( ln_timing ) CALL tic_tac(.FALSE.)
! !
! North fold boundary condition ! Apply the North pole folding
! !
DO jf = 1, ipf ij2 = 0
DO jf = 1, kfld
! !
SELECT CASE ( cd_nat(jf) ) ! which grid number? SELECT CASE ( cd_nat(jf) ) ! which grid number?
CASE ('T','W') ; iig = 1 ! T-, W-point CASE ('T','W') ; iig = 1 ! T-, W-point
...@@ -206,76 +170,43 @@ ...@@ -206,76 +170,43 @@
CASE ('F') ; iig = 4 ! F-point CASE ('F') ; iig = 4 ! F-point
END SELECT END SELECT
! !
DO jl = 1, ipl ; DO jk = 1, ipk ihsz = ihls(jf) ! shorter name
! iisht = nn_hls - ihsz
! if T point with F-point pivot : must be done first iFT = COUNT( (/ ihsz > 0 .AND. c_NFtype == 'F' .AND. cd_nat(jf) == 'T' /) ) ! F-folding and T grid
! --> specific correction of 3 points near the 2 pivots (to be clean, usually masked -> so useless) !
IF( c_NFtype == 'F' .AND. iig == 1 ) THEN DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf)
ij1 = jpj - khls ! j-index in the receiving array
ij2 = 1 ! only 1 line in the buffer
DO ji = mi0(khls), mi1(khls) ! change because of EW periodicity as we also change jpiglo-khls
iib = nfd_jisnd(mi0( khls),iig) ! i-index in the buffer
iin = nfd_rksnd(mi0( khls),iig) ! neigbhour-index in the buffer
IF( nfd_rknei(iin) == -1 .AND. kfillmode == jpfillnothing ) CYCLE
ptab(jf)%pt4d(ji,ij1,jk,jl) = zbufr(iib,ij2,jk,jl,iin) ! no psgn(jf)
END DO
DO ji = mi0(jpiglo/2+1), mi1(jpiglo/2+1)
iib = nfd_jisnd(mi0( jpiglo/2+1),iig) ! i-index in the buffer
iin = nfd_rksnd(mi0( jpiglo/2+1),iig) ! neigbhour-index in the buffer
IF( nfd_rknei(iin) == -1 .AND. kfillmode == jpfillnothing ) CYCLE
ptab(jf)%pt4d(ji,ij1,jk,jl) = zbufr(iib,ij2,jk,jl,iin) ! no psgn(jf)
END DO
DO ji = mi0(jpiglo-khls), mi1(jpiglo-khls)
iib = nfd_jisnd(mi0(jpiglo-khls),iig) ! i-index in the buffer
iin = nfd_rksnd(mi0(jpiglo-khls),iig) ! neigbhour-index in the buffer
IF( nfd_rknei(iin) == -1 .AND. kfillmode == jpfillnothing ) CYCLE
ptab(jf)%pt4d(ji,ij1,jk,jl) = zbufr(iib,ij2,jk,jl,iin) ! no psgn(jf)
END DO
ENDIF
! !
! Apply the North pole folding. DO jj = 1,ihsz ! NP folding for the last ihls(jf) lines of this field
DO jj = 1, ipjfld(jf) ! for all lines to be exchanged for this field ij1 = ipj(jf) - jj + 1 ! j-index in the receiving array (from the top -> reverse order for jj)
ij1 = ijrcv(jj,jf) ! j-index in the receiving array ij2 = ij2 + 1
ij2 = ijbuf(jj,jf) ! j-index in the buffer ij3 = ihsz+1 - jj + 1
iis = ii1st(jj,jf) ! stating i-index in the receiving array DO ji = 1, ipi(jf)
iie = iiend(jj,jf) ! ending i-index in the receiving array ii1 = ji + iisht
DO ji = iis, iie inei = nfd_rksnd(ii1,ij3,iig) ! neigbhour-index in the buffer
iib = nfd_jisnd(ji,iig) ! i-index in the buffer IF( nfd_rknei(inei) == -1 .AND. kfillmode == jpfillnothing ) CYCLE ! no neighbourg and do nothing to fill
iin = nfd_rksnd(ji,iig) ! neigbhour-index in the buffer ii2 = nfd_jisnd(ii1,ij3,iig) ! i-index in the buffer, starts at 1 in the inner domain
IF( nfd_rknei(iin) == -1 .AND. kfillmode == jpfillnothing ) CYCLE ptab(jf)%pt4d(ji,ij1,jk,jl) = psgn(jf) * zbufr(ii2,ij2,inei)
ptab(jf)%pt4d(ji,ij1,jk,jl) = psgn(jf) * zbufr(iib,ij2,jk,jl,iin)
END DO END DO
END DO END DO
! DO jj = ihsz+1, ipjfld(jf)+iFT ! NP folding for line ipj-ihsz that can be partially modified
! re-apply periodocity when we modified the eastern side of the inner domain (and not the full line) ij1 = ipj(jf) - jj + 1 ! j-index in the receiving array (from the top -> reverse order for jj)
IF( c_NFtype == 'T' ) THEN ! * North fold T-point pivot ij2 = ij2 + 1 - iFT
IF( iig <= 2 ) THEN ; iis = mi0(1) ; iie = mi1(khls) ! 'T','W','U': update west halo ij3 = 1
ELSE ; iis = 1 ; iie = 0 ! 'V','F' : full line already exchanged DO ji = 1, ipi(jf)
ENDIF ii1 = ji + iisht
ENDIF IF( lnfd_same(ii1,iig) ) CYCLE ! do nothing if should not be modified
IF( c_NFtype == 'F' ) THEN ! * North fold F-point pivot inei = nfd_rksnd(ii1,ij3,iig) ! neigbhour-index in the buffer
IF( iig <= 2 ) THEN ; iis = 1 ; iie = 0 ! 'T','W','U': nothing to do IF( nfd_rknei(inei) == -1 .AND. kfillmode == jpfillnothing ) CYCLE ! no neighbourg and do nothing to fill
ELSEIF( iig == 3 ) THEN ; iis = mi0(1) ; iie = mi1(khls) ! 'V' : update west halo ii2 = nfd_jisnd(ii1,ij3,iig) ! i-index in the buffer, starts at 1 in the inner domain
ELSEIF( khls > 1 ) THEN ; iis = mi0(1) ; iie = mi1(khls-1) ! 'F' and khls > 1 ptab(jf)%pt4d(ji,ij1,jk,jl) = psgn(jf) * zbufr(ii2,ij2,inei)
ELSE ; iis = 1 ; iie = 0 ! 'F' and khls == 1 : nothing to do END DO
ENDIF
ENDIF
jj = ipjfld(jf) ! only for the last line of this field
ij1 = ijrcv(jj,jf) ! j-index in the receiving array
ij2 = ijbuf(jj,jf) ! j-index in the buffer
DO ji = iis, iie
iib = nfd_jisnd(ji,iig) ! i-index in the buffer
iin = nfd_rksnd(ji,iig) ! neigbhour-index in the buffer
IF( nfd_rknei(iin) == -1 .AND. kfillmode == jpfillnothing ) CYCLE
ptab(jf)%pt4d(ji,ij1,jk,jl) = psgn(jf) * zbufr(iib,ij2,jk,jl,iin)
END DO END DO
! !
END DO ; END DO ! ipl ; ipk END DO ; END DO ! jk ; jl
! !
END DO ! ipf END DO ! jf
! !
DEALLOCATE( zbufr, ireq_r, ijsnd, ijbuf, ijrcv, ii1st, iiend, ipjfld ) DEALLOCATE( zbufr, ireq_r, ipjfld )
! !
CALL mpi_waitall(nfd_nbnei, ireq_s, MPI_STATUSES_IGNORE, ierr) ! wait for all Isend CALL mpi_waitall(nfd_nbnei, ireq_s, MPI_STATUSES_IGNORE, ierr) ! wait for all Isend
! !
...@@ -283,114 +214,128 @@ ...@@ -283,114 +214,128 @@
! !
ELSE !== allgather exchanges ==! ELSE !== allgather exchanges ==!
! !
! how many lines do we exchange at max? -> ipj (no further optimizations in this case...) DO jf = 1, kfld
ipj = khls + 2 ! how many lines do we send for each field?
! how many lines do we need at max? -> ipj2 (no further optimizations in this case...) ipj1(jf) = ihls(jf) + COUNT( (/ c_NFtype == 'T' .OR. cd_nat(jf) == 'V' .OR. cd_nat(jf) == 'F' /) ) &
ipj2 = 2 * khls + 2 & + COUNT( (/ c_NFtype == 'F' .AND. cd_nat(jf) == 'T' .AND. ihls(jf) == 0 /) )
! how many lines do we need for each field?
ipj2(jf) = 2 * ihls(jf) + COUNT( (/ c_NFtype == 'T' /) ) + COUNT( (/ cd_nat(jf) == 'V' .OR. cd_nat(jf) == 'F' /) ) &
& + COUNT( (/ c_NFtype == 'F' .AND. cd_nat(jf) == 'T' .AND. ihls(jf) == 0 /) )
END DO
! !
i0max = jpimax - 2 * khls i0max = MAXVAL( nfni_0, mask = nfproc /= -1 ) ! largest value of Ni_0 among processors (we are not sending halos)
ibuffsize = i0max * ipj * ipk * ipl * ipf ibuffsize = i0max * SUM( ipj1(:) * ipk(:) * ipl(:) ) ! use i0max because each proc must have the same buffer size
ALLOCATE( znorthloc(i0max,ipj,ipk,ipl,ipf), znorthglo(i0max,ipj,ipk,ipl,ipf,ndim_rank_north) ) ALLOCATE( znorthloc(i0max, ibuffsize/i0max), znorthall(i0max, ibuffsize/i0max, ndim_rank_north) )
! !
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk ! put in znorthloc ipj j-lines of ptab ij1 = 0 ! initalize line index
DO jj = 1, ipj DO jf = 1, kfld ; DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf)
ij2 = jpj - ipj2 + jj ! the first ipj lines of the last ipj2 lines DO jj = 1, ipj1(jf) ! put in znorthloc ipj1(jf) j-lines of ptab
ij2 = ipj(jf) - ipj2(jf) + jj ! the first ipj1 lines of the last ipj2 lines
ij1 = ij1 + 1
DO ji = 1, Ni_0 DO ji = 1, Ni_0
ii2 = Nis0 - 1 + ji ! inner domain: Nis0 to Nie0 ii2 = ihls(jf) + ji ! copy only the inner domain
znorthloc(ji,jj,jk,jl,jf) = ptab(jf)%pt4d(ii2,ij2,jk,jl) znorthloc(ji,ij1) = ptab(jf)%pt4d(ii2,ij2,jk,jl)
END DO END DO
DO ji = Ni_0+1, i0max DO ji = Ni_0+1, i0max ! avoid to send uninitialized values
znorthloc(ji,jj,jk,jl,jf) = HUGE(0._/**/PRECISION) ! avoid sending uninitialized values (make sure we don't use it) znorthloc(ji,ij1) = zhuge ! and make sure we don't use it
END DO END DO
END DO END DO
END DO ; END DO ; END DO END DO ; END DO ; END DO
! !
! start waiting time measurement ! start waiting time measurement
IF( ln_timing ) CALL tic_tac(.TRUE.)
#if ! defined key_mpi_off #if ! defined key_mpi_off
CALL MPI_ALLGATHER( znorthloc, ibuffsize, MPI_TYPE, znorthglo, ibuffsize, MPI_TYPE, ncomm_north, ierr ) IF( ln_timing ) CALL tic_tac( .TRUE.) ! start waiting time measurement
! fill znorthall with the znorthloc of each northern process
CALL MPI_ALLGATHER( znorthloc, ibuffsize, MPI_TYPE, znorthall, ibuffsize, MPI_TYPE, ncomm_north, ierr )
IF( ln_timing ) CALL tic_tac(.FALSE.) ! stop waiting time measurement
#endif #endif
! stop waiting time measurement DEALLOCATE( znorthloc ) ! no more need of znorthloc
IF( ln_timing ) CALL tic_tac(.FALSE.)
DEALLOCATE( znorthloc )
ALLOCATE( ztabglo(ipf) )
DO jf = 1, ipf
ALLOCATE( ztabglo(jf)%pt4d(jpiglo,ipj2,ipk,ipl) )
END DO
! !
! need to fill only the first ipj lines of ztabglo as lbc_nfd don't use the last khls lines DO jf = 1, kfld
ijnr = 0 !
DO jr = 1, jpni ! recover the global north array ihsz = ihls(jf) ! shorter name
iproc = nfproc(jr) iisht = nn_hls - ihsz
impp = nfimpp(jr) ALLOCATE( ztabglo(1)%pt4d(Ni0glo+2*ihsz,ipj2(jf),ipk(jf),ipl(jf)) )
ipi = nfjpi( jr) - 2 * khls ! corresponds to Ni_0 but for subdomain iproc !
IF( iproc == -1 ) THEN ! No neighbour (land proc that was suppressed) iFU = COUNT( (/ c_NFtype == 'F' .AND. cd_nat(jf) == 'U' /) ) ! F-folding and U grid
! IF( iFU == 0 ) ztabglo(1)%pt4d(:,ipj2(jf)-ihsz,:,:) = zhuge ! flag off the line that is not fully modified
SELECT CASE ( kfillmode ) !
CASE ( jpfillnothing ) ! no filling ! need to fill only the first ipj1(j) lines of ztabglo as lbc_nfd don't use the last ihsz lines
CALL ctl_stop( 'STOP', 'mpp_nfd_generic : cannot use jpfillnothing with ln_nnogather = F') ijnr = 0
CASE ( jpfillcopy ) ! filling with inner domain values DO jr = 1, jpni ! recover the global north array using each northern process
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk iproc = nfproc(jr) ! process number
DO jj = 1, ipj impp = nfimpp(jr) + ihsz ! ( = +nn_hls-iisht) ! inner domain position (without halos) of subdomain iproc
ij2 = jpj - ipj2 + jj ! the first ipj lines of the last ipj2 lines ipi0 = nfni_0(jr) ! Ni_0 but for subdomain iproc
DO ji = 1, ipi !
ii1 = impp + khls + ji - 1 ! corresponds to mig(khls + ji) but for subdomain iproc IF( iproc == -1 ) THEN ! No neighbour (land proc that was suppressed)
ztabglo(jf)%pt4d(ii1,jj,jk,jl) = ptab(jf)%pt4d(Nis0,ij2,jk,jl) ! chose to take the 1st inner domain point !
SELECT CASE ( kfillmode )
CASE ( jpfillnothing ) ! no filling
CALL ctl_stop( 'STOP', 'mpp_nfd_generic : cannot use jpfillnothing with ln_nnogather = F')
CASE ( jpfillcopy ) ! filling with inner domain values
DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf)
DO jj = 1, ipj1(jf)
ij2 = ipj(jf) - ipj2(jf) + jj ! the first ipj1(jf) lines of the last ipj2(jf) lines
DO ji = 1, ipi0
ii1 = impp + ji - 1 ! inner iproc-subdomain in the global domain with ihsz halos
ztabglo(1)%pt4d(ii1,jj,jk,jl) = ptab(jf)%pt4d(ihsz+1,ij2,jk,jl) ! take the 1st inner domain point
END DO
END DO END DO
END DO END DO ; END DO
END DO ; END DO ; END DO CASE ( jpfillcst ) ! filling with constant value
CASE ( jpfillcst ) ! filling with constant value DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf)
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk DO jj = 1, ipj1(jf)
DO jj = 1, ipj DO ji = 1, ipi0
DO ji = 1, ipi ii1 = impp + ji - 1 ! inner iproc-subdomain in the global domain with ihsz halos
ii1 = impp + khls + ji - 1 ! corresponds to mig(khls + ji) but for subdomain iproc ztabglo(1)%pt4d(ii1,jj,jk,jl) = pfillval
ztabglo(jf)%pt4d(ii1,jj,jk,jl) = pfillval END DO
END DO
END DO ; END DO
END SELECT
!
ELSE ! use neighbour values
ijnr = ijnr + 1
ij1 = SUM( ipj1(1:jf-1) * ipk(1:jf-1) * ipl(1:jf-1) ) ! reset line offset, return 0 if jf = 1
DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf)
DO jj = 1, ipj1(jf)
ij1 = ij1 + 1
DO ji = 1, ipi0
ii1 = impp + ji - 1 ! inner iproc-subdomain in the global domain with ihsz halos
ztabglo(1)%pt4d(ii1,jj,jk,jl) = znorthall(ji, ij1, ijnr)
END DO END DO
END DO END DO
END DO ; END DO ; END DO END DO ; END DO
END SELECT ENDIF
! !
ELSE END DO ! jpni
ijnr = ijnr + 1
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk
DO jj = 1, ipj
DO ji = 1, ipi
ii1 = impp + khls + ji - 1 ! corresponds to mig(khls + ji) but for subdomain iproc
ztabglo(jf)%pt4d(ii1,jj,jk,jl) = znorthglo(ji,jj,jk,jl,jf,ijnr)
END DO
END DO
END DO ; END DO ; END DO
ENDIF
! !
END DO ! jpni CALL lbc_nfd( ztabglo, cd_nat(jf:jf), psgn(jf:jf), 1 ) ! North fold boundary condition
DEALLOCATE( znorthglo ) !
! DO jl = 1, ipl(jf) ; DO jk = 1, ipk(jf) ! Scatter back to ARRAY_IN
DO jf = 1, ipf DO jj = 0, ihsz-1
CALL lbc_nfd( ztabglo(jf:jf), cd_nat(jf:jf), psgn(jf:jf), khls, 1 ) ! North fold boundary condition ij1 = ipj( jf) - jj ! last ihsz lines
DO jl = 1, ipl ; DO jk = 1, ipk ! e-w periodicity ij2 = ipj2(jf) - jj ! last ihsz lines
DO jj = 1, khls + 1 DO ji= 1, ipi(jf)
ij1 = ipj2 - (khls + 1) + jj ! need only the last khls + 1 lines until ipj2 ii2 = mig(ji+iisht,ihsz) ! warning, mig is expecting local domain indices related to nn_hls
ztabglo(jf)%pt4d( 1: khls,ij1,jk,jl) = ztabglo(jf)%pt4d(jpiglo-2*khls+1:jpiglo-khls,ij1,jk,jl) ptab(jf)%pt4d(ji,ij1,jk,jl) = ztabglo(1)%pt4d(ii2,ij2,jk,jl)
ztabglo(jf)%pt4d(jpiglo-khls+1:jpiglo,ij1,jk,jl) = ztabglo(jf)%pt4d( khls+1: 2*khls,ij1,jk,jl) END DO
END DO END DO
END DO ; END DO DO jj = ihsz, ihsz - iFU
END DO ij1 = ipj( jf) - jj ! last ihsz+1 line
! ij2 = ipj2(jf) - jj ! last ihsz+1 line
DO jf = 1, ipf ; DO jl = 1, ipl ; DO jk = 1, ipk ! Scatter back to ARRAY_IN DO ji= 1, ipi(jf)
DO jj = 1, khls + 1 ii2 = mig(ji+iisht,ihsz) ! warning, mig is expecting local domain indices related to nn_hls
ij1 = jpj - (khls + 1) + jj ! last khls + 1 lines until jpj zztmp = ztabglo(1)%pt4d(ii2,ij2,jk,jl)
ij2 = ipj2 - (khls + 1) + jj ! last khls + 1 lines until ipj2 IF( zztmp /= zhuge ) ptab(jf)%pt4d(ji,ij1,jk,jl) = zztmp ! apply it only if it was modified by lbc_nfd
DO ji= 1, jpi END DO
ii2 = mig(ji)
ptab(jf)%pt4d(ji,ij1,jk,jl) = ztabglo(jf)%pt4d(ii2,ij2,jk,jl)
END DO END DO
END DO END DO ; END DO
END DO ; END DO ; END DO !
DEALLOCATE( ztabglo(1)%pt4d )
!
END DO ! jf
! !
DO jf = 1, ipf DEALLOCATE( znorthall )
DEALLOCATE( ztabglo(jf)%pt4d )
END DO
DEALLOCATE( ztabglo )
! !
ENDIF ! ln_nnogather ENDIF ! ln_nnogather
! !
......
src/OCE/LBC/mppini.F90
View file @ 3efbcb17
...@@ -70,7 +70,8 @@ CONTAINS ...@@ -70,7 +70,8 @@ CONTAINS
jpi = jpiglo jpi = jpiglo
jpj = jpjglo jpj = jpjglo
jpk = MAX( 2, jpkglo ) jpk = MAX( 2, jpkglo )
jpij = jpi*jpj !jpij = jpi*jpj
jpij = Ni0glo*Nj0glo
jpni = 1 jpni = 1
jpnj = 1 jpnj = 1
jpnij = jpni*jpnj jpnij = jpni*jpnj
...@@ -165,6 +166,10 @@ CONTAINS ...@@ -165,6 +166,10 @@ CONTAINS
902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nammpp in configuration namelist' ) 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nammpp in configuration namelist' )
! !
nn_hls = MAX(1, nn_hls) ! nn_hls must be > 0 nn_hls = MAX(1, nn_hls) ! nn_hls must be > 0
# if defined key_mpi2
WRITE(numout,*) ' use key_mpi2, we force nn_comm = 1'
nn_comm = 1
# endif
IF(lwp) THEN IF(lwp) THEN
WRITE(numout,*) ' Namelist nammpp' WRITE(numout,*) ' Namelist nammpp'
IF( jpni < 1 .OR. jpnj < 1 ) THEN IF( jpni < 1 .OR. jpnj < 1 ) THEN
...@@ -303,7 +308,7 @@ CONTAINS ...@@ -303,7 +308,7 @@ CONTAINS
9002 FORMAT (a, i4, a) 9002 FORMAT (a, i4, a)
9003 FORMAT (a, i5) 9003 FORMAT (a, i5)
ALLOCATE( nfimpp(jpni), nfproc(jpni), nfjpi(jpni), & ALLOCATE( nfimpp(jpni), nfproc(jpni), nfjpi(jpni), nfni_0(jpni), &
& iin(jpnij), ijn(jpnij), & & iin(jpnij), ijn(jpnij), &
& iimppt(jpni,jpnj), ijmppt(jpni,jpnj), ijpi(jpni,jpnj), ijpj(jpni,jpnj), ipproc(jpni,jpnj), & & iimppt(jpni,jpnj), ijmppt(jpni,jpnj), ijpi(jpni,jpnj), ijpj(jpni,jpnj), ipproc(jpni,jpnj), &
& inei(8,jpni,jpnj), llnei(8,jpni,jpnj), & & inei(8,jpni,jpnj), llnei(8,jpni,jpnj), &
...@@ -327,7 +332,8 @@ CONTAINS ...@@ -327,7 +332,8 @@ CONTAINS
jpi = ijpi(ii,ij) jpi = ijpi(ii,ij)
jpj = ijpj(ii,ij) jpj = ijpj(ii,ij)
jpk = MAX( 2, jpkglo ) jpk = MAX( 2, jpkglo )
jpij = jpi*jpj !jpij = jpi*jpj
jpij = (jpi-2*nn_hls)*(jpj-2*nn_hls)
nimpp = iimppt(ii,ij) nimpp = iimppt(ii,ij)
njmpp = ijmppt(ii,ij) njmpp = ijmppt(ii,ij)
! !
...@@ -375,7 +381,8 @@ CONTAINS ...@@ -375,7 +381,8 @@ CONTAINS
! Store informations for the north pole folding communications ! Store informations for the north pole folding communications
nfproc(:) = ipproc(:,jpnj) nfproc(:) = ipproc(:,jpnj)
nfimpp(:) = iimppt(:,jpnj) nfimpp(:) = iimppt(:,jpnj)
nfjpi (:) = ijpi(:,jpnj) nfjpi (:) = ijpi(:,jpnj) ! needed only for mpp_lbc_north_icb_generic.h90
nfni_0(:) = ijpi(:,jpnj) - 2 * nn_hls
! !
! 3. Define Western, Eastern, Southern and Northern neighbors + corners in the subdomain grid reference ! 3. Define Western, Eastern, Southern and Northern neighbors + corners in the subdomain grid reference
! ------------------------------------------------------------------------------------------------------ ! ------------------------------------------------------------------------------------------------------
...@@ -489,7 +496,9 @@ CONTAINS ...@@ -489,7 +496,9 @@ CONTAINS
! ----------------------------------------- ! -----------------------------------------
! !
! set default neighbours ! set default neighbours
mpinei(:) = impi(:,narea) mpinei(:) = impi(:,narea) ! should be just local but is still used in icblbc and mpp_lnk_icb_generic.h90...
mpiSnei(0,:) = -1 ! no comm if no halo (but still need to call the NP Folding that may modify the last line)
mpiRnei(0,:) = -1
DO jh = 1, n_hlsmax DO jh = 1, n_hlsmax
mpiSnei(jh,:) = impi(:,narea) ! default definition mpiSnei(jh,:) = impi(:,narea) ! default definition
mpiRnei(jh,:) = impi(:,narea) mpiRnei(jh,:) = impi(:,narea)
...@@ -516,6 +525,7 @@ CONTAINS ...@@ -516,6 +525,7 @@ CONTAINS
! !
! ! Prepare mpp north fold ! ! Prepare mpp north fold
! !
l_NFold = l_NFold .AND. ANY( nfproc /= -1 ) ! make sure that we kept at least 1 proc along the last line
llmpiNFold = jpni > 1 .AND. l_NFold ! is the North fold done with an MPI communication? llmpiNFold = jpni > 1 .AND. l_NFold ! is the North fold done with an MPI communication?
l_IdoNFold = ijn(narea) == jpnj .AND. l_NFold ! is this process doing North fold? l_IdoNFold = ijn(narea) == jpnj .AND. l_NFold ! is this process doing North fold?
! !
...@@ -1179,7 +1189,7 @@ CONTAINS ...@@ -1179,7 +1189,7 @@ CONTAINS
! !
ALLOCATE( zmsk0(ipi,ipj), zmsk(ipi,ipj) ) ALLOCATE( zmsk0(ipi,ipj), zmsk(ipi,ipj) )
zmsk0(jh+1:jh+Ni_0,jh+1:jh+Nj_0) = REAL(COUNT(lloce, dim = 3), wp) ! define inner domain -> need REAL to use lbclnk zmsk0(jh+1:jh+Ni_0,jh+1:jh+Nj_0) = REAL(COUNT(lloce, dim = 3), wp) ! define inner domain -> need REAL to use lbclnk
CALL lbc_lnk('mppini', zmsk0, 'T', 1._wp, khls = jh) ! fill halos CALL lbc_lnk( ' mppini', zmsk0, 'T', 1._wp ) ! fill halos
! Beware about the mask we must use here : ! Beware about the mask we must use here :
DO jj = jh+1, jh+Nj_0 DO jj = jh+1, jh+Nj_0
DO ji = jh+1, jh+Ni_0 DO ji = jh+1, jh+Ni_0
...@@ -1192,7 +1202,7 @@ CONTAINS ...@@ -1192,7 +1202,7 @@ CONTAINS
& + zmsk0(ji+1,jj) + zmsk0(ji,jj+1) + zmsk0(ji+1,jj+1) & + zmsk0(ji+1,jj) + zmsk0(ji,jj+1) + zmsk0(ji+1,jj+1)
END DO END DO
END DO END DO
CALL lbc_lnk('mppini', zmsk, 'T', 1._wp, khls = jh) ! fill halos again! CALL lbc_lnk( 'mppini', zmsk, 'T', 1._wp ) ! fill halos again!
! !
iiwe = jh ; iiea = Ni_0 ! bottom-left corner - 1 of the sent data iiwe = jh ; iiea = Ni_0 ! bottom-left corner - 1 of the sent data
ijso = jh ; ijno = Nj_0 ijso = jh ; ijno = Nj_0
...@@ -1265,7 +1275,7 @@ ENDIF ...@@ -1265,7 +1275,7 @@ ENDIF
! used in IOM. This works even if jpnij .ne. jpni*jpnj. ! used in IOM. This works even if jpnij .ne. jpni*jpnj.
iglo( :) = (/ Ni0glo, Nj0glo /) iglo( :) = (/ Ni0glo, Nj0glo /)
iloc( :) = (/ Ni_0 , Nj_0 /) iloc( :) = (/ Ni_0 , Nj_0 /)
iabsf(:) = (/ Nis0 , Njs0 /) + (/ nimpp, njmpp /) - 1 - nn_hls ! corresponds to mig0(Nis0) but mig0 is not yet defined! iabsf(:) = (/ Nis0 , Njs0 /) + (/ nimpp, njmpp /) - 1 - nn_hls ! corresponds to mig(Nis0,0) but mig is not yet defined!
iabsl(:) = iabsf(:) + iloc(:) - 1 iabsl(:) = iabsf(:) + iloc(:) - 1
ihals(:) = (/ 0 , 0 /) ihals(:) = (/ 0 , 0 /)
ihale(:) = (/ 0 , 0 /) ihale(:) = (/ 0 , 0 /)
...@@ -1300,10 +1310,10 @@ ENDIF ...@@ -1300,10 +1310,10 @@ ENDIF
INTEGER, INTENT(in ) :: knum ! layout.dat unit INTEGER, INTENT(in ) :: knum ! layout.dat unit
! !
REAL(wp), DIMENSION(jpi,jpj,2,4) :: zinfo REAL(wp), DIMENSION(jpi,jpj,2,4) :: zinfo
INTEGER , DIMENSION(10) :: irknei ! too many elements but safe... INTEGER , DIMENSION(0:10) :: irknei ! too many elements but safe...
INTEGER :: ji, jj, jg, jn ! dummy loop indices INTEGER :: ji, jj, jg, jn ! dummy loop indices
INTEGER :: iitmp INTEGER :: iitmp
LOGICAL :: lnew LOGICAL :: llnew
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF (lwp) THEN IF (lwp) THEN
...@@ -1317,29 +1327,28 @@ ENDIF ...@@ -1317,29 +1327,28 @@ ENDIF
WRITE(knum,*) WRITE(knum,*)
WRITE(knum,*) WRITE(knum,*)
WRITE(knum,*) 'Number of subdomains located along the north fold : ', ndim_rank_north WRITE(knum,*) 'Number of subdomains located along the north fold : ', ndim_rank_north
WRITE(knum,*) 'Rank of the subdomains located along the north fold : ', ndim_rank_north WRITE(knum,*) 'Rank of the subdomains located along the north fold : '
DO jn = 1, ndim_rank_north, 5 DO jn = 1, ndim_rank_north, 5
WRITE(knum,*) nrank_north( jn:MINVAL( (/jn+4,ndim_rank_north/) ) ) WRITE(knum,*) nrank_north( jn:MINVAL( (/jn+4,ndim_rank_north/) ) )
END DO END DO
ENDIF ENDIF
nfd_nbnei = 0 ! defaul def (useless?) nfd_nbnei = 0 ! default def (useless?)
IF( ln_nnogather ) THEN IF( ln_nnogather ) THEN
! !
! Use the "gather nfd" to know how to do the nfd: for ji point, which process send data from which of its ji-index? ! Use the "gather nfd" to know how to do the nfd: for ji point, which process send data from which of its ji-index?
! Note that nfd is perfectly symetric: I receive data from X <=> I send data to X (-> no deadlock) ! Note that nfd is perfectly symetric: I receive data from X <=> I send data to X (-> no deadlock)
! !
zinfo(:,:,:,:) = HUGE(0._wp) ! default def to make sure we don't use the halos DO jg = 1, 4 ! grid type: T, U, V, F
DO jg = 1, 4 ! grid type: T, U, V, F
DO jj = nn_hls+1, jpj-nn_hls ! inner domain (warning do_loop_substitute not yet defined) DO jj = nn_hls+1, jpj-nn_hls ! inner domain (warning do_loop_substitute not yet defined)
DO ji = nn_hls+1, jpi-nn_hls ! inner domain (warning do_loop_substitute not yet defined) DO ji = nn_hls+1, jpi-nn_hls ! inner domain (warning do_loop_substitute not yet defined)
zinfo(ji,jj,1,jg) = REAL(narea, wp) ! mpi_rank + 1 (as default lbc_lnk fill is 0 zinfo(ji,jj,1,jg) = REAL(narea, wp) ! mpi_rank + 1 (note: lbc_lnk will put 0 if no neighbour)
zinfo(ji,jj,2,jg) = REAL(ji, wp) ! ji of this proc zinfo(ji,jj,2,jg) = REAL(ji, wp) ! ji of this proc
END DO END DO
END DO END DO
END DO END DO
! !
ln_nnogather = .FALSE. ! force "classical" North pole folding to fill all halos -> should be no more HUGE values... ln_nnogather = .FALSE. ! force "classical" North pole folding to fill all halos
CALL lbc_lnk( 'mppini', zinfo(:,:,:,1), 'T', 1._wp ) ! Do 4 calls instead of 1 to save memory as the nogather version CALL lbc_lnk( 'mppini', zinfo(:,:,:,1), 'T', 1._wp ) ! Do 4 calls instead of 1 to save memory as the nogather version
CALL lbc_lnk( 'mppini', zinfo(:,:,:,2), 'U', 1._wp ) ! creates buffer arrays with jpiglo as the first dimension CALL lbc_lnk( 'mppini', zinfo(:,:,:,2), 'U', 1._wp ) ! creates buffer arrays with jpiglo as the first dimension
CALL lbc_lnk( 'mppini', zinfo(:,:,:,3), 'V', 1._wp ) ! CALL lbc_lnk( 'mppini', zinfo(:,:,:,3), 'V', 1._wp ) !
...@@ -1348,24 +1357,52 @@ ENDIF ...@@ -1348,24 +1357,52 @@ ENDIF
IF( l_IdoNFold ) THEN ! only the procs involed in the NFD must take care of this IF( l_IdoNFold ) THEN ! only the procs involed in the NFD must take care of this
ALLOCATE( nfd_rksnd(jpi,4), nfd_jisnd(jpi,4) ) ! neighbour rand and remote ji-index for each grid (T, U, V, F) ALLOCATE( nfd_rksnd(jpi,nn_hls+1,4), nfd_jisnd(jpi,nn_hls+1,4), lnfd_same(jpi,4) )
nfd_rksnd(:,:) = NINT( zinfo(:, jpj, 1, :) ) - 1 ! neighbour MPI rank nfd_rksnd(:,:,:) = NINT( zinfo(:,jpj-nn_hls:jpj,1,:) ) - 1 ! neighbour MPI rank (-1 means no neighbour)
nfd_jisnd(:,:) = NINT( zinfo(:, jpj, 2, :) ) - nn_hls ! neighbour ji index (shifted as we don't send the halos) ! Use some tricks for mpp_nfd_generic.h90:
WHERE( nfd_rksnd == -1 ) nfd_jisnd = 1 ! use ji=1 if no neighbour, see mpp_nfd_generic.h90 ! 1) neighbour ji index (shifted as we don't send the halos)
nfd_jisnd(:,:,:) = NINT( zinfo(:,jpj-nn_hls:jpj,2,:) ) - nn_hls
nfd_nbnei = 1 ! Number of neighbour sending data for the nfd. We have at least 1 neighbour! ! 2) use ji=1 if no neighbour
irknei(1) = nfd_rksnd(1,1) ! which is the 1st one (I can be neighbour of myself, exclude land-proc are also ok) WHERE( nfd_rksnd == -1 ) nfd_jisnd = 1
! 3) control which points must be modified by the NP folding on line jpjglo-nn_hls
lnfd_same(:,:) = .TRUE.
IF( c_NFtype == 'T' ) THEN
lnfd_same(mi0(jpiglo/2+2,nn_hls):mi1(jpiglo-nn_hls,nn_hls), 1) = .FALSE.
lnfd_same(mi0(jpiglo/2+1,nn_hls):mi1(jpiglo-nn_hls,nn_hls), 2) = .FALSE.
lnfd_same(mi0( nn_hls+1,nn_hls):mi1(jpiglo-nn_hls,nn_hls),3:4) = .FALSE.
IF( l_Iperio ) THEN ! in case the ew-periodicity was done before calling the NP folding
lnfd_same(mi0( 1,nn_hls):mi1(nn_hls,nn_hls),1:4) = .FALSE.
lnfd_same(mi0(jpiglo-nn_hls+1,nn_hls):mi1(jpiglo,nn_hls),3:4) = .FALSE.
ENDIF
ELSEIF( c_NFtype == 'F' ) THEN
lnfd_same(mi0(jpiglo/2+1 ,nn_hls):mi1(jpiglo/2+1 ,nn_hls),1) = .FALSE.
lnfd_same(mi0(jpiglo-nn_hls,nn_hls):mi1(jpiglo-nn_hls ,nn_hls),1) = .FALSE.
lnfd_same(mi0(jpiglo/2+1 ,nn_hls):mi1(jpiglo-nn_hls ,nn_hls),3) = .FALSE.
lnfd_same(mi0(jpiglo/2+1 ,nn_hls):mi1(jpiglo-nn_hls-1,nn_hls),4) = .FALSE.
IF( l_Iperio ) THEN ! in case the ew-periodicity was done before calling the NP folding
lnfd_same(mi0(nn_hls,nn_hls):mi1(nn_hls ,nn_hls),1) = .FALSE.
lnfd_same(mi0( 1,nn_hls):mi1(nn_hls ,nn_hls),3) = .FALSE.
IF( nn_hls > 1 ) lnfd_same(mi0( 1,nn_hls):mi1(nn_hls-1,nn_hls),4) = .FALSE.
ENDIF
ENDIF
WHERE( lnfd_same ) nfd_jisnd(:,1,:) = HUGE(0) ! make sure we dont use it
nfd_nbnei = 0
irknei(0) = HUGE(0)
DO jg = 1, 4 DO jg = 1, 4
DO ji = 1, jpi ! we must be able to fill the full line including halos DO jj = 1, nn_hls+1
lnew = .TRUE. ! new neighbour? DO ji = 1, jpi ! we must be able to fill the full line including halos
DO jn = 1, nfd_nbnei IF( jj == 1 .AND. lnfd_same(ji,jg) ) CYCLE
IF( irknei(jn) == nfd_rksnd(ji,jg) ) lnew = .FALSE. ! already found llnew = .TRUE. ! new neighbour?
DO jn = 0, nfd_nbnei
IF( irknei(jn) == nfd_rksnd(ji,jj,jg) ) llnew = .FALSE. ! already found
END DO
IF( llnew ) THEN
jn = nfd_nbnei + 1
nfd_nbnei = jn
irknei(jn) = nfd_rksnd(ji,jj,jg)
ENDIF
END DO END DO
IF( lnew ) THEN
jn = nfd_nbnei + 1
nfd_nbnei = jn
irknei(jn) = nfd_rksnd(ji,jg)
ENDIF
END DO END DO
END DO END DO
...@@ -1373,14 +1410,20 @@ ENDIF ...@@ -1373,14 +1410,20 @@ ENDIF
nfd_rknei(:) = irknei(1:nfd_nbnei) nfd_rknei(:) = irknei(1:nfd_nbnei)
! re-number nfd_rksnd according to the indexes of nfd_rknei ! re-number nfd_rksnd according to the indexes of nfd_rknei
DO jg = 1, 4 DO jg = 1, 4
DO ji = 1, jpi DO jj = 1, nn_hls+1
iitmp = nfd_rksnd(ji,jg) ! must store a copy of nfd_rksnd(ji,jg) to make sure we don't change it twice DO ji = 1, jpi
DO jn = 1, nfd_nbnei IF( jj == 1 .AND. lnfd_same(ji,jg) ) THEN
IF( iitmp == nfd_rknei(jn) ) nfd_rksnd(ji,jg) = jn nfd_rksnd(ji,jj,jg) = HUGE(0) ! make sure we don't use it
ELSE
iitmp = nfd_rksnd(ji,jj,jg) ! must store a copy of nfd_rksnd(ji,jj,jg) so we don't change it twice
DO jn = 1, nfd_nbnei
IF( iitmp == nfd_rknei(jn) ) nfd_rksnd(ji,jj,jg) = jn
END DO
ENDIF
END DO END DO
END DO END DO
END DO END DO
IF( ldwrtlay ) THEN IF( ldwrtlay ) THEN
WRITE(knum,*) WRITE(knum,*)
WRITE(knum,*) 'north fold exchanges with explicit point-to-point messaging :' WRITE(knum,*) 'north fold exchanges with explicit point-to-point messaging :'
...@@ -1411,7 +1454,18 @@ ENDIF ...@@ -1411,7 +1454,18 @@ ENDIF
Ni_0 = Nie0 - Nis0 + 1 Ni_0 = Nie0 - Nis0 + 1
Nj_0 = Nje0 - Njs0 + 1 Nj_0 = Nje0 - Njs0 + 1
! !
jpkm1 = jpk-1 ! " " jpkm1 = jpk-1
!
ntile = 0 ! Initialise "no tile" by default
nijtile = 1
ntsi = Nis0
ntsj = Njs0
ntei = Nie0
ntej = Nje0
nthl = 0
nthr = 0
nthb = 0
ntht = 0
! !
END SUBROUTINE init_doloop END SUBROUTINE init_doloop
...@@ -1424,38 +1478,74 @@ ENDIF ...@@ -1424,38 +1478,74 @@ ENDIF
!! !!
!! ** Method : !! ** Method :
!! !!
!! ** Action : - mig , mjg : local domain indices ==> global domain, including halos, indices !! Local domain indices: Same values for the same point, different upper/lower bounds
!! - mig0, mjg0: local domain indices ==> global domain, excluding halos, indices !! e.g. with nn_hls = 2
!! jh = 0 x,x,3,...,jpi-2, x, x
!! jh = 1 x,2,3,...,jpi-2,jpi-1, x
!! jh = 2 1,2,3,...,jpi-2,jpi-1,jpi
!!
!! or jh = 0 x,x,3,...,Ni_0+2, x, x
!! jh = 1 x,2,3,...,Ni_0+2,Ni_0+3, x
!! jh = 2 1,2,3,...,Ni_0+2,Ni_0+3,Ni_0+4
!!
!! Global domain indices: different values for the same point, all starts at 1
!! e.g. with nn_hls = 2
!! jh = 0 1,2,3, ...,jpiglo-4, x, x,x,x
!! jh = 1 1,2,3, ...,jpiglo-4,jpiglo-3,jpiglo-2, x,x
!! jh = 2 1,2,3,...,jpiglo-4,jpiglo-3,jpiglo-2,jpiglo-1,jpiglo
!!
!! or jh = 0 1,2,3, ...,Ni0glo , x, x,x,x
!! jh = 1 1,2,3, ...,Ni0glo ,Ni0glo+1,Ni0glo+2, x,x
!! jh = 2 1,2,3,...,Ni0glo,Ni0glo+1,Ni0glo+2,Ni0glo+3,Ni0glo+4
!! ^
!! |
!! |
!! iimpp
!!
!! ** Action : - mig , mjg : local domain indices ==> global domain indices
!! - mi0 , mi1 : global domain indices ==> local domain indices !! - mi0 , mi1 : global domain indices ==> local domain indices
!! - mj0 , mj1 (if global point not in the local domain ==> mi0>mi1 and/or mj0>mj1) !! - mj0 , mj1 (if global point not in the local domain ==> mi0>mi1 and/or mj0>mj1)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER :: ji, jj ! dummy loop argument INTEGER :: ji, jj, jh ! dummy loop argument
INTEGER :: ipi, ipj, ipiglo, ipjglo, iimpp, ijmpp, ishft
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
ALLOCATE( mig(jpi), mjg(jpj), mig0(jpi), mjg0(jpj) ) ALLOCATE( mig(jpi , 0:nn_hls), mjg(jpj , 0:nn_hls) )
ALLOCATE( mi0(jpiglo), mi1(jpiglo), mj0(jpjglo), mj1(jpjglo) ) ALLOCATE( mi0(jpiglo, 0:nn_hls), mi1(jpiglo, 0:nn_hls), mj0(jpjglo, 0:nn_hls), mj1(jpjglo, 0:nn_hls) )
! !
DO ji = 1, jpi ! local domain indices ==> global domain indices, including halos DO jh = 0, nn_hls
mig(ji) = ji + nimpp - 1 !
END DO ishft = nn_hls - jh
DO jj = 1, jpj !
mjg(jj) = jj + njmpp - 1 ipi = Ni_0 + 2*jh ; ipj = Nj_0 + 2*jh
END DO ipiglo = Ni0glo + 2*jh ; ipjglo = Nj0glo + 2*jh
! ! local domain indices ==> global domain indices, excluding halos iimpp = nimpp - ishft ; ijmpp = njmpp - ishft
! !
mig0(:) = mig(:) - nn_hls ! local domain indices ==> global domain indices, including jh halos
mjg0(:) = mjg(:) - nn_hls !
! ! global domain, including halos, indices ==> local domain indices DO ji = ishft + 1, ishft + ipi
! ! (return (m.0,m.1)=(1,0) if data domain gridpoint is to the west/south of the mig(ji,jh) = ji + iimpp - 1
! ! local domain, or (m.0,m.1)=(jp.+1,jp.) to the east/north of local domain. END DO
DO ji = 1, jpiglo !
mi0(ji) = MAX( 1 , MIN( ji - nimpp + 1, jpi+1 ) ) DO jj = ishft + 1, ishft + ipj
mi1(ji) = MAX( 0 , MIN( ji - nimpp + 1, jpi ) ) mjg(jj,jh) = jj + ijmpp - 1
END DO END DO
DO jj = 1, jpjglo !
mj0(jj) = MAX( 1 , MIN( jj - njmpp + 1, jpj+1 ) ) ! global domain, including jh halos, indices ==> local domain indices
mj1(jj) = MAX( 0 , MIN( jj - njmpp + 1, jpj ) ) ! return (m.0,m.1)=(1,0) if data domain gridpoint is to the west/south of the
END DO ! local domain, or (m.0,m.1)=(jp.+1,jp.) to the east/north of local domain.
!
DO ji = 1, ipiglo
mi0(ji,jh) = MAX( 1 , MIN( ji - iimpp + 1, ipi+ishft+1 ) )
mi1(ji,jh) = MAX( 0 , MIN( ji - iimpp + 1, ipi+ishft ) )
END DO
!
DO jj = 1, ipjglo
mj0(jj,jh) = MAX( 1 , MIN( jj - ijmpp + 1, ipj+ishft+1 ) )
mj1(jj,jh) = MAX( 0 , MIN( jj - ijmpp + 1, ipj+ishft ) )
END DO
!
END DO ! jh
! !
END SUBROUTINE init_locglo END SUBROUTINE init_locglo
......
src/OCE/LDF/ldfdyn.F90
View file @ 3efbcb17
...@@ -321,10 +321,10 @@ CONTAINS ...@@ -321,10 +321,10 @@ CONTAINS
l_ldfdyn_time = .TRUE. ! will be calculated by call to ldf_dyn routine in step.F90 l_ldfdyn_time = .TRUE. ! will be calculated by call to ldf_dyn routine in step.F90
! !
! ! allocate arrays used in ldf_dyn. ! ! allocate arrays used in ldf_dyn.
ALLOCATE( dtensq(jpi,jpj,jpk) , dshesq(jpi,jpj,jpk) , esqt(jpi,jpj) , esqf(jpi,jpj) , STAT=ierr ) ALLOCATE( dtensq(A2D(1),jpk) , dshesq(A2D(1),jpk) , esqt(A2D(0)) , esqf(A2D(0)) , STAT=ierr )
IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'ldf_dyn_init: failed to allocate Smagorinsky arrays') IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'ldf_dyn_init: failed to allocate Smagorinsky arrays')
! !
DO_2D( 1, 1, 1, 1 ) ! Set local gridscale values DO_2D( 0, 0, 0, 0 ) ! Set local gridscale values
esqt(ji,jj) = ( 2._wp * e1e2t(ji,jj) / ( e1t(ji,jj) + e2t(ji,jj) ) )**2 esqt(ji,jj) = ( 2._wp * e1e2t(ji,jj) / ( e1t(ji,jj) + e2t(ji,jj) ) )**2
esqf(ji,jj) = ( 2._wp * e1e2f(ji,jj) / ( e1f(ji,jj) + e2f(ji,jj) ) )**2 esqf(ji,jj) = ( 2._wp * e1e2f(ji,jj) / ( e1f(ji,jj) + e2f(ji,jj) ) )**2
END_2D END_2D
...@@ -419,7 +419,7 @@ CONTAINS ...@@ -419,7 +419,7 @@ CONTAINS
! ! of |U|L^3/16 in blp case ! ! of |U|L^3/16 in blp case
DO jk = 1, jpkm1 DO jk = 1, jpkm1
! !
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 1, 0, 1 )
zdb = ( uu(ji,jj,jk,Kbb) * r1_e2u(ji,jj) - uu(ji-1,jj,jk,Kbb) * r1_e2u(ji-1,jj) ) & zdb = ( uu(ji,jj,jk,Kbb) * r1_e2u(ji,jj) - uu(ji-1,jj,jk,Kbb) * r1_e2u(ji-1,jj) ) &
& * r1_e1t(ji,jj) * e2t(ji,jj) & & * r1_e1t(ji,jj) * e2t(ji,jj) &
& - ( vv(ji,jj,jk,Kbb) * r1_e1v(ji,jj) - vv(ji,jj-1,jk,Kbb) * r1_e1v(ji,jj-1) ) & & - ( vv(ji,jj,jk,Kbb) * r1_e1v(ji,jj) - vv(ji,jj-1,jk,Kbb) * r1_e1v(ji,jj-1) ) &
...@@ -437,8 +437,6 @@ CONTAINS ...@@ -437,8 +437,6 @@ CONTAINS
! !
END DO END DO
! !
CALL lbc_lnk( 'ldfdyn', dtensq, 'T', 1.0_wp ) ! lbc_lnk on dshesq not needed
!
DO jk = 1, jpkm1 DO jk = 1, jpkm1
! !
DO_2D( 0, 0, 0, 0 ) ! T-point value DO_2D( 0, 0, 0, 0 ) ! T-point value
...@@ -455,7 +453,7 @@ CONTAINS ...@@ -455,7 +453,7 @@ CONTAINS
! !
END_2D END_2D
! !
DO_2D( 1, 0, 1, 0 ) ! F-point value DO_2D( 0, 0, 0, 0 ) ! F-point value
! !
zu2pv2_ij_p1 = uu(ji ,jj+1,jk, kbb) * uu(ji ,jj+1,jk, kbb) + vv(ji+1,jj ,jk, kbb) * vv(ji+1,jj ,jk, kbb) zu2pv2_ij_p1 = uu(ji ,jj+1,jk, kbb) * uu(ji ,jj+1,jk, kbb) + vv(ji+1,jj ,jk, kbb) * vv(ji+1,jj ,jk, kbb)
zu2pv2_ij = uu(ji ,jj ,jk, kbb) * uu(ji ,jj ,jk, kbb) + vv(ji ,jj ,jk, kbb) * vv(ji ,jj ,jk, kbb) zu2pv2_ij = uu(ji ,jj ,jk, kbb) * uu(ji ,jj ,jk, kbb) + vv(ji ,jj ,jk, kbb) * vv(ji ,jj ,jk, kbb)
......
src/OCE/LDF/ldfslp.F90
View file @ 3efbcb17
...@@ -113,16 +113,17 @@ CONTAINS ...@@ -113,16 +113,17 @@ CONTAINS
REAL(wp), INTENT(in), DIMENSION(:,:,:) :: pn2 ! Brunt-Vaisala frequency (locally ref.) REAL(wp), INTENT(in), DIMENSION(:,:,:) :: pn2 ! Brunt-Vaisala frequency (locally ref.)
!! !!
INTEGER :: ji , jj , jk ! dummy loop indices INTEGER :: ji , jj , jk ! dummy loop indices
INTEGER :: ii0, ii1 ! temporary integer
INTEGER :: ij0, ij1 ! temporary integer
REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16, zcofw, z1_slpmax ! local scalars REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16, zcofw, z1_slpmax ! local scalars
REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - -
REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - - REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - -
REAL(wp) :: zck, zfk, zbw ! - - REAL(wp) :: zck, zfk, zbw ! - -
REAL(wp) :: zdepu, zdepv ! - - REAL(wp) :: zdepu, zdepv ! - -
REAL(wp), DIMENSION(jpi,jpj) :: zslpml_hmlpu, zslpml_hmlpv REAL(wp) :: zuslp, zvslp ! - -
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgru, zwz, zdzr REAL(wp) :: zwslpi, zwslpj ! - -
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrv, zww REAL(wp), DIMENSION(A2D(0)) :: zslpml_hmlpu, zslpml_hmlpv
REAL(wp), DIMENSION(A2D(1),jpk) :: zdzr
REAL(wp), DIMENSION(A2D(1),jpk) :: zgru, zgrv
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwz, zww
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( ln_timing ) CALL timing_start('ldf_slp') IF( ln_timing ) CALL timing_start('ldf_slp')
...@@ -154,15 +155,15 @@ CONTAINS ...@@ -154,15 +155,15 @@ CONTAINS
ENDIF ENDIF
! !
zdzr(:,:,1) = 0._wp !== Local vertical density gradient at T-point == ! (evaluated from N^2) zdzr(:,:,1) = 0._wp !== Local vertical density gradient at T-point == ! (evaluated from N^2)
DO jk = 2, jpkm1 DO_3D( 1, 1, 1, 1, 2, jpkm1 )
! ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point ! ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point
! ! trick: tmask(ik ) = 0 => all pn2 = 0 => zdzr = 0 ! ! trick: tmask(ik ) = 0 => all pn2 = 0 => zdzr = 0
! ! else tmask(ik+1) = 0 => pn2(ik+1) = 0 => zdzr divides by 1 ! ! else tmask(ik+1) = 0 => pn2(ik+1) = 0 => zdzr divides by 1
! ! umask(ik+1) /= 0 => all pn2 /= 0 => zdzr divides by 2 ! ! umask(ik+1) /= 0 => all pn2 /= 0 => zdzr divides by 2
! ! NB: 1/(tmask+1) = (1-.5*tmask) substitute a / by a * ==> faster ! ! NB: 1/(tmask+1) = (1-.5*tmask) substitute a / by a * ==> faster
zdzr(:,:,jk) = zm1_g * ( prd(:,:,jk) + 1._wp ) & zdzr(ji,jj,jk) = zm1_g * ( prd(ji,jj,jk) + 1._wp ) &
& * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) & * ( pn2(ji,jj,jk) + pn2(ji,jj,jk+1) ) * ( 1._wp - 0.5_wp * tmask(ji,jj,jk+1) )
END DO END_3D
! !
! !== Slopes just below the mixed layer ==! ! !== Slopes just below the mixed layer ==!
CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr, Kmm ) ! output: uslpml, vslpml, wslpiml, wslpjml CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr, Kmm ) ! output: uslpml, vslpml, wslpiml, wslpjml
...@@ -231,28 +232,23 @@ CONTAINS ...@@ -231,28 +232,23 @@ CONTAINS
CALL lbc_lnk( 'ldfslp', zwz, 'U', -1.0_wp, zww, 'V', -1.0_wp ) ! lateral boundary conditions CALL lbc_lnk( 'ldfslp', zwz, 'U', -1.0_wp, zww, 'V', -1.0_wp ) ! lateral boundary conditions
! !
! !* horizontal Shapiro filter ! !* horizontal Shapiro filter
DO jk = 2, jpkm1 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! rows jj=2 and =jpjm1 only
DO_2D( 0, 0, 0, 0 ) ! rows jj=2 and =jpjm1 only zuslp = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) &
uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) &
& + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) &
& + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) &
& + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & & + 4.* zwz(ji ,jj ,jk) )
& + 4.* zwz(ji ,jj ,jk) ) zvslp = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) &
vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) &
& + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) &
& + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) &
& + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & & + 4.* zww(ji,jj ,jk) )
& + 4.* zww(ji,jj ,jk) )
END_2D
! !* decrease along coastal boundaries ! !* decrease along coastal boundaries
DO_2D( 0, 0, 0, 0 ) uslp(ji,jj,jk) = zuslp * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk ) ) * 0.5_wp &
uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk ) ) * 0.5_wp & & * ( umask(ji,jj ,jk) + umask(ji,jj ,jk+1) ) * 0.5_wp
& * ( umask(ji,jj ,jk) + umask(ji,jj ,jk+1) ) * 0.5_wp vslp(ji,jj,jk) = zvslp * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk ) ) * 0.5_wp &
vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk ) ) * 0.5_wp & & * ( vmask(ji ,jj,jk) + vmask(ji ,jj,jk+1) ) * 0.5_wp
& * ( vmask(ji ,jj,jk) + vmask(ji ,jj,jk+1) ) * 0.5_wp END_3D
END_2D
END DO
! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd ) ! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd )
! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd )
...@@ -291,29 +287,25 @@ CONTAINS ...@@ -291,29 +287,25 @@ CONTAINS
CALL lbc_lnk( 'ldfslp', zwz, 'T', -1.0_wp, zww, 'T', -1.0_wp ) ! lateral boundary conditions CALL lbc_lnk( 'ldfslp', zwz, 'T', -1.0_wp, zww, 'T', -1.0_wp ) ! lateral boundary conditions
! !
! !* horizontal Shapiro filter ! !* horizontal Shapiro filter
DO jk = 2, jpkm1 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! rows jj=2 and =jpjm1 only
DO_2D( 0, 0, 0, 0 ) ! rows jj=2 and =jpjm1 only zcofw = wmask(ji,jj,jk) * z1_16
zcofw = wmask(ji,jj,jk) * z1_16 zwslpi = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) &
wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) &
& + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) &
& + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) &
& + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & & + 4.* zwz(ji ,jj ,jk) ) * zcofw
& + 4.* zwz(ji ,jj ,jk) ) * zcofw
zwslpj = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) &
wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) &
& + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) &
& + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) &
& + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & & + 4.* zww(ji ,jj ,jk) ) * zcofw
& + 4.* zww(ji ,jj ,jk) ) * zcofw ! !* decrease in vicinity of topography
END_2D zck = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) &
! !* decrease in vicinity of topography & * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25
DO_2D( 0, 0, 0, 0 ) wslpi(ji,jj,jk) = zwslpi * zck
zck = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) & wslpj(ji,jj,jk) = zwslpj * zck
& * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25 END_3D
wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck
wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck
END_2D
END DO
! IV. Lateral boundary conditions ! IV. Lateral boundary conditions
! =============================== ! ===============================
...@@ -574,8 +566,8 @@ CONTAINS ...@@ -574,8 +566,8 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
REAL(wp), DIMENSION(:,:,:), INTENT(in) :: prd ! in situ density REAL(wp), DIMENSION(:,:,:), INTENT(in) :: prd ! in situ density
REAL(wp), DIMENSION(:,:,:), INTENT(in) :: pn2 ! Brunt-Vaisala frequency (locally ref.) REAL(wp), DIMENSION(:,:,:), INTENT(in) :: pn2 ! Brunt-Vaisala frequency (locally ref.)
REAL(wp), DIMENSION(:,:,:), INTENT(in) :: p_gru, p_grv ! i- & j-gradient of density (u- & v-pts) REAL(wp), DIMENSION(A2D(1),jpk), INTENT(in) :: p_gru, p_grv ! i- & j-gradient of density (u- & v-pts)
REAL(wp), DIMENSION(:,:,:), INTENT(in) :: p_dzr ! z-gradient of density (T-point) REAL(wp), DIMENSION(A2D(1),jpk), INTENT(in) :: p_dzr ! z-gradient of density (T-point)
INTEGER , INTENT(in) :: Kmm ! ocean time level indices INTEGER , INTENT(in) :: Kmm ! ocean time level indices
!! !!
INTEGER :: ji , jj , jk ! dummy loop indices INTEGER :: ji , jj , jk ! dummy loop indices
...@@ -591,11 +583,6 @@ CONTAINS ...@@ -591,11 +583,6 @@ CONTAINS
zm1_2g = -0.5_wp / grav zm1_2g = -0.5_wp / grav
z1_slpmax = 1._wp / rn_slpmax z1_slpmax = 1._wp / rn_slpmax
! !
uslpml (1,:) = 0._wp ; uslpml (jpi,:) = 0._wp
vslpml (1,:) = 0._wp ; vslpml (jpi,:) = 0._wp
wslpiml(1,:) = 0._wp ; wslpiml(jpi,:) = 0._wp
wslpjml(1,:) = 0._wp ; wslpjml(jpi,:) = 0._wp
!
! !== surface mixed layer mask ! ! !== surface mixed layer mask !
DO_3D( 1, 1, 1, 1, 1, jpk ) ! =1 inside the mixed layer, =0 otherwise DO_3D( 1, 1, 1, 1, 1, jpk ) ! =1 inside the mixed layer, =0 otherwise
ik = nmln(ji,jj) - 1 ik = nmln(ji,jj) - 1
...@@ -657,8 +644,6 @@ CONTAINS ...@@ -657,8 +644,6 @@ CONTAINS
wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik) wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik)
wslpjml(ji,jj) = zaj / ( zbj - zeps ) * tmask (ji,jj,ik) wslpjml(ji,jj) = zaj / ( zbj - zeps ) * tmask (ji,jj,ik)
END_2D END_2D
!!gm this lbc_lnk should be useless....
CALL lbc_lnk( 'ldfslp', uslpml , 'U', -1.0_wp , vslpml , 'V', -1.0_wp , wslpiml, 'W', -1.0_wp , wslpjml, 'W', -1.0_wp )
! !
END SUBROUTINE ldf_slp_mxl END SUBROUTINE ldf_slp_mxl
...@@ -700,8 +685,8 @@ CONTAINS ...@@ -700,8 +685,8 @@ CONTAINS
ELSE ! Madec operator : slopes at u-, v-, and w-points ELSE ! Madec operator : slopes at u-, v-, and w-points
IF(lwp) WRITE(numout,*) ' ==>>> iso operator (Madec)' IF(lwp) WRITE(numout,*) ' ==>>> iso operator (Madec)'
ALLOCATE( omlmask(jpi,jpj,jpk) , & ALLOCATE( omlmask(jpi,jpj,jpk) , &
& uslp(jpi,jpj,jpk) , uslpml(jpi,jpj) , wslpi(jpi,jpj,jpk) , wslpiml(jpi,jpj) , & & uslp(jpi,jpj,jpk) , uslpml(A2D(0)) , wslpi(jpi,jpj,jpk) , wslpiml(A2D(0)) , &
& vslp(jpi,jpj,jpk) , vslpml(jpi,jpj) , wslpj(jpi,jpj,jpk) , wslpjml(jpi,jpj) , STAT=ierr ) & vslp(jpi,jpj,jpk) , vslpml(A2D(0)) , wslpj(jpi,jpj,jpk) , wslpjml(A2D(0)) , STAT=ierr )
IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'ldf_slp_init : unable to allocate Madec operator slope ' ) IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'ldf_slp_init : unable to allocate Madec operator slope ' )
! Direction of lateral diffusion (tracers and/or momentum) ! Direction of lateral diffusion (tracers and/or momentum)
......
src/OCE/LDF/ldftra.F90
View file @ 3efbcb17
...@@ -794,8 +794,8 @@ CONTAINS ...@@ -794,8 +794,8 @@ CONTAINS
! !
INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zztmp ! local scalar REAL(wp) :: zztmp ! local scalar
REAL(wp), DIMENSION(A2D(nn_hls)) :: zw2d ! 2D workspace REAL(wp), DIMENSION(A2D(0)) :: zw2d ! 2D workspace
REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zw3d ! 3D workspace REAL(wp), DIMENSION(A2D(0),jpk) :: zw3d ! 3D workspace
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
!!gm I don't like this routine.... Crazy way of doing things, not optimal at all... !!gm I don't like this routine.... Crazy way of doing things, not optimal at all...
...@@ -806,42 +806,53 @@ CONTAINS ...@@ -806,42 +806,53 @@ CONTAINS
! !
! !== eiv velocities: calculate and output ==! ! !== eiv velocities: calculate and output ==!
! !
zw3d(:,:,jpk) = 0._wp ! bottom value always 0 zw3d(:,:,jpk) = 0._wp ! bottom value always 0
! !
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e2u e3u u_eiv = -dk[psi_uw] IF( iom_use('uoce_eiv') ) THEN
zw3d(ji,jj,jk) = ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk) ) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm) ) DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e2u e3u u_eiv = -dk[psi_uw]
END_3D zw3d(ji,jj,jk) = ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk) ) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm) )
CALL iom_put( "uoce_eiv", zw3d ) END_3D
CALL iom_put( "uoce_eiv", zw3d )
ENDIF
! !
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e1v e3v v_eiv = -dk[psi_vw] IF( iom_use('ueiv_masstr') ) THEN
zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk) ) / ( e1v(ji,jj) * e3v(ji,jj,jk,Kmm) ) DO_3D( 0, 0, 0, 0, 1, jpkm1 )
END_3D zw3d(ji,jj,jk) = rho0 * ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk) )
CALL iom_put( "voce_eiv", zw3d ) END_3D
CALL iom_put( "ueiv_masstr", zw3d ) ! mass transport in i-direction
ENDIF
! !
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e1 e2 w_eiv = dk[psix] + dk[psix] IF( iom_use('voce_eiv') ) THEN
zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk) - psi_vw(ji ,jj-1,jk) & DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e1v e3v v_eiv = -dk[psi_vw]
& + psi_uw(ji,jj,jk) - psi_uw(ji-1,jj ,jk) ) / e1e2t(ji,jj) zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk) ) / ( e1v(ji,jj) * e3v(ji,jj,jk,Kmm) )
END_3D END_3D
CALL iom_put( "woce_eiv", zw3d ) CALL iom_put( "voce_eiv", zw3d )
ENDIF
! !
IF( iom_use('weiv_masstr') ) THEN ! vertical mass transport & its square value IF( iom_use('veiv_masstr') ) THEN
DO_2D( 0, 0, 0, 0 ) DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zw2d(ji,jj) = rho0 * e1e2t(ji,jj) zw3d(ji,jj,jk) = rho0 * ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk) )
END_2D END_3D
DO jk = 1, jpk CALL iom_put( "veiv_masstr", zw3d ) ! mass transport in j-direction
zw3d(:,:,jk) = zw3d(:,:,jk) * zw2d(:,:) ENDIF
END DO !
CALL iom_put( "weiv_masstr" , zw3d ) IF( iom_use('woce_eiv') ) THEN
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e1 e2 w_eiv = dk[psix] + dk[psix]
zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk) - psi_vw(ji ,jj-1,jk) &
& + psi_uw(ji,jj,jk) - psi_uw(ji-1,jj ,jk) ) / e1e2t(ji,jj)
END_3D
CALL iom_put( "woce_eiv", zw3d )
ENDIF ENDIF
! !
IF( iom_use('ueiv_masstr') ) THEN IF( iom_use('weiv_masstr') ) THEN
zw3d(:,:,:) = 0.e0 DO_3D( 0, 0, 0, 0, 1, jpkm1 )
DO jk = 1, jpkm1 zw3d(ji,jj,jk) = rho0 * ( psi_vw(ji,jj,jk) - psi_vw(ji ,jj-1,jk) &
zw3d(:,:,jk) = rho0 * ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) ) & + psi_uw(ji,jj,jk) - psi_uw(ji-1,jj ,jk) )
END DO END_3D
CALL iom_put( "ueiv_masstr", zw3d ) ! mass transport in i-direction CALL iom_put( "weiv_masstr" , zw3d ) ! mass transport in z-direction
ENDIF ENDIF
! !
!
zztmp = 0.5_wp * rho0 * rcp zztmp = 0.5_wp * rho0 * rcp
IF( iom_use('ueiv_heattr') .OR. iom_use('ueiv_heattr3d') ) THEN IF( iom_use('ueiv_heattr') .OR. iom_use('ueiv_heattr3d') ) THEN
zw2d(:,:) = 0._wp zw2d(:,:) = 0._wp
...@@ -855,49 +866,47 @@ CONTAINS ...@@ -855,49 +866,47 @@ CONTAINS
CALL iom_put( "ueiv_heattr3d", zztmp * zw3d ) ! heat transport in i-direction CALL iom_put( "ueiv_heattr3d", zztmp * zw3d ) ! heat transport in i-direction
ENDIF ENDIF
! !
IF( iom_use('veiv_masstr') ) THEN IF( iom_use('veiv_heattr') .OR. iom_use('veiv_heattr3d') .OR. iom_use('sophteiv') ) THEN
zw3d(:,:,:) = 0.e0 zw2d(:,:) = 0._wp
DO jk = 1, jpkm1 zw3d(:,:,:) = 0._wp
zw3d(:,:,jk) = rho0 * ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) ) DO_3D( 0, 0, 0, 0, 1, jpkm1 )
END DO zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj ,jk) ) &
CALL iom_put( "veiv_masstr", zw3d ) ! mass transport in i-direction & * ( ts (ji,jj,jk,jp_tem,Kmm) + ts (ji,jj+1,jk,jp_tem,Kmm) )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk)
END_3D
CALL iom_put( "veiv_heattr" , zztmp * zw2d ) ! heat transport in j-direction
CALL iom_put( "veiv_heattr3d", zztmp * zw3d ) ! heat transport in j-direction
!
IF( iom_use( 'sophteiv' ) .AND. l_diaptr ) CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d )
ENDIF ENDIF
! !
zw2d(:,:) = 0._wp
zw3d(:,:,:) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj ,jk) ) &
& * ( ts (ji,jj,jk,jp_tem,Kmm) + ts (ji,jj+1,jk,jp_tem,Kmm) )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk)
END_3D
CALL iom_put( "veiv_heattr" , zztmp * zw2d ) ! heat transport in j-direction
CALL iom_put( "veiv_heattr3d", zztmp * zw3d ) ! heat transport in j-direction
!
IF( iom_use( 'sophteiv' ) .AND. l_diaptr ) CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d )
! !
zztmp = 0.5_wp * 0.5 zztmp = 0.5_wp * 0.5
IF( iom_use('ueiv_salttr') .OR. iom_use('ueiv_salttr3d')) THEN IF( iom_use('ueiv_salttr') .OR. iom_use('ueiv_salttr3d') ) THEN
zw2d(:,:) = 0._wp zw2d(:,:) = 0._wp
zw3d(:,:,:) = 0._wp zw3d(:,:,:) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zw3d(ji,jj,jk) = zw3d(ji,jj,jk) * ( psi_uw(ji,jj,jk+1) - psi_uw(ji ,jj,jk) ) & zw3d(ji,jj,jk) = zw3d(ji,jj,jk) * ( psi_uw(ji,jj,jk+1) - psi_uw(ji ,jj,jk) ) &
& * ( ts (ji,jj,jk,jp_sal,Kmm) + ts (ji+1,jj,jk,jp_sal,Kmm) ) & * ( ts (ji,jj,jk,jp_sal,Kmm) + ts (ji+1,jj,jk,jp_sal,Kmm) )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk)
END_3D END_3D
CALL iom_put( "ueiv_salttr", zztmp * zw2d ) ! salt transport in i-direction CALL iom_put( "ueiv_salttr", zztmp * zw2d ) ! salt transport in i-direction
CALL iom_put( "ueiv_salttr3d", zztmp * zw3d ) ! salt transport in i-direction CALL iom_put( "ueiv_salttr3d", zztmp * zw3d ) ! salt transport in i-direction
ENDIF ENDIF
zw2d(:,:) = 0._wp
zw3d(:,:,:) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj ,jk) ) &
& * ( ts (ji,jj,jk,jp_sal,Kmm) + ts (ji,jj+1,jk,jp_sal,Kmm) )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk)
END_3D
CALL iom_put( "veiv_salttr" , zztmp * zw2d ) ! salt transport in j-direction
CALL iom_put( "veiv_salttr3d", zztmp * zw3d ) ! salt transport in j-direction
! !
IF( iom_use( 'sopsteiv' ) .AND. l_diaptr ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d ) IF( iom_use('veiv_salttr') .OR. iom_use('veiv_salttr3d') .OR. iom_use('sopsteiv') ) THEN
zw2d(:,:) = 0._wp
zw3d(:,:,:) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj ,jk) ) &
& * ( ts (ji,jj,jk,jp_sal,Kmm) + ts (ji,jj+1,jk,jp_sal,Kmm) )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk)
END_3D
CALL iom_put( "veiv_salttr" , zztmp * zw2d ) ! salt transport in j-direction
CALL iom_put( "veiv_salttr3d", zztmp * zw3d ) ! salt transport in j-direction
!
IF( iom_use( 'sopsteiv' ) .AND. l_diaptr ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d )
ENDIF
! !
! !
END SUBROUTINE ldf_eiv_dia END SUBROUTINE ldf_eiv_dia
......
src/OCE/OBS/mpp_map.F90
View file @ 3efbcb17
...@@ -10,8 +10,8 @@ MODULE mpp_map ...@@ -10,8 +10,8 @@ MODULE mpp_map
!! mppmap_init : Initialize mppmap. !! mppmap_init : Initialize mppmap.
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
USE par_kind, ONLY : wp ! Precision variables USE par_kind, ONLY : wp ! Precision variables
USE par_oce , ONLY : jpi, jpj, Nis0, Nie0, Njs0, Nje0 ! Ocean parameters USE par_oce , ONLY : jpi, jpj, Nis0, Nie0, Njs0, Nje0, nn_hls ! Ocean parameters
USE dom_oce , ONLY : mig, mjg, narea ! Ocean space and time domain variables USE dom_oce , ONLY : mig, mjg, narea ! Ocean space and time domain variables
#if ! defined key_mpi_off #if ! defined key_mpi_off
USE lib_mpp , ONLY : mpi_comm_oce ! MPP library USE lib_mpp , ONLY : mpi_comm_oce ! MPP library
#endif #endif
...@@ -64,7 +64,7 @@ INCLUDE 'mpif.h' ...@@ -64,7 +64,7 @@ INCLUDE 'mpif.h'
imppmap(:,:) = 0 imppmap(:,:) = 0
! ! Setup local grid points ! ! Setup local grid points
imppmap(mig(1):mig(jpi),mjg(1):mjg(jpj)) = narea imppmap(mig(1,nn_hls):mig(jpi,nn_hls),mjg(1,nn_hls):mjg(jpj,nn_hls)) = narea
! Get global data ! Get global data
......
src/OCE/OBS/obs_grd_bruteforce.h90
View file @ 3efbcb17
...@@ -111,9 +111,9 @@ ...@@ -111,9 +111,9 @@
zmskg(:,:) = -1.e+10 zmskg(:,:) = -1.e+10
DO jj = kldj, klej DO jj = kldj, klej
DO ji = kldi, klei DO ji = kldi, klei
zlamg(mig(ji),mjg(jj)) = pglam(ji,jj) zlamg(mig(ji,nn_hls),mjg(jj,nn_hls)) = pglam(ji,jj)
zphig(mig(ji),mjg(jj)) = pgphi(ji,jj) zphig(mig(ji,nn_hls),mjg(jj,nn_hls)) = pgphi(ji,jj)
zmskg(mig(ji),mjg(jj)) = pmask(ji,jj) zmskg(mig(ji,nn_hls),mjg(jj,nn_hls)) = pmask(ji,jj)
END DO END DO
END DO END DO
CALL mpp_global_max( zlamg ) CALL mpp_global_max( zlamg )
......
src/OCE/OBS/obs_grid.F90
View file @ 3efbcb17
...@@ -280,9 +280,9 @@ CONTAINS ...@@ -280,9 +280,9 @@ CONTAINS
! Add various grids here. ! Add various grids here.
DO jj = 1, jpj DO jj = 1, jpj
DO ji = 1, jpi DO ji = 1, jpi
zlamg(mig(ji),mjg(jj)) = glamt(ji,jj) zlamg(mig(ji,nn_hls),mjg(jj,nn_hls)) = glamt(ji,jj)
zphig(mig(ji),mjg(jj)) = gphit(ji,jj) zphig(mig(ji,nn_hls),mjg(jj,nn_hls)) = gphit(ji,jj)
zmskg(mig(ji),mjg(jj)) = tmask(ji,jj,1) zmskg(mig(ji,nn_hls),mjg(jj,nn_hls)) = tmask(ji,jj,1)
END DO END DO
END DO END DO
CALL mpp_global_max( zlamg ) CALL mpp_global_max( zlamg )
......
src/OCE/OBS/obs_inter_sup.F90
View file @ 3efbcb17
...@@ -279,8 +279,8 @@ CONTAINS ...@@ -279,8 +279,8 @@ CONTAINS
! Pack interpolation data to be sent ! Pack interpolation data to be sent
DO ji = 1, itot DO ji = 1, itot
ii = mi1(igrdij_recv(2*ji-1)) ii = mi1(igrdij_recv(2*ji-1),nn_hls)
ij = mj1(igrdij_recv(2*ji)) ij = mj1(igrdij_recv(2*ji ),nn_hls)
DO jk = 1, kpk DO jk = 1, kpk
zsend(jk,ji) = pval(ii,ij,jk) zsend(jk,ji) = pval(ii,ij,jk)
END DO END DO
......
src/OCE/OBS/obs_write.F90
View file @ 3efbcb17
...@@ -245,8 +245,8 @@ CONTAINS ...@@ -245,8 +245,8 @@ CONTAINS
fbdata%iobsi(jo,jvar) = profdata%mi(jo,jvar) fbdata%iobsi(jo,jvar) = profdata%mi(jo,jvar)
fbdata%iobsj(jo,jvar) = profdata%mj(jo,jvar) fbdata%iobsj(jo,jvar) = profdata%mj(jo,jvar)
ELSE ELSE
fbdata%iobsi(jo,jvar) = mig(profdata%mi(jo,jvar)) fbdata%iobsi(jo,jvar) = mig(profdata%mi(jo,jvar),nn_hls)
fbdata%iobsj(jo,jvar) = mjg(profdata%mj(jo,jvar)) fbdata%iobsj(jo,jvar) = mjg(profdata%mj(jo,jvar),nn_hls)
ENDIF ENDIF
END DO END DO
CALL greg2jul( 0, & CALL greg2jul( 0, &
...@@ -511,8 +511,8 @@ CONTAINS ...@@ -511,8 +511,8 @@ CONTAINS
fbdata%iobsi(jo,1) = surfdata%mi(jo) fbdata%iobsi(jo,1) = surfdata%mi(jo)
fbdata%iobsj(jo,1) = surfdata%mj(jo) fbdata%iobsj(jo,1) = surfdata%mj(jo)
ELSE ELSE
fbdata%iobsi(jo,1) = mig(surfdata%mi(jo)) fbdata%iobsi(jo,1) = mig(surfdata%mi(jo),nn_hls)
fbdata%iobsj(jo,1) = mjg(surfdata%mj(jo)) fbdata%iobsj(jo,1) = mjg(surfdata%mj(jo),nn_hls)
ENDIF ENDIF
CALL greg2jul( 0, & CALL greg2jul( 0, &
& surfdata%nmin(jo), & & surfdata%nmin(jo), &
......
src/OCE/SBC/cpl_oasis3.F90
View file @ 3efbcb17
...@@ -64,7 +64,6 @@ MODULE cpl_oasis3 ...@@ -64,7 +64,6 @@ MODULE cpl_oasis3
INTEGER :: nrcv ! total number of fields received INTEGER :: nrcv ! total number of fields received
INTEGER :: nsnd ! total number of fields sent INTEGER :: nsnd ! total number of fields sent
INTEGER :: ncplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data INTEGER :: ncplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data
INTEGER, PUBLIC, PARAMETER :: nmaxfld=62 ! Maximum number of coupling fields
INTEGER, PUBLIC, PARAMETER :: nmaxcat=5 ! Maximum number of coupling fields INTEGER, PUBLIC, PARAMETER :: nmaxcat=5 ! Maximum number of coupling fields
INTEGER, PUBLIC, PARAMETER :: nmaxcpl=5 ! Maximum number of coupling fields INTEGER, PUBLIC, PARAMETER :: nmaxcpl=5 ! Maximum number of coupling fields
...@@ -78,7 +77,7 @@ MODULE cpl_oasis3 ...@@ -78,7 +77,7 @@ MODULE cpl_oasis3
INTEGER :: ncplmodel ! Maximum number of models to/from which this variable may be sent/received INTEGER :: ncplmodel ! Maximum number of models to/from which this variable may be sent/received
END TYPE FLD_CPL END TYPE FLD_CPL
TYPE(FLD_CPL), DIMENSION(nmaxfld), PUBLIC :: srcv, ssnd !: Coupling fields TYPE(FLD_CPL), DIMENSION(:), ALLOCATABLE, PUBLIC :: srcv, ssnd !: Coupling fields
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: exfld ! Temporary buffer for receiving REAL(wp), DIMENSION(:,:), ALLOCATABLE :: exfld ! Temporary buffer for receiving
...@@ -153,15 +152,6 @@ CONTAINS ...@@ -153,15 +152,6 @@ CONTAINS
CALL oasis_abort ( ncomp_id, 'cpl_define', 'ncplmodel is larger than nmaxcpl, increase nmaxcpl') ; RETURN CALL oasis_abort ( ncomp_id, 'cpl_define', 'ncplmodel is larger than nmaxcpl, increase nmaxcpl') ; RETURN
ENDIF ENDIF
nrcv = krcv
IF( nrcv > nmaxfld ) THEN
CALL oasis_abort ( ncomp_id, 'cpl_define', 'nrcv is larger than nmaxfld, increase nmaxfld') ; RETURN
ENDIF
nsnd = ksnd
IF( nsnd > nmaxfld ) THEN
CALL oasis_abort ( ncomp_id, 'cpl_define', 'nsnd is larger than nmaxfld, increase nmaxfld') ; RETURN
ENDIF
! !
! ... Define the shape for the area that excludes the halo as we don't want them to be "seen" by oasis ! ... Define the shape for the area that excludes the halo as we don't want them to be "seen" by oasis
! !
...@@ -181,11 +171,11 @@ CONTAINS ...@@ -181,11 +171,11 @@ CONTAINS
! ... Define the partition, excluding halos as we don't want them to be "seen" by oasis ! ... Define the partition, excluding halos as we don't want them to be "seen" by oasis
! ----------------------------------------------------------------- ! -----------------------------------------------------------------
paral(1) = 2 ! box partitioning paral(1) = 2 ! box partitioning
paral(2) = Ni0glo * mjg0(nn_hls) + mig0(nn_hls) ! NEMO lower left corner global offset, without halos paral(2) = Ni0glo * mjg(nn_hls,0) + mig(nn_hls,0) ! NEMO lower left corner global offset, without halos
paral(3) = Ni_0 ! local extent in i, excluding halos paral(3) = Ni_0 ! local extent in i, excluding halos
paral(4) = Nj_0 ! local extent in j, excluding halos paral(4) = Nj_0 ! local extent in j, excluding halos
paral(5) = Ni0glo ! global extent in x, excluding halos paral(5) = Ni0glo ! global extent in x, excluding halos
IF( sn_cfctl%l_oasout ) THEN IF( sn_cfctl%l_oasout ) THEN
WRITE(numout,*) ' multiexchg: paral (1:5)', paral WRITE(numout,*) ' multiexchg: paral (1:5)', paral
...@@ -419,6 +409,7 @@ CONTAINS ...@@ -419,6 +409,7 @@ CONTAINS
IF( .NOT. ll_1st ) THEN IF( .NOT. ll_1st ) THEN
CALL lbc_lnk( 'cpl_oasis3', pdata(:,:,jc), srcv(kid)%clgrid, srcv(kid)%nsgn ) CALL lbc_lnk( 'cpl_oasis3', pdata(:,:,jc), srcv(kid)%clgrid, srcv(kid)%nsgn )
ENDIF ENDIF
!!clem: mettre T instead of clgrid
ENDDO ENDDO
! !
......
src/OCE/SBC/cyclone.F90
View file @ 3efbcb17
...@@ -54,9 +54,9 @@ CONTAINS ...@@ -54,9 +54,9 @@ CONTAINS
!! ** Action : - open TC data, find TCs for the current timestep !! ** Action : - open TC data, find TCs for the current timestep
!! - for each potential TC, add the winds on the grid !! - for each potential TC, add the winds on the grid
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER , INTENT(in) :: kt ! time step index INTEGER , INTENT(in) :: kt ! time step index
REAL(wp), INTENT(out), DIMENSION(jpi,jpj) :: pwnd_i ! wind speed i-components at T-point ORCA direction REAL(wp), INTENT(out), DIMENSION(A2D(0)) :: pwnd_i ! wind speed i-components at T-point ORCA direction
REAL(wp), INTENT(out), DIMENSION(jpi,jpj) :: pwnd_j ! wind speed j-components at T-point ORCA direction REAL(wp), INTENT(out), DIMENSION(A2D(0)) :: pwnd_j ! wind speed j-components at T-point ORCA direction
! !
!! !!
INTEGER :: ji, jj , jtc ! loop arguments INTEGER :: ji, jj , jtc ! loop arguments
...@@ -79,7 +79,7 @@ CONTAINS ...@@ -79,7 +79,7 @@ CONTAINS
REAL(wp) :: zwnd_r, zwnd_t ! radial and tangential components of the wind REAL(wp) :: zwnd_r, zwnd_t ! radial and tangential components of the wind
REAL(wp) :: zvmax ! timestep interpolated vmax REAL(wp) :: zvmax ! timestep interpolated vmax
REAL(wp) :: zrlon, zrlat ! temporary REAL(wp) :: zrlon, zrlat ! temporary
REAL(wp), DIMENSION(jpi,jpj) :: zwnd_x, zwnd_y ! zonal and meridional components of the wind REAL(wp), DIMENSION(A2D(0)) :: zwnd_x, zwnd_y ! zonal and meridional components of the wind
REAL(wp), DIMENSION(14,5) :: ztct ! tropical cyclone track data at kt REAL(wp), DIMENSION(14,5) :: ztct ! tropical cyclone track data at kt
! !
CHARACTER(len=100) :: cn_dir ! Root directory for location of files CHARACTER(len=100) :: cn_dir ! Root directory for location of files
...@@ -146,7 +146,7 @@ CONTAINS ...@@ -146,7 +146,7 @@ CONTAINS
! fitted B parameter (Willoughby 2004) ! fitted B parameter (Willoughby 2004)
zb = 2. zb = 2.
DO_2D( 1, 1, 1, 1 ) DO_2D( 0, 0, 0, 0 )
! calc distance between TC center and any point following great circle ! calc distance between TC center and any point following great circle
! source : http://www.movable-type.co.uk/scripts/latlong.html ! source : http://www.movable-type.co.uk/scripts/latlong.html
...@@ -207,7 +207,7 @@ CONTAINS ...@@ -207,7 +207,7 @@ CONTAINS
zA=0 zA=0
ENDIF ENDIF
DO_2D( 1, 1, 1, 1 ) DO_2D( 0, 0, 0, 0 )
zzrglam = rad * glamt(ji,jj) - zrlon zzrglam = rad * glamt(ji,jj) - zrlon
zzrgphi = rad * gphit(ji,jj) zzrgphi = rad * gphit(ji,jj)
......
src/OCE/SBC/fldread.F90
View file @ 3efbcb17
...@@ -699,9 +699,10 @@ CONTAINS ...@@ -699,9 +699,10 @@ CONTAINS
INTEGER :: imf ! size of the structure sd INTEGER :: imf ! size of the structure sd
INTEGER :: ill ! character length INTEGER :: ill ! character length
INTEGER :: iv ! indice of V component INTEGER :: iv ! indice of V component
INTEGER :: ipi, ipj
CHARACTER (LEN=100) :: clcomp ! dummy weight name CHARACTER (LEN=100) :: clcomp ! dummy weight name
REAL(wp), DIMENSION(jpi,jpj) :: utmp, vtmp ! temporary arrays for vector rotation REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utmp, vtmp ! temporary arrays for vector rotation
REAL(wp), DIMENSION(:,:,:), POINTER :: dta_u, dta_v ! short cut REAL(wp), DIMENSION(:,:,:), POINTER :: dta_u, dta_v ! short cut
!!--------------------------------------------------------------------- !!---------------------------------------------------------------------
! !
!! (sga: following code should be modified so that pairs arent searched for each time !! (sga: following code should be modified so that pairs arent searched for each time
...@@ -724,11 +725,14 @@ CONTAINS ...@@ -724,11 +725,14 @@ CONTAINS
IF( sd(ju)%ln_tint ) THEN ; dta_u => sd(ju)%fdta(:,:,:,jn) ; dta_v => sd(iv)%fdta(:,:,:,jn) IF( sd(ju)%ln_tint ) THEN ; dta_u => sd(ju)%fdta(:,:,:,jn) ; dta_v => sd(iv)%fdta(:,:,:,jn)
ELSE ; dta_u => sd(ju)%fnow(:,:,: ) ; dta_v => sd(iv)%fnow(:,:,: ) ELSE ; dta_u => sd(ju)%fnow(:,:,: ) ; dta_v => sd(iv)%fnow(:,:,: )
ENDIF ENDIF
ipi = SIZE(dta_u,1) ; ipj = SIZE(dta_u,2)
ALLOCATE( utmp(ipi,ipj), vtmp(ipi,ipj) )
DO jk = 1, SIZE( sd(ju)%fnow, 3 ) DO jk = 1, SIZE( sd(ju)%fnow, 3 )
CALL rot_rep( dta_u(:,:,jk), dta_v(:,:,jk), 'T', 'en->i', utmp(:,:) ) CALL rot_rep( dta_u(:,:,jk), dta_v(:,:,jk), 'T', 'en->i', utmp(:,:) )
CALL rot_rep( dta_u(:,:,jk), dta_v(:,:,jk), 'T', 'en->j', vtmp(:,:) ) CALL rot_rep( dta_u(:,:,jk), dta_v(:,:,jk), 'T', 'en->j', vtmp(:,:) )
dta_u(:,:,jk) = utmp(:,:) ; dta_v(:,:,jk) = vtmp(:,:) dta_u(:,:,jk) = utmp(:,:) ; dta_v(:,:,jk) = vtmp(:,:)
END DO END DO
DEALLOCATE( utmp, vtmp )
sd(ju)%rotn(jn) = .TRUE. ! vector was rotated sd(ju)%rotn(jn) = .TRUE. ! vector was rotated
IF( lwp .AND. kt == nit000 ) WRITE(numout,*) & IF( lwp .AND. kt == nit000 ) WRITE(numout,*) &
& 'fld_read: vector pair ('//TRIM(sd(ju)%clvar)//', '//TRIM(sd(iv)%clvar)//') rotated on to model grid' & 'fld_read: vector pair ('//TRIM(sd(ju)%clvar)//', '//TRIM(sd(iv)%clvar)//') rotated on to model grid'
...@@ -1100,7 +1104,7 @@ CONTAINS ...@@ -1100,7 +1104,7 @@ CONTAINS
CHARACTER (len=5) :: clname ! CHARACTER (len=5) :: clname !
INTEGER , DIMENSION(4) :: ddims INTEGER , DIMENSION(4) :: ddims
INTEGER :: isrc INTEGER :: isrc
REAL(wp), DIMENSION(jpi,jpj) :: data_tmp REAL(wp), DIMENSION(A2D(0)) :: data_tmp
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( nxt_wgt > tot_wgts ) THEN IF( nxt_wgt > tot_wgts ) THEN
...@@ -1155,9 +1159,9 @@ CONTAINS ...@@ -1155,9 +1159,9 @@ CONTAINS
ELSE ; ref_wgts(nxt_wgt)%numwgt = 16 ELSE ; ref_wgts(nxt_wgt)%numwgt = 16
ENDIF ENDIF
ALLOCATE( ref_wgts(nxt_wgt)%data_jpi(Nis0:Nie0,Njs0:Nje0,4) ) ALLOCATE( ref_wgts(nxt_wgt)%data_jpi(A2D(0),4) )
ALLOCATE( ref_wgts(nxt_wgt)%data_jpj(Nis0:Nie0,Njs0:Nje0,4) ) ALLOCATE( ref_wgts(nxt_wgt)%data_jpj(A2D(0),4) )
ALLOCATE( ref_wgts(nxt_wgt)%data_wgt(Nis0:Nie0,Njs0:Nje0,ref_wgts(nxt_wgt)%numwgt) ) ALLOCATE( ref_wgts(nxt_wgt)%data_wgt(A2D(0),ref_wgts(nxt_wgt)%numwgt) )
DO jn = 1,4 DO jn = 1,4
WRITE(clname,'(a3,i2.2)') 'src',jn WRITE(clname,'(a3,i2.2)') 'src',jn
...@@ -1332,7 +1336,9 @@ CONTAINS ...@@ -1332,7 +1336,9 @@ CONTAINS
INTEGER, DIMENSION(3) :: rec1_lsm, recn_lsm ! temporary arrays for start and length in case of seaoverland INTEGER, DIMENSION(3) :: rec1_lsm, recn_lsm ! temporary arrays for start and length in case of seaoverland
INTEGER :: ii_lsm1,ii_lsm2,ij_lsm1,ij_lsm2 ! temporary indices INTEGER :: ii_lsm1,ii_lsm2,ij_lsm1,ij_lsm2 ! temporary indices
INTEGER :: ji, jj, jk, jn, jir, jjr ! loop counters INTEGER :: ji, jj, jk, jn, jir, jjr ! loop counters
INTEGER :: ipk INTEGER :: ipi, ipj, ipk
INTEGER :: iisht, ijsht
INTEGER :: ii, ij
INTEGER :: ni, nj ! lengths INTEGER :: ni, nj ! lengths
INTEGER :: jpimin,jpiwid ! temporary indices INTEGER :: jpimin,jpiwid ! temporary indices
INTEGER :: jpimin_lsm,jpiwid_lsm ! temporary indices INTEGER :: jpimin_lsm,jpiwid_lsm ! temporary indices
...@@ -1343,7 +1349,11 @@ CONTAINS ...@@ -1343,7 +1349,11 @@ CONTAINS
INTEGER :: itmpi,itmpj,itmpz ! lengths INTEGER :: itmpi,itmpj,itmpz ! lengths
REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: ztmp_fly_dta ! local array of values on input grid REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: ztmp_fly_dta ! local array of values on input grid
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
ipi = SIZE(dta, 1)
ipj = SIZE(dta, 2)
ipk = SIZE(dta, 3) ipk = SIZE(dta, 3)
iisht = ( jpi - ipi ) / 2
ijsht = ( jpj - ipj ) / 2
! !
!! for weighted interpolation we have weights at four corners of a box surrounding !! for weighted interpolation we have weights at four corners of a box surrounding
!! a model grid point, each weight is multiplied by a grid value (bilinear case) !! a model grid point, each weight is multiplied by a grid value (bilinear case)
...@@ -1438,7 +1448,9 @@ CONTAINS ...@@ -1438,7 +1448,9 @@ CONTAINS
DO_3D( 0, 0, 0, 0, 1,ipk ) DO_3D( 0, 0, 0, 0, 1,ipk )
ni = ref_wgts(kw)%data_jpi(ji,jj,jn) + 1 ni = ref_wgts(kw)%data_jpi(ji,jj,jn) + 1
nj = ref_wgts(kw)%data_jpj(ji,jj,jn) + 1 nj = ref_wgts(kw)%data_jpj(ji,jj,jn) + 1
dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn) * ref_wgts(kw)%fly_dta(ni,nj,jk) ii = ji - iisht
ij = jj - ijsht
dta(ii,ij,jk) = dta(ii,ij,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn) * ref_wgts(kw)%fly_dta(ni,nj,jk)
END_3D END_3D
END DO END DO
...@@ -1482,7 +1494,9 @@ CONTAINS ...@@ -1482,7 +1494,9 @@ CONTAINS
!!$ DO_3D( 0, 0, 0, 0, 1,ipk ) !!$ DO_3D( 0, 0, 0, 0, 1,ipk )
!!$ ni = ref_wgts(kw)%data_jpi(ji,jj,jn) + 1 !!$ ni = ref_wgts(kw)%data_jpi(ji,jj,jn) + 1
!!$ nj = ref_wgts(kw)%data_jpj(ji,jj,jn) + 1 !!$ nj = ref_wgts(kw)%data_jpj(ji,jj,jn) + 1
!!$ dta(ji,jj,jk) = dta(ji,jj,jk) & !!$ ii = ji - iisht
!!$ ij = jj - ijsht
!!$ dta(ii,ij,jk) = dta(ii,ij,jk) &
!!$ ! gradient in the i direction !!$ ! gradient in the i direction
!!$ & + ref_wgts(kw)%data_wgt(ji,jj,jn+4) * 0.5_wp * & !!$ & + ref_wgts(kw)%data_wgt(ji,jj,jn+4) * 0.5_wp * &
!!$ & (ref_wgts(kw)%fly_dta(ni+1,nj ,jk) - ref_wgts(kw)%fly_dta(ni-1,nj ,jk)) & !!$ & (ref_wgts(kw)%fly_dta(ni+1,nj ,jk) - ref_wgts(kw)%fly_dta(ni-1,nj ,jk)) &
...@@ -1500,8 +1514,10 @@ CONTAINS ...@@ -1500,8 +1514,10 @@ CONTAINS
DO_3D( 0, 0, 0, 0, 1,ipk ) DO_3D( 0, 0, 0, 0, 1,ipk )
ni = ref_wgts(kw)%data_jpi(ji,jj,jn) ni = ref_wgts(kw)%data_jpi(ji,jj,jn)
nj = ref_wgts(kw)%data_jpj(ji,jj,jn) nj = ref_wgts(kw)%data_jpj(ji,jj,jn)
ii = ji - iisht
ij = jj - ijsht
! gradient in the i direction ! gradient in the i direction
dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+4) * 0.5_wp * & dta(ii,ij,jk) = dta(ii,ij,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+4) * 0.5_wp * &
& (ref_wgts(kw)%fly_dta(ni+2,nj+1,jk) - ref_wgts(kw)%fly_dta(ni ,nj+1,jk)) & (ref_wgts(kw)%fly_dta(ni+2,nj+1,jk) - ref_wgts(kw)%fly_dta(ni ,nj+1,jk))
END_3D END_3D
END DO END DO
...@@ -1509,8 +1525,10 @@ CONTAINS ...@@ -1509,8 +1525,10 @@ CONTAINS
DO_3D( 0, 0, 0, 0, 1,ipk ) DO_3D( 0, 0, 0, 0, 1,ipk )
ni = ref_wgts(kw)%data_jpi(ji,jj,jn) ni = ref_wgts(kw)%data_jpi(ji,jj,jn)
nj = ref_wgts(kw)%data_jpj(ji,jj,jn) nj = ref_wgts(kw)%data_jpj(ji,jj,jn)
ii = ji - iisht
ij = jj - ijsht
! gradient in the j direction ! gradient in the j direction
dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+8) * 0.5_wp * & dta(ii,ij,jk) = dta(ii,ij,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+8) * 0.5_wp * &
& (ref_wgts(kw)%fly_dta(ni+1,nj+2,jk) - ref_wgts(kw)%fly_dta(ni+1,nj ,jk)) & (ref_wgts(kw)%fly_dta(ni+1,nj+2,jk) - ref_wgts(kw)%fly_dta(ni+1,nj ,jk))
END_3D END_3D
END DO END DO
...@@ -1518,8 +1536,10 @@ CONTAINS ...@@ -1518,8 +1536,10 @@ CONTAINS
DO_3D( 0, 0, 0, 0, 1,ipk ) DO_3D( 0, 0, 0, 0, 1,ipk )
ni = ref_wgts(kw)%data_jpi(ji,jj,jn) ni = ref_wgts(kw)%data_jpi(ji,jj,jn)
nj = ref_wgts(kw)%data_jpj(ji,jj,jn) nj = ref_wgts(kw)%data_jpj(ji,jj,jn)
ii = ji - iisht
ij = jj - ijsht
! gradient in the ij direction ! gradient in the ij direction
dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+12) * 0.25_wp * ( & dta(ii,ij,jk) = dta(ii,ij,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+12) * 0.25_wp * ( &
& (ref_wgts(kw)%fly_dta(ni+2,nj+2,jk) - ref_wgts(kw)%fly_dta(ni ,nj+2,jk)) - & & (ref_wgts(kw)%fly_dta(ni+2,nj+2,jk) - ref_wgts(kw)%fly_dta(ni ,nj+2,jk)) - &
& (ref_wgts(kw)%fly_dta(ni+2,nj ,jk) - ref_wgts(kw)%fly_dta(ni ,nj ,jk))) & (ref_wgts(kw)%fly_dta(ni+2,nj ,jk) - ref_wgts(kw)%fly_dta(ni ,nj ,jk)))
END_3D END_3D
......
src/OCE/SBC/geo2ocean.F90
View file @ 3efbcb17
...@@ -57,15 +57,23 @@ CONTAINS ...@@ -57,15 +57,23 @@ CONTAINS
!! ** Purpose : Rotate the Repere: Change vector componantes between !! ** Purpose : Rotate the Repere: Change vector componantes between
!! geographic grid <--> stretched coordinates grid. !! geographic grid <--> stretched coordinates grid.
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pxin, pyin ! vector componantes REAL(wp), DIMENSION(:,:), INTENT(in ) :: pxin, pyin ! vector componantes
CHARACTER(len=1), INTENT(in ) :: cd_type ! define the nature of pt2d array grid-points CHARACTER(len=1), INTENT(in ) :: cd_type ! define the nature of pt2d array grid-points
CHARACTER(len=5), INTENT(in ) :: cdtodo ! type of transpormation: CHARACTER(len=5), INTENT(in ) :: cdtodo ! type of transpormation:
! ! 'en->i' = east-north to i-component ! ! 'en->i' = east-north to i-component
! ! 'en->j' = east-north to j-component ! ! 'en->j' = east-north to j-component
! ! 'ij->e' = (i,j) components to east ! ! 'ij->e' = (i,j) components to east
! ! 'ij->n' = (i,j) components to north ! ! 'ij->n' = (i,j) components to north
REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: prot REAL(wp), DIMENSION(:,:), INTENT( out) :: prot
!
INTEGER :: ipi, ipj, iipi, ijpj
INTEGER :: iisht, ijsht
INTEGER :: ii, ij, ii1, ij1
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
ipi = SIZE(pxin, 1) ; ipj = SIZE(pxin, 2)
iisht = ( jpi - ipi ) / 2 ; ijsht = ( jpj - ipj ) / 2
ii1 = 1 + iisht ; ij1 = 1 + iisht
iipi = ipi + iisht ; ijpj = ipj + ijsht
! !
IF( lmust_init ) THEN ! at 1st call only: set gsin. & gcos. IF( lmust_init ) THEN ! at 1st call only: set gsin. & gcos.
IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*)
...@@ -80,34 +88,50 @@ CONTAINS ...@@ -80,34 +88,50 @@ CONTAINS
! !
CASE( 'en->i' ) ! east-north to i-component CASE( 'en->i' ) ! east-north to i-component
SELECT CASE (cd_type) SELECT CASE (cd_type)
CASE ('T') ; prot(:,:) = pxin(:,:) * gcost(:,:) + pyin(:,:) * gsint(:,:) CASE ('T') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcost(ii1:iipi,ij1:ijpj) &
CASE ('U') ; prot(:,:) = pxin(:,:) * gcosu(:,:) + pyin(:,:) * gsinu(:,:) & + pyin(1:ipi,1:ipj) * gsint(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(:,:) = pxin(:,:) * gcosv(:,:) + pyin(:,:) * gsinv(:,:) CASE ('U') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcosu(ii1:iipi,ij1:ijpj) &
CASE ('F') ; prot(:,:) = pxin(:,:) * gcosf(:,:) + pyin(:,:) * gsinf(:,:) & + pyin(1:ipi,1:ipj) * gsinu(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcosv(ii1:iipi,ij1:ijpj) &
& + pyin(1:ipi,1:ipj) * gsinv(ii1:iipi,ij1:ijpj)
CASE ('F') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcosf(ii1:iipi,ij1:ijpj) &
& + pyin(1:ipi,1:ipj) * gsinf(ii1:iipi,ij1:ijpj)
CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' ) CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' )
END SELECT END SELECT
CASE ('en->j') ! east-north to j-component CASE ('en->j') ! east-north to j-component
SELECT CASE (cd_type) SELECT CASE (cd_type)
CASE ('T') ; prot(:,:) = pyin(:,:) * gcost(:,:) - pxin(:,:) * gsint(:,:) CASE ('T') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcost(ii1:iipi,ij1:ijpj) &
CASE ('U') ; prot(:,:) = pyin(:,:) * gcosu(:,:) - pxin(:,:) * gsinu(:,:) & - pxin(1:ipi,1:ipj) * gsint(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(:,:) = pyin(:,:) * gcosv(:,:) - pxin(:,:) * gsinv(:,:) CASE ('U') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcosu(ii1:iipi,ij1:ijpj) &
CASE ('F') ; prot(:,:) = pyin(:,:) * gcosf(:,:) - pxin(:,:) * gsinf(:,:) & - pxin(1:ipi,1:ipj) * gsinu(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcosv(ii1:iipi,ij1:ijpj) &
& - pxin(1:ipi,1:ipj) * gsinv(ii1:iipi,ij1:ijpj)
CASE ('F') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcosf(ii1:iipi,ij1:ijpj) &
& - pxin(1:ipi,1:ipj) * gsinf(ii1:iipi,ij1:ijpj)
CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' ) CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' )
END SELECT END SELECT
CASE ('ij->e') ! (i,j)-components to east CASE ('ij->e') ! (i,j)-components to east
SELECT CASE (cd_type) SELECT CASE (cd_type)
CASE ('T') ; prot(:,:) = pxin(:,:) * gcost(:,:) - pyin(:,:) * gsint(:,:) CASE ('T') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcost(ii1:iipi,ij1:ijpj) &
CASE ('U') ; prot(:,:) = pxin(:,:) * gcosu(:,:) - pyin(:,:) * gsinu(:,:) & - pyin(1:ipi,1:ipj) * gsint(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(:,:) = pxin(:,:) * gcosv(:,:) - pyin(:,:) * gsinv(:,:) CASE ('U') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcosu(ii1:iipi,ij1:ijpj) &
CASE ('F') ; prot(:,:) = pxin(:,:) * gcosf(:,:) - pyin(:,:) * gsinf(:,:) & - pyin(1:ipi,1:ipj) * gsinu(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcosv(ii1:iipi,ij1:ijpj) &
& - pyin(1:ipi,1:ipj) * gsinv(ii1:iipi,ij1:ijpj)
CASE ('F') ; prot(1:ipi,1:ipj) = pxin(1:ipi,1:ipj) * gcosf(ii1:iipi,ij1:ijpj) &
& - pyin(1:ipi,1:ipj) * gsinf(ii1:iipi,ij1:ijpj)
CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' ) CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' )
END SELECT END SELECT
CASE ('ij->n') ! (i,j)-components to north CASE ('ij->n') ! (i,j)-components to north
SELECT CASE (cd_type) SELECT CASE (cd_type)
CASE ('T') ; prot(:,:) = pyin(:,:) * gcost(:,:) + pxin(:,:) * gsint(:,:) CASE ('T') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcost(ii1:iipi,ij1:ijpj) &
CASE ('U') ; prot(:,:) = pyin(:,:) * gcosu(:,:) + pxin(:,:) * gsinu(:,:) & + pxin(1:ipi,1:ipj) * gsint(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(:,:) = pyin(:,:) * gcosv(:,:) + pxin(:,:) * gsinv(:,:) CASE ('U') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcosu(ii1:iipi,ij1:ijpj) &
CASE ('F') ; prot(:,:) = pyin(:,:) * gcosf(:,:) + pxin(:,:) * gsinf(:,:) & + pxin(1:ipi,1:ipj) * gsinu(ii1:iipi,ij1:ijpj)
CASE ('V') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcosv(ii1:iipi,ij1:ijpj) &
& + pxin(1:ipi,1:ipj) * gsinv(ii1:iipi,ij1:ijpj)
CASE ('F') ; prot(1:ipi,1:ipj) = pyin(1:ipi,1:ipj) * gcosf(ii1:iipi,ij1:ijpj) &
& + pxin(1:ipi,1:ipj) * gsinf(ii1:iipi,ij1:ijpj)
CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' ) CASE DEFAULT ; CALL ctl_stop( 'Only T, U, V and F grid points are coded' )
END SELECT END SELECT
CASE DEFAULT ; CALL ctl_stop( 'rot_rep: Syntax Error in the definition of cdtodo' ) CASE DEFAULT ; CALL ctl_stop( 'rot_rep: Syntax Error in the definition of cdtodo' )
...@@ -286,14 +310,17 @@ CONTAINS ...@@ -286,14 +310,17 @@ CONTAINS
!! ** Method : Change a vector from geocentric to east/north !! ** Method : Change a vector from geocentric to east/north
!! !!
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pxx, pyy, pzz REAL(wp), DIMENSION(:,:), INTENT(in ) :: pxx, pyy, pzz
CHARACTER(len=1) , INTENT(in ) :: cgrid CHARACTER(len=1) , INTENT(in ) :: cgrid
REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pte, ptn REAL(wp), DIMENSION(:,:), INTENT( out) :: pte, ptn
! !
REAL(wp), PARAMETER :: rpi = 3.141592653e0 REAL(wp), PARAMETER :: rpi = 3.141592653e0
REAL(wp), PARAMETER :: rad = rpi / 180.e0 REAL(wp), PARAMETER :: rad = rpi / 180.e0
INTEGER :: ig ! INTEGER :: ig !
INTEGER :: ierr ! local integer INTEGER :: ierr ! local integer
INTEGER :: ipi, ipj, iipi, ijpj
INTEGER :: iisht, ijsht
INTEGER :: ii, ij, ii1, ij1
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( .NOT. ALLOCATED( gsinlon ) ) THEN IF( .NOT. ALLOCATED( gsinlon ) ) THEN
...@@ -345,10 +372,16 @@ CONTAINS ...@@ -345,10 +372,16 @@ CONTAINS
CALL ctl_stop( ctmp1 ) CALL ctl_stop( ctmp1 )
END SELECT END SELECT
! !
pte = - gsinlon(:,:,ig) * pxx + gcoslon(:,:,ig) * pyy ipi = SIZE(pxx, 1) ; ipj = SIZE(pxx, 2)
ptn = - gcoslon(:,:,ig) * gsinlat(:,:,ig) * pxx & iisht = ( jpi - ipi ) / 2 ; ijsht = ( jpj - ipj ) / 2
& - gsinlon(:,:,ig) * gsinlat(:,:,ig) * pyy & ii1 = 1 + iisht ; ij1 = 1 + iisht
& + gcoslat(:,:,ig) * pzz iipi = ipi + iisht ; ijpj = ipj + ijsht
!
pte(1:ipi,1:ipj) = - gsinlon(ii1:iipi,ij1:ijpj,ig) * pxx(1:ipi,1:ipj) &
& + gcoslon(ii1:iipi,ij1:ijpj,ig) * pyy(1:ipi,1:ipj)
ptn(1:ipi,1:ipj) = - gcoslon(ii1:iipi,ij1:ijpj,ig) * gsinlat(ii1:iipi,ij1:ijpj,ig) * pxx(1:ipi,1:ipj) &
& - gsinlon(ii1:iipi,ij1:ijpj,ig) * gsinlat(ii1:iipi,ij1:ijpj,ig) * pyy(1:ipi,1:ipj) &
& + gcoslat(ii1:iipi,ij1:ijpj,ig) * pzz(1:ipi,1:ipj)
! !
END SUBROUTINE geo2oce END SUBROUTINE geo2oce
...@@ -371,6 +404,9 @@ CONTAINS ...@@ -371,6 +404,9 @@ CONTAINS
REAL(wp), PARAMETER :: rad = rpi / 180.e0_wp REAL(wp), PARAMETER :: rad = rpi / 180.e0_wp
INTEGER :: ig ! INTEGER :: ig !
INTEGER :: ierr ! local integer INTEGER :: ierr ! local integer
INTEGER :: ipi, ipj, iipi, ijpj
INTEGER :: iisht, ijsht
INTEGER :: ii, ij, ii1, ij1
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
IF( .NOT. ALLOCATED( gsinlon ) ) THEN IF( .NOT. ALLOCATED( gsinlon ) ) THEN
...@@ -422,9 +458,16 @@ CONTAINS ...@@ -422,9 +458,16 @@ CONTAINS
CALL ctl_stop( ctmp1 ) CALL ctl_stop( ctmp1 )
END SELECT END SELECT
! !
pxx = - gsinlon(:,:,ig) * pte - gcoslon(:,:,ig) * gsinlat(:,:,ig) * ptn ipi = SIZE(pte, 1) ; ipj = SIZE(pte, 2)
pyy = gcoslon(:,:,ig) * pte - gsinlon(:,:,ig) * gsinlat(:,:,ig) * ptn iisht = ( jpi - ipi ) / 2 ; ijsht = ( jpj - ipj ) / 2
pzz = gcoslat(:,:,ig) * ptn ii1 = 1 + iisht ; ij1 = 1 + iisht
iipi = ipi + iisht ; ijpj = ipj + ijsht
!
pxx(1:ipi,1:ipj) = - gsinlon(ii1:iipi,ij1:ijpj,ig) * pte(1:ipi,1:ipj) &
& - gcoslon(ii1:iipi,ij1:ijpj,ig) * gsinlat(ii1:iipi,ij1:ijpj,ig) * ptn(1:ipi,1:ipj)
pyy(1:ipi,1:ipj) = gcoslon(ii1:iipi,ij1:ijpj,ig) * pte(1:ipi,1:ipj) &
& - gsinlon(ii1:iipi,ij1:ijpj,ig) * gsinlat(ii1:iipi,ij1:ijpj,ig) * ptn(1:ipi,1:ipj)
pzz(1:ipi,1:ipj) = gcoslat(ii1:iipi,ij1:ijpj,ig) * ptn(1:ipi,1:ipj)
! !
END SUBROUTINE oce2geo END SUBROUTINE oce2geo
......
src/OCE/SBC/ocealb.F90
View file @ 3efbcb17
...@@ -17,7 +17,9 @@ MODULE ocealb ...@@ -17,7 +17,9 @@ MODULE ocealb
PRIVATE PRIVATE
PUBLIC oce_alb ! routine called by sbccpl PUBLIC oce_alb ! routine called by sbccpl
!! Substitution
# include "do_loop_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: ocealb.F90 10069 2018-08-28 14:12:24Z nicolasmartin $ !! $Id: ocealb.F90 10069 2018-08-28 14:12:24Z nicolasmartin $
...@@ -31,8 +33,8 @@ CONTAINS ...@@ -31,8 +33,8 @@ CONTAINS
!! !!
!! ** Purpose : Computation of the albedo of the ocean !! ** Purpose : Computation of the albedo of the ocean
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
REAL(wp), DIMENSION(:,:), INTENT(out) :: palb_os ! albedo of ocean under overcast sky REAL(wp), DIMENSION(A2D(0)), INTENT(out) :: palb_os ! albedo of ocean under overcast sky
REAL(wp), DIMENSION(:,:), INTENT(out) :: palb_cs ! albedo of ocean under clear sky REAL(wp), DIMENSION(A2D(0)), INTENT(out) :: palb_cs ! albedo of ocean under clear sky
!! !!
REAL(wp) :: zcoef REAL(wp) :: zcoef
REAL(wp) :: rmue = 0.40 ! cosine of local solar altitude REAL(wp) :: rmue = 0.40 ! cosine of local solar altitude
......