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
No results found
Show changes
Hide whitespace changes
Inline Side-by-side
Showing
with 174 additions and 74 deletions
src/OCE/ISF/isfload.F90
View file @ 672ff863
......@@ -26,6 +26,7 @@ MODULE isfload
!
!! * Substitutions
# include "do_loop_substitute.h90"
# include "single_precision_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
......@@ -87,7 +88,7 @@ CONTAINS
#if defined key_qco && key_isf
CALL eos( zts_top(:,:,:), gdept_0(:,:,jk), zrhd(:,:,jk) )
#else
CALL eos( zts_top(:,:,:), gdept(:,:,jk,Kmm), zrhd(:,:,jk) )
CALL eos( zts_top(:,:,:), CASTSP(gdept(:,:,jk,Kmm)), zrhd(:,:,jk) )
#endif
END DO
!
......
src/OCE/ISF/isfparmlt.F90
View file @ 672ff863
......@@ -30,6 +30,7 @@ MODULE isfparmlt
!! * Substitutions
# include "single_precision_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
......@@ -110,9 +111,9 @@ CONTAINS
! compute ptfrz
! 1. ------------Mean freezing point
DO jk = 1,jpk
CALL eos_fzp(ts(:,:,jk,jp_sal,Kmm), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
CALL eos_fzp(CASTSP(ts(:,:,jk,jp_sal,Kmm)), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
END DO
CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
CALL isf_tbl(Kmm, CASTDP(ztfrz3d), ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
!
pqfwf(:,:) = sf_isfpar_fwf(1)%fnow(:,:,1) ! fresh water flux from the isf (fwfisf <0 mean melting) ( > 0 from isf to oce)
pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocean/ice shelf flux assume to be equal to latent heat flux ( > 0 from isf to oce)
......@@ -151,9 +152,9 @@ CONTAINS
!
! 0. ------------Mean freezing point
DO jk = 1,jpk
CALL eos_fzp(ts(:,:,jk,jp_sal,Kmm), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
CALL eos_fzp(CASTSP(ts(:,:,jk,jp_sal,Kmm)), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
END DO
CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
CALL isf_tbl(Kmm, CASTDP(ztfrz3d), ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
!
! 1. ------------Mean temperature
CALL isf_tbl(Kmm, ts(:,:,:,jp_tem,Kmm), ztavg, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
......@@ -204,9 +205,9 @@ CONTAINS
!
! 1. ------------Mean freezing point (needed for heat content flux)
DO jk = 1,jpk
CALL eos_fzp(ts(:,:,jk,jp_sal,Kmm), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
CALL eos_fzp(CASTSP(ts(:,:,jk,jp_sal,Kmm)), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
END DO
CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
CALL isf_tbl(Kmm, CASTDP(ztfrz3d), ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
!
! 2. ------------Scale isf melt pattern with total amount from oasis
! ice shelf 2d map of fwf from isf to oce
......
src/OCE/ISF/isfstp.F90
View file @ 672ff863
......@@ -35,6 +35,7 @@ MODULE isfstp
PUBLIC isf_stp, isf_init, isf_nam ! routine called in sbcmod and divhor
!! * Substitutions
# include "single_precision_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
......@@ -86,7 +87,7 @@ CONTAINS
DO jk = 1, jpk
ze3t(:,:,jk) = e3t(:,:,jk,Kmm)
END DO
CALL isf_tbl_lvl( ht(:,:), ze3t , misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav )
CALL isf_tbl_lvl( CASTSP(ht(:,:)), ze3t , misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav )
#else
CALL isf_tbl_lvl( ht(:,:), e3t(:,:,:,Kmm), misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav )
#endif
......@@ -115,7 +116,7 @@ CONTAINS
DO jk = 1, jpk
ze3t(:,:,jk) = e3t(:,:,jk,Kmm)
END DO
CALL isf_tbl_lvl( ht(:,:), ze3t , misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par )
CALL isf_tbl_lvl( CASTSP(ht(:,:)), ze3t , misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par )
#else
CALL isf_tbl_lvl( ht(:,:), e3t(:,:,:,Kmm), misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par )
#endif
......
src/OCE/ISF/isftbl.F90
View file @ 672ff863
......@@ -44,7 +44,7 @@ CONTAINS
!!-------------------------- IN -------------------------------------
INTEGER , INTENT(in ) :: Kmm ! ocean time level index
CHARACTER(len=1) , INTENT(in ) :: cd_ptin ! point of variable in/out
REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pvarin ! 3d variable to average over the tbl
REAL(dp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pvarin ! 3d variable to average over the tbl
INTEGER, DIMENSION(jpi,jpj) , INTENT(in ) :: ktop ! top level
REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: phtbl ! tbl thickness
!!-------------------------- IN OPTIONAL -----------------------------
......@@ -130,7 +130,7 @@ CONTAINS
INTEGER, DIMENSION(jpi,jpj) , INTENT(in ) :: ktop, kbot ! top and bottom level of the top boundary layer
REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: phtbl, pfrac ! fraction of bottom level to be affected by the tbl
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pe3 ! vertical scale factor
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pvarin ! tbl property to average between ktop, kbot over phtbl
REAL(dp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pvarin ! tbl property to average between ktop, kbot over phtbl
!!--------------------------------------------------------------------
INTEGER :: ji,jj,jk ! loop indices
INTEGER :: ikt, ikb ! top and bottom levels
......
src/OCE/ISF/isfutils.F90
View file @ 672ff863
......@@ -10,12 +10,12 @@ MODULE isfutils
!! isfutils : - read_2dcstdta to read a constant input file with iom_get
!! - debug to print array sum, min, max in ocean.output
!!----------------------------------------------------------------------
USE par_kind
USE iom , ONLY: iom_open, iom_get, iom_close, jpdom_global ! read input file
USE lib_fortran , ONLY: glob_sum, glob_min, glob_max ! compute global value
USE par_oce , ONLY: jpi,jpj,jpk, jpnij, Nis0, Nie0, Njs0, Nje0 ! domain size
USE dom_oce , ONLY: narea ! local domain
USE in_out_manager, ONLY: i8, wp, lwp, numout ! miscelenious
USE in_out_manager, ONLY: lwp, numout ! miscelenious
USE lib_mpp
IMPLICIT NONE
......@@ -28,6 +28,7 @@ MODULE isfutils
PUBLIC read_2dcstdta, debug
# include "single_precision_substitute.h90"
CONTAINS
SUBROUTINE read_2dcstdta(cdfile, cdvar, pvar)
......@@ -70,9 +71,9 @@ CONTAINS
!!--------------------------------------------------------------------
!
! global min/max/sum to check data range and NaN
zsum = glob_sum( 'debug', pvar(:,:) )
zmin = glob_min( 'debug', pvar(:,:) )
zmax = glob_max( 'debug', pvar(:,:) )
zsum =glob_sum( 'debug', CASTDP(pvar(:,:)) )
zmin = glob_min( 'debug', CASTDP(pvar(:,:)) )
zmax = glob_max( 'debug', CASTDP(pvar(:,:)) )
!
! basic check sum to check reproducibility
! TRANSFER function find out the integer corresponding to pvar(i,j) bit pattern
......@@ -112,7 +113,7 @@ CONTAINS
!!
!!-------------------------- IN -------------------------------------
CHARACTER(LEN=*) , INTENT(in ) :: cdtxt
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pvar
REAL(dp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pvar
!!--------------------------------------------------------------------
REAL(wp) :: zmin, zmax, zsum
INTEGER(i8) :: imodd, ip
......
src/OCE/LBC/lbc_lnk_neicoll_generic.h90
View file @ 672ff863
SUBROUTINE lbc_lnk_neicoll_/**/PRECISION( cdname, ptab, cd_nat, psgn, kfld, kfillmode, pfillval, khls, lsend, lrecv, ld4only )
SUBROUTINE lbc_lnk_neicoll_/**/PRECISION( cdname, ptab, cd_nat, psgn, kfld, kfillmode, pfillval, khls, lsend, lrecv, ld4only )
CHARACTER(len=*) , INTENT(in ) :: cdname ! name of the calling subroutine
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
......@@ -270,4 +269,3 @@
ENDIF
!
END SUBROUTINE lbc_lnk_neicoll_/**/PRECISION
src/OCE/LBC/lbc_lnk_pt2pt_generic.h90
View file @ 672ff863
#if ! defined BLOCK_ISEND && ! defined BLOCK_FILL
SUBROUTINE lbc_lnk_pt2pt_/**/PRECISION( cdname, ptab, cd_nat, psgn, kfld, kfillmode, pfillval, khls, lsend, lrecv, ld4only )
CHARACTER(len=*) , INTENT(in ) :: cdname ! name of the calling subroutine
......
src/OCE/LBC/lbc_nfd_ext_generic.h90
View file @ 672ff863
SUBROUTINE lbc_nfd_ext_/**/PRECISION( ptab, cd_nat, psgn, kextj )
SUBROUTINE lbc_nfd_ext_/**/PRECISION( ptab, cd_nat, psgn, kextj )
!!----------------------------------------------------------------------
REAL(PRECISION), DIMENSION(:,1-kextj:),INTENT(inout) :: ptab
CHARACTER(len=1), INTENT(in ) :: cd_nat ! nature of array grid-points
......@@ -117,4 +116,3 @@
ENDIF ! c_NFtype == 'F'
!
END SUBROUTINE lbc_nfd_ext_/**/PRECISION
src/OCE/LBC/lbc_nfd_generic.h90
View file @ 672ff863
SUBROUTINE lbc_nfd_/**/PRECISION( ptab, cd_nat, psgn, khls, kfld )
SUBROUTINE lbc_nfd_/**/PRECISION( ptab, cd_nat, psgn, khls, kfld )
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
REAL(PRECISION), DIMENSION(:), INTENT(in ) :: psgn ! sign used across the north fold boundary
......@@ -388,4 +387,3 @@
END DO ! ipf
!
END SUBROUTINE lbc_nfd_/**/PRECISION
src/OCE/LBC/mpp_lnk_icb_generic.h90
View file @ 672ff863
......@@ -180,4 +180,4 @@
# undef LBCNORTH
# undef PRECISION
# undef SENDROUTINE
# undef RECVROUTINE
# undef RECVROUTINE
src/OCE/LBC/mpp_loc_generic.h90
View file @ 672ff863
!== IN: ptab is an array ==!
!== IN: ptab is an array ==!
# if defined SINGLE_PRECISION
# define ARRAY_TYPE(i,j,k) REAL(sp) , INTENT(in ) :: ARRAY_IN(i,j,k)
# define ARRAY_TYPE(i,j,k) REAL(dp) , INTENT(in ) :: ARRAY_IN(i,j,k)
#if ! defined key_mpi_off
# define MPI_TYPE MPI_2REAL
#endif
......
src/OCE/LBC/mpp_nfd_generic.h90
View file @ 672ff863
SUBROUTINE mpp_nfd_/**/PRECISION( ptab, cd_nat, psgn, kfillmode, pfillval, khls, kfld )
SUBROUTINE mpp_nfd_/**/PRECISION( ptab, cd_nat, psgn, kfillmode, pfillval, khls, kfld )
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
REAL(PRECISION), DIMENSION(:), INTENT(in ) :: psgn ! sign used across the north fold boundary
......@@ -395,4 +394,3 @@
ENDIF ! ln_nnogather
!
END SUBROUTINE mpp_nfd_/**/PRECISION
src/OCE/LBC/mppini.F90
View file @ 672ff863
......@@ -166,6 +166,7 @@ CONTAINS
902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nammpp in configuration namelist' )
!
nn_hls = MAX(1, nn_hls) ! nn_hls must be > 0
IF( nn_hls > 1 ) CALL ctl_warn( 'mpp_init', 'you use nn_hls > 1, this may significantly slow down NEMO performances' )
IF(lwp) THEN
WRITE(numout,*) ' Namelist nammpp'
IF( jpni < 1 .OR. jpnj < 1 ) THEN
......@@ -332,8 +333,8 @@ CONTAINS
nimpp = iimppt(ii,ij)
njmpp = ijmppt(ii,ij)
!
CALL init_doloop ! set start/end indices of do-loop, depending on the halo width value (nn_hls)
CALL init_locglo ! define now functions needed to convert indices from/to global to/from local domains
CALL init_doloop ! set start/end indices of do-loop, depending on the halo width value (nn_hls)
!
IF(lwp) THEN
WRITE(numout,*)
......@@ -1408,9 +1409,46 @@ ENDIF
Njs0 = 1+nn_hls
Nie0 = jpi-nn_hls
Nje0 = jpj-nn_hls
!
Ni_0 = Nie0 - Nis0 + 1
Nj_0 = Nje0 - Njs0 + 1
! Start and end indexes for actual coupling fields on extended grid
Nis0_ext = Nis0
Nie0_ext = Nie0
Njs0_ext = Njs0
Nje0_ext = Nje0
IF (mig(1) == 1 ) THEN
! Drag start column 1 place to the left/west
Nis0_ext=Nis0_ext-1
ENDIF
IF (mig(jpi) == jpiglo ) THEN
! Drag end column 1 place to the right/east
Nie0_ext=Nie0_ext+1
ENDIF
! RSRH we don't adjust anything in the N-S (j) direction
! since the content of the N-fold is catared for by populating values
! using lbc_lnk rather than relying on getting them from the coupler.
! This is rather confused by the fact that jpjglo has a value BIGGER
! than it did at pre 4.2... e.g. for ORCA1 it's set to 333 instead of 332
! which is rather baffling and confuses some dimensioning calculations
! which previously worked fine pre 4.2.
! Set up dimensions for old style coupling exchanges on extended grid
Ni_0_ext = Ni_0
Nj_0_ext = Nj_0
IF (mig(1) == 1 .OR. mig(jpi) == jpiglo ) THEN
! We're at the extreme left or right edge of the grid so need to cater
! for an extra column
Ni_0_ext = Ni_0_ext + 1
ENDIF
!
jpkm1 = jpk-1 ! " "
!
......
src/OCE/LDF/ldfdyn.F90
View file @ 672ff863
......@@ -102,7 +102,7 @@ CONTAINS
!! = 20 = F(i,j) = F(e1,e2) or F(e1^3,e2^3) (lap or bilap case)
!! =-30 => = F(i,j,k) = shape read in 'eddy_viscosity_3D.nc' file
!! = 30 = F(i,j,k) = 2D (case 20) + decrease with depth (case 10)
!! = 31 = F(i,j,k,t) = F(local velocity) ( |u|e /12 laplacian operator
!! = 31 = F(i,j,k,t) = F(local velocity) ( |u|e / 2 laplacian operator
!! or |u|e^3/12 bilaplacian operator )
!! = 32 = F(i,j,k,t) = F(local deformation rate and gridscale) (D and L) (Smagorinsky)
!! ( L^2|D| laplacian operator
......@@ -357,7 +357,7 @@ CONTAINS
!! ** Method : time varying eddy viscosity coefficients:
!!
!! nn_ahm_ijk_t = 31 ahmt, ahmf = F(i,j,k,t) = F(local velocity)
!! ( |u|e /12 or |u|e^3/12 for laplacian or bilaplacian operator )
!! ( |u|e / 2 or |u|e^3/2 for laplacian, u|e^3/12 for bilaplacian operator )
!!
!! nn_ahm_ijk_t = 32 ahmt, ahmf = F(i,j,k,t) = F(local deformation rate and gridscale) (D and L) (Smagorinsky)
!! ( L^2|D| or L^4|D|/8 for laplacian or bilaplacian operator )
......@@ -379,20 +379,20 @@ CONTAINS
!
CASE( 31 ) !== time varying 3D field ==! = F( local velocity )
!
IF( ln_dynldf_lap ) THEN ! laplacian operator : |u| e /12 = |u/144| e
DO jk = 1, jpkm1
IF( ln_dynldf_lap ) THEN ! laplacian operator : |u| e / 2 = |u/ 4| e
DO jk = 1, jpkm1 ! r1_8 = 1 / (2*2 * 2)
DO_2D( 0, 0, 0, 0 )
zu2pv2_ij = uu(ji ,jj ,jk,Kbb) * uu(ji ,jj ,jk,Kbb) + vv(ji ,jj ,jk,Kbb) * vv(ji ,jj ,jk,Kbb)
zu2pv2_ij_m1 = uu(ji-1,jj ,jk,Kbb) * uu(ji-1,jj ,jk,Kbb) + vv(ji ,jj-1,jk,Kbb) * vv(ji ,jj-1,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)
zemax = MAX( e1t(ji,jj) , e2t(ji,jj) )
ahmt(ji,jj,jk) = SQRT( (zu2pv2_ij + zu2pv2_ij_m1) * r1_288 ) * zemax * tmask(ji,jj,jk) ! 288= 12*12 * 2
ahmt(ji,jj,jk) = SQRT( (zu2pv2_ij + zu2pv2_ij_m1) * r1_8 ) * zemax * tmask(ji,jj,jk)
zemax = MAX( e1f(ji,jj) , e2f(ji,jj) )
ahmf(ji,jj,jk) = SQRT( (zu2pv2_ij + zu2pv2_ij_p1) * r1_288 ) * zemax * fmask(ji,jj,jk) ! 288= 12*12 * 2
ahmf(ji,jj,jk) = SQRT( (zu2pv2_ij + zu2pv2_ij_p1) * r1_8 ) * zemax * fmask(ji,jj,jk)
END_2D
END DO
ELSEIF( ln_dynldf_blp ) THEN ! bilaplacian operator : sqrt( |u| e^3 /12 ) = sqrt( |u/144| e ) * e
DO jk = 1, jpkm1
DO jk = 1, jpkm1 ! r1_288 = 1 / (12*12 * 2)
DO_2D( 0, 0, 0, 0 )
zu2pv2_ij = uu(ji ,jj ,jk,Kbb) * uu(ji ,jj ,jk,Kbb) + vv(ji ,jj ,jk,Kbb) * vv(ji ,jj ,jk,Kbb)
zu2pv2_ij_m1 = uu(ji-1,jj ,jk,Kbb) * uu(ji-1,jj ,jk,Kbb) + vv(ji ,jj-1,jk,Kbb) * vv(ji ,jj-1,jk,Kbb)
......
src/OCE/LDF/ldfslp.F90
View file @ 672ff863
......@@ -74,6 +74,7 @@ MODULE ldfslp
!! * Substitutions
# include "do_loop_substitute.h90"
# include "single_precision_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
......@@ -320,8 +321,8 @@ CONTAINS
CALL lbc_lnk( 'ldfslp', uslp , 'U', -1.0_wp , vslp , 'V', -1.0_wp , wslpi, 'W', -1.0_wp, wslpj, 'W', -1.0_wp )
IF(sn_cfctl%l_prtctl) THEN
CALL prt_ctl(tab3d_1=uslp , clinfo1=' slp - u : ', tab3d_2=vslp, clinfo2=' v : ')
CALL prt_ctl(tab3d_1=wslpi, clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ')
CALL prt_ctl(tab3d_1=CASTDP(uslp), clinfo1=' slp - u : ', tab3d_2=CASTDP(vslp), clinfo2=' v : ')
CALL prt_ctl(tab3d_1=CASTDP(wslpi), clinfo1=' slp - wi: ', tab3d_2=CASTDP(wslpj), clinfo2=' wj: ')
ENDIF
!
IF( ln_timing ) CALL timing_stop('ldf_slp')
......
src/OCE/LDF/ldftra.F90
View file @ 672ff863
......@@ -71,6 +71,8 @@ MODULE ldftra
INTEGER , PUBLIC :: nn_aei_ijk_t !: choice of time/space variation of the eiv coeff.
REAL(wp), PUBLIC :: rn_Ue !: lateral diffusive velocity [m/s]
REAL(wp), PUBLIC :: rn_Le !: lateral diffusive length [m]
INTEGER, PUBLIC :: nn_ldfeiv_shape !: shape of bounding coefficient (Treguier et al formulation only)
! ! Flag to control the type of lateral diffusive operator
INTEGER, PARAMETER, PUBLIC :: np_ERROR =-10 ! error in specification of lateral diffusion
......@@ -389,7 +391,7 @@ CONTAINS
!! with a reduction to 0 in vicinity of the Equator
!! nn_aht_ijk_t = 21 ahtu, ahtv = F(i,j, t) = F(growth rate of baroclinic instability)
!!
!! = 31 ahtu, ahtv = F(i,j,k,t) = F(local velocity) ( |u|e /12 laplacian operator
!! = 31 ahtu, ahtv = F(i,j,k,t) = F(local velocity) ( |u|e / 2 laplacian operator
!! or |u|e^3/12 bilaplacian operator )
!!
!! * time varying EIV coefficients: call to ldf_eiv routine
......@@ -440,10 +442,10 @@ CONTAINS
END DO
!
CASE( 31 ) !== time varying 3D field ==! = F( local velocity )
IF( ln_traldf_lap ) THEN ! laplacian operator |u| e /12
IF( ln_traldf_lap ) THEN ! laplacian operator |u| e / 2
DO jk = 1, jpkm1
ahtu(:,:,jk) = ABS( uu(:,:,jk,Kbb) ) * e1u(:,:) * r1_12 ! n.b. uu,vv are masked
ahtv(:,:,jk) = ABS( vv(:,:,jk,Kbb) ) * e2v(:,:) * r1_12
ahtu(:,:,jk) = ABS( uu(:,:,jk,Kbb) ) * e1u(:,:) * r1_2 ! n.b. uu,vv are masked
ahtv(:,:,jk) = ABS( vv(:,:,jk,Kbb) ) * e2v(:,:) * r1_2
END DO
ELSEIF( ln_traldf_blp ) THEN ! bilaplacian operator sqrt( |u| e^3 /12 ) = sqrt( |u| e /12 ) * e
DO jk = 1, jpkm1
......@@ -495,7 +497,8 @@ CONTAINS
REAL(wp) :: zah_max, zUfac ! - scalar
!!
NAMELIST/namtra_eiv/ ln_ldfeiv , ln_ldfeiv_dia, & ! eddy induced velocity (eiv)
& nn_aei_ijk_t, rn_Ue, rn_Le ! eiv coefficient
& nn_aei_ijk_t, rn_Ue, rn_Le, & ! eiv coefficient
& nn_ldfeiv_shape
!!----------------------------------------------------------------------
!
IF(lwp) THEN ! control print
......@@ -588,6 +591,7 @@ CONTAINS
IF(lwp) WRITE(numout,*) ' ==>>> eddy induced velocity coef. = F( latitude, longitude, time )'
IF(lwp) WRITE(numout,*) ' = F( growth rate of baroclinic instability )'
IF(lwp) WRITE(numout,*) ' maximum allowed value: aei0 = ', aei0, ' m2/s'
IF(lwp) WRITE(numout,*) ' shape of bounding coefficient : ',nn_ldfeiv_shape
!
l_ldfeiv_time = .TRUE. ! will be calculated by call to ldf_tra routine in step.F90
!
......@@ -636,14 +640,21 @@ CONTAINS
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: paeiu, paeiv ! eiv coefficient [m2/s]
!
INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zfw, ze3w, zn2, z1_f20, zzaei ! local scalars
REAL(wp), DIMENSION(jpi,jpj) :: zn, zah, zhw, zRo, zaeiw ! 2D workspace
REAL(wp) :: zfw, ze3w, zn2, z1_f20, zzaei, z2_3 ! local scalars
REAL(wp), DIMENSION(jpi,jpj) :: zn, zah, zhw, zRo, zRo_lim, zTclinic_recip, zaeiw, zratio ! 2D workspace
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmodslp ! 3D workspace
!!----------------------------------------------------------------------
!
zn (:,:) = 0._wp ! Local initialization
zmodslp(:,:,:) = 0._wp
zhw(:,:) = 5._wp
zah(:,:) = 0._wp
zRo(:,:) = 0._wp
zRo_lim(:,:) = 0._wp
zTclinic_recip(:,:) = 0._wp
zratio(:,:) = 0._wp
zaeiw(:,:) = 0._wp
! ! Compute lateral diffusive coefficient at T-point
IF( ln_traldf_triad ) THEN
DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpk )
......@@ -670,8 +681,9 @@ CONTAINS
! eddies using the isopycnal slopes calculated in ldfslp.F :
! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w))
ze3w = e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)
zah(ji,jj) = zah(ji,jj) + zn2 * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk) &
& + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) * ze3w
zmodslp(ji,jj,jk) = wslpi(ji,jj,jk) * wslpi(ji,jj,jk) &
& + wslpj(ji,jj,jk) * wslpj(ji,jj,jk)
zah(ji,jj) = zah(ji,jj) + zn2 * zmodslp(ji,jj,jk) * ze3w
zhw(ji,jj) = zhw(ji,jj) + ze3w
END_3D
ENDIF
......@@ -679,17 +691,69 @@ CONTAINS
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zfw = MAX( ABS( 2. * omega * SIN( rad * gphit(ji,jj) ) ) , 1.e-10 )
! Rossby radius at w-point taken betwenn 2 km and 40km
zRo(ji,jj) = MAX( 2.e3 , MIN( .4 * zn(ji,jj) / zfw, 40.e3 ) )
zRo(ji,jj) = .4 * zn(ji,jj) / zfw
zRo_lim(ji,jj) = MAX( 2.e3 , MIN( zRo(ji,jj), 40.e3 ) )
! Compute aeiw by multiplying Ro^2 and T^-1
zaeiw(ji,jj) = zRo(ji,jj) * zRo(ji,jj) * SQRT( zah(ji,jj) / zhw(ji,jj) ) * tmask(ji,jj,1)
zTclinic_recip(ji,jj) = SQRT( MAX(zah(ji,jj),0._wp) / zhw(ji,jj) ) * tmask(ji,jj,1)
zaeiw(ji,jj) = zRo_lim(ji,jj) * zRo_lim(ji,jj) * zTclinic_recip(ji,jj)
END_2D
IF( iom_use('N_2d') ) CALL iom_put('N_2d',zn(:,:)/zhw(:,:))
IF( iom_use('modslp') ) CALL iom_put('modslp',SQRT(zmodslp(:,:,:)) )
CALL iom_put('RossRad',zRo)
CALL iom_put('RossRadlim',zRo_lim)
CALL iom_put('Tclinic_recip',zTclinic_recip)
! !== Bound on eiv coeff. ==!
z1_f20 = 1._wp / ( 2._wp * omega * sin( rad * 20._wp ) )
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zzaei = MIN( 1._wp, ABS( ff_t(ji,jj) * z1_f20 ) ) * zaeiw(ji,jj) ! tropical decrease
zaeiw(ji,jj) = MIN( zzaei , paei0 ) ! Max value = paei0
END_2D
z2_3 = 2._wp/3._wp
SELECT CASE(nn_ldfeiv_shape)
CASE(0) !! Standard shape applied - decrease in tropics and cap.
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zzaei = MIN( 1._wp, ABS( ff_t(ji,jj) * z1_f20 ) ) * zaeiw(ji,jj) ! tropical decrease
zaeiw(ji,jj) = MIN( zzaei, paei0 )
END_2D
CASE(1) !! Abrupt cut-off on Rossby radius:
! JD : modifications here to introduce scaling by local rossby radius of deformation vs local grid scale
! arbitrary decision that GM is de-activated if local rossy radius larger than 2 times local grid scale
! based on Hallberg (2013)
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
IF ( zRo(ji,jj) >= ( 2._wp * MIN( e1t(ji,jj), e2t(ji,jj) ) ) ) THEN
! TODO : use a version of zRo that integrates over a few time steps ?
zaeiw(ji,jj) = 0._wp
ELSE
zaeiw(ji,jj) = MIN( zaeiw(ji,jj), paei0 )
ENDIF
END_2D
CASE(2) !! Rossby radius ramp type 1:
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zratio(ji,jj) = zRo(ji,jj)/MIN(e1t(ji,jj),e2t(ji,jj))
zaeiw(ji,jj) = MIN( zaeiw(ji,jj), MAX( 0._wp, MIN( 1._wp, z2_3*(2._wp - zratio(ji,jj)) ) ) * paei0 )
END_2D
CALL iom_put('RR_GS',zratio)
CASE(3) !! Rossby radius ramp type 2:
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zratio(ji,jj) = MIN(e1t(ji,jj),e2t(ji,jj))/zRo(ji,jj)
zaeiw(ji,jj) = MIN( zaeiw(ji,jj), MAX( 0._wp, MIN( 1._wp, z2_3*( zratio(ji,jj) - 0.5_wp ) ) ) * paei0 )
END_2D
CASE(4) !! bathymetry ramp:
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zaeiw(ji,jj) = MIN( zaeiw(ji,jj), MAX( 0._wp, MIN( 1._wp, 0.001*(ht_0(ji,jj) - 2000._wp) ) ) * paei0 )
END_2D
CASE(5) !! Rossby radius ramp type 1 applied to Treguier et al coefficient rather than cap:
!! Note the ramp is RR/GS=[2.0,1.0] (not [2.0,0.5] as for cases 2,3) and we ramp up
!! to 5% of the Treguier et al coefficient, aiming for peak values of around 100m2/s
!! at high latitudes rather than 2000m2/s which is what you get in eORCA025 with an
!! uncapped coefficient.
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zratio(ji,jj) = zRo(ji,jj)/MIN(e1t(ji,jj),e2t(ji,jj))
zaeiw(ji,jj) = MAX( 0._wp, MIN( 1._wp, 2._wp - zratio(ji,jj) ) ) * 0.05 * zaeiw(ji,jj)
zaeiw(ji,jj) = MIN( zaeiw(ji,jj), paei0 )
END_2D
CALL iom_put('RR_GS',zratio)
CASE DEFAULT
CALL ctl_stop('ldf_eiv: Unrecognised option for nn_ldfeiv_shape.')
END SELECT
IF( nn_hls == 1 ) CALL lbc_lnk( 'ldftra', zaeiw(:,:), 'W', 1.0_wp ) ! lateral boundary condition
!
DO_2D( 0, 0, 0, 0 )
......@@ -730,7 +794,7 @@ CONTAINS
CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
! TEMP: [tiling] Can be A2D(nn_hls) after all lbc_lnks removed in the nn_hls = 2 case in tra_adv_fct
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu ! in : 3 ocean transport components [m3/s]
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pv ! out: 3 ocean transport components [m3/s]
REAL(dp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pv ! out: 3 ocean transport components [m3/s]
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pw ! increased by the eiv [m3/s]
!!
INTEGER :: ji, jj, jk ! dummy loop indices
......@@ -867,7 +931,7 @@ CONTAINS
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' ) ) CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d )
IF( iom_use( 'sophteiv' ) ) CALL dia_ptr_hst( jp_tem, 'eiv', 0.5_wp * zw3d )
!
zztmp = 0.5_wp * 0.5
IF( iom_use('ueiv_salttr') .OR. iom_use('ueiv_salttr3d')) THEN
......@@ -891,7 +955,7 @@ CONTAINS
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' ) ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d )
IF( iom_use( 'sopsteiv' ) ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5_wp * zw3d )
!
!
END SUBROUTINE ldf_eiv_dia
......
src/OCE/OBS/README.rst deleted 120000 → 0
View file @ c29e0862
../../../doc/rst/source/da.rst
\ No newline at end of file
src/OCE/OBS/ddatetoymdhms.h90
View file @ 672ff863
!!----------------------------------------------------------------------
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: ddatetoymdhms.h90 13226 2020-07-02 14:24:31Z orioltp $
!! Software governed by the CeCILL license (see ./LICENSE)
......@@ -20,7 +20,7 @@
!! * Modules used
!! * Arguments
real(dp), INTENT(IN) :: ddate
real(wp), INTENT(IN) :: ddate
INTEGER, INTENT(OUT) :: kyea
INTEGER, INTENT(OUT) :: kmon
INTEGER, INTENT(OUT) :: kday
......
src/OCE/OBS/diaobs.F90
View file @ 672ff863
......@@ -161,8 +161,8 @@ CONTAINS
LOGICAL, DIMENSION(:), ALLOCATABLE :: llvar ! Logical for profile variable read
LOGICAL, DIMENSION(jpmaxnfiles) :: lmask ! Used for finding number of sstbias files
!
REAL(dp) :: rn_dobsini ! Obs window start date YYYYMMDD.HHMMSS
REAL(dp) :: rn_dobsend ! Obs window end date YYYYMMDD.HHMMSS
REAL(wp) :: rn_dobsini ! Obs window start date YYYYMMDD.HHMMSS
REAL(wp) :: rn_dobsend ! Obs window end date YYYYMMDD.HHMMSS
REAL(wp) :: ztype_avglamscl ! Local version of rn_*_avglamscl
REAL(wp) :: ztype_avgphiscl ! Local version of rn_*_avgphiscl
REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zglam ! Model longitudes for profile variables
......@@ -404,7 +404,7 @@ CONTAINS
ltype_night = .FALSE.
ENDIF
CALL obs_setinterpopts( nsurftypes, jtype, TRIM(cobstypessurf(jtype)), &
CALL obs_setinterpopts( nsurftypes, jtype, cobstypessurf(jtype), &
& nn_2dint_default, n2dint_type, &
& ztype_avglamscl, ztype_avgphiscl, &
& ltype_fp_indegs, ltype_night, &
......@@ -426,6 +426,9 @@ CONTAINS
ENDIF
!
IF( ln_grid_global ) THEN
IF( jpnij < jpni * jpnj ) THEN
CALL ctl_stop( 'STOP', 'dia_obs_init: ln_grid_global=T is not available when land subdomains are suppressed' )
END IF
CALL ctl_warn( 'dia_obs_init: ln_grid_global=T may cause memory issues when used with a large number of processors' )
ENDIF
!
......@@ -902,8 +905,8 @@ CONTAINS
IMPLICIT NONE
!! * Arguments
REAL(KIND=dp), INTENT(OUT) :: ddobs ! Date in YYYYMMDD.HHMMSS
INTEGER :: kstp
REAL(KIND=wp), INTENT(OUT) :: ddobs ! Date in YYYYMMDD.HHMMSS
INTEGER, INTENT(IN) :: kstp
!! * Local declarations
INTEGER :: iyea ! date - (year, month, day, hour, minute)
......@@ -986,7 +989,7 @@ CONTAINS
IMPLICIT NONE
!! * Arguments
REAL(KIND=dp), INTENT(OUT) :: ddobsini ! Initial date in YYYYMMDD.HHMMSS
REAL(KIND=wp), INTENT(OUT) :: ddobsini ! Initial date in YYYYMMDD.HHMMSS
CALL calc_date( nit000 - 1, ddobsini )
......@@ -1013,7 +1016,7 @@ CONTAINS
IMPLICIT NONE
!! * Arguments
REAL(dp), INTENT(OUT) :: ddobsfin ! Final date in YYYYMMDD.HHMMSS
REAL(wp), INTENT(OUT) :: ddobsfin ! Final date in YYYYMMDD.HHMMSS
CALL calc_date( nitend, ddobsfin )
......@@ -1074,7 +1077,7 @@ CONTAINS
& ravgphiscl_type !N/S diameter of obs footprint for this type
LOGICAL, INTENT(IN) :: lfp_indegs_type !T=> footprint in degrees, F=> in metres
LOGICAL, INTENT(IN) :: lavnight_type !T=> obs represent night time average
CHARACTER(len=8), INTENT(IN) :: ctypein
CHARACTER(len=lca), INTENT(IN) :: ctypein
INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: &
& n2dint
......
src/OCE/OBS/find_obs_proc.h90
View file @ 672ff863
!!----------------------------------------------------------------------
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: find_obs_proc.h90 13286 2020-07-09 15:48:29Z smasson $
!! Software governed by the CeCILL license (see ./LICENSE)
......