Changeset 14053 for NEMO/trunk/src/OCE
- Timestamp:
- 2020-12-03T14:48:38+01:00 (3 years ago)
- Location:
- NEMO/trunk/src/OCE
- Files:
-
- 1 deleted
- 36 edited
- 1 copied
-
DIA/diawri.F90 (modified) (7 diffs)
-
DOM/dom_oce.F90 (modified) (7 diffs)
-
DOM/domain.F90 (modified) (18 diffs)
-
DOM/dommsk.F90 (modified) (1 diff)
-
DOM/domqco.F90 (modified) (11 diffs)
-
DOM/domvvl.F90 (modified) (3 diffs)
-
DOM/domzgr_substitute.h90 (modified) (2 diffs)
-
DOM/istate.F90 (modified) (5 diffs)
-
DOM/phycst.F90 (modified) (1 diff)
-
DYN/dynadv.F90 (modified) (1 diff)
-
DYN/dynatf_qco.F90 (modified) (4 diffs)
-
DYN/dynldf_lap_blp.F90 (modified) (5 diffs)
-
DYN/dynspg_ts.F90 (modified) (16 diffs)
-
DYN/dynvor.F90 (modified) (31 diffs)
-
DYN/sshwzv.F90 (modified) (6 diffs)
-
IOM/iom.F90 (modified) (5 diffs)
-
IOM/restart.F90 (modified) (5 diffs)
-
ISF/isfcpl.F90 (modified) (10 diffs)
-
ISF/isfdynatf.F90 (modified) (1 diff)
-
ISF/isfrst.F90 (modified) (4 diffs)
-
LBC/mppini.F90 (modified) (1 diff)
-
LDF/ldfdyn.F90 (modified) (6 diffs)
-
SBC/sbcapr.F90 (modified) (1 diff)
-
SBC/sbcice_cice.F90 (modified) (2 diffs)
-
SBC/sbcmod.F90 (modified) (2 diffs)
-
SBC/sbcrnf.F90 (modified) (2 diffs)
-
TRA/traatf_qco.F90 (modified) (4 diffs)
-
TRA/traqsr.F90 (modified) (2 diffs)
-
TRA/trasbc.F90 (modified) (3 diffs)
-
USR/usrdef_istate.F90 (modified) (6 diffs)
-
ZDF/zdfddm.F90 (modified) (1 diff)
-
nemogcm.F90 (modified) (11 diffs)
-
oce.F90 (modified) (3 diffs)
-
step.F90 (modified) (1 diff)
-
step_oce.F90 (modified) (7 diffs)
-
stpMLF.F90 (deleted)
-
stpctl.F90 (modified) (3 diffs)
-
stpmlf.F90 (copied) (copied from NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/OCE/stpmlf.F90)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/DIA/diawri.F90
r14045 r14053 19 19 !! 3.7 ! 2014-01 (G. Madec) remove eddy induced velocity from no-IOM output 20 20 !! ! change name of output variables in dia_wri_state 21 !! 4.0 ! 2020-10 (A. Nasser, S. Techene) add diagnostic for SWE 21 22 !!---------------------------------------------------------------------- 22 23 … … 119 120 INTEGER :: ji, jj, jk ! dummy loop indices 120 121 INTEGER :: ikbot ! local integer 121 REAL(wp):: ze3122 122 REAL(wp):: zztmp , zztmpx ! local scalar 123 123 REAL(wp):: zztmp2, zztmpy ! - - 124 REAL(wp):: ze3 124 125 REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace 125 126 REAL(wp), DIMENSION(jpi,jpj,jpk) :: z3d ! 3D workspace … … 138 139 CALL iom_put("e3u_0", e3u_0(:,:,:) ) 139 140 CALL iom_put("e3v_0", e3v_0(:,:,:) ) 141 CALL iom_put("e3f_0", e3f_0(:,:,:) ) 140 142 ! 141 143 IF ( iom_use("e3t") .OR. iom_use("e3tdef") ) THEN ! time-varying e3t … … 164 166 CALL iom_put( "e3w" , z3d(:,:,:) ) 165 167 ENDIF 168 IF ( iom_use("e3f") ) THEN ! time-varying e3f caution here at Kaa 169 DO jk = 1, jpk 170 z3d(:,:,jk) = e3f(:,:,jk) 171 END DO 172 CALL iom_put( "e3f" , z3d(:,:,:) ) 173 ENDIF 166 174 167 175 IF( ll_wd ) THEN ! sea surface height (brought back to the reference used for wetting and drying) 168 CALL iom_put( "ssh" , (ssh(:,:,Kmm)+ssh_ref)* tmask(:,:,1) )176 CALL iom_put( "ssh" , (ssh(:,:,Kmm)+ssh_ref)*ssmask(:,:) ) 169 177 ELSE 170 178 CALL iom_put( "ssh" , ssh(:,:,Kmm) ) ! sea surface height 171 179 ENDIF 172 180 173 IF( iom_use("wetdep") ) & ! wet depth 174 CALL iom_put( "wetdep" , ht_0(:,:) + ssh(:,:,Kmm) ) 181 IF( iom_use("wetdep") ) CALL iom_put( "wetdep" , ht_0(:,:) + ssh(:,:,Kmm) ) ! wet depth 182 183 #if defined key_qco 184 IF( iom_use("ht") ) CALL iom_put( "ht" , ht(:,:) ) ! water column at t-point 185 IF( iom_use("hu") ) CALL iom_put( "hu" , hu(:,:,Kmm) ) ! water column at u-point 186 IF( iom_use("hv") ) CALL iom_put( "hv" , hv(:,:,Kmm) ) ! water column at v-point 187 IF( iom_use("hf") ) CALL iom_put( "hf" , hf_0(:,:)*( 1._wp + r3f(:,:) ) ) ! water column at f-point (caution here at Naa) 188 #endif 175 189 176 190 CALL iom_put( "toce", ts(:,:,:,jp_tem,Kmm) ) ! 3D temperature … … 326 340 ENDIF 327 341 ! 342 IF ( iom_use("sKE") ) THEN ! surface kinetic energy at T point 343 z2d(:,:) = 0._wp 344 DO_2D( 0, 0, 0, 0 ) 345 z2d(ji,jj) = 0.25_wp * ( uu(ji ,jj,1,Kmm) * uu(ji ,jj,1,Kmm) * e1e2u(ji ,jj) * e3u(ji ,jj,1,Kmm) & 346 & + uu(ji-1,jj,1,Kmm) * uu(ji-1,jj,1,Kmm) * e1e2u(ji-1,jj) * e3u(ji-1,jj,1,Kmm) & 347 & + vv(ji,jj ,1,Kmm) * vv(ji,jj ,1,Kmm) * e1e2v(ji,jj ) * e3v(ji,jj ,1,Kmm) & 348 & + vv(ji,jj-1,1,Kmm) * vv(ji,jj-1,1,Kmm) * e1e2v(ji,jj-1) * e3v(ji,jj-1,1,Kmm) ) & 349 & * r1_e1e2t(ji,jj) / e3t(ji,jj,1,Kmm) * ssmask(ji,jj) 350 END_2D 351 CALL lbc_lnk( 'diawri', z2d, 'T', 1. ) 352 IF ( iom_use("sKE" ) ) CALL iom_put( "sKE" , z2d ) 353 ENDIF 354 ! 355 IF ( iom_use("sKEf") ) THEN ! surface kinetic energy at F point 356 z2d(:,:) = 0._wp ! CAUTION : only valid in SWE, not with bathymetry 357 DO_2D( 0, 0, 0, 0 ) 358 z2d(ji,jj) = 0.25_wp * ( uu(ji,jj ,1,Kmm) * uu(ji,jj ,1,Kmm) * e1e2u(ji,jj ) * e3u(ji,jj ,1,Kmm) & 359 & + uu(ji,jj+1,1,Kmm) * uu(ji,jj+1,1,Kmm) * e1e2u(ji,jj+1) * e3u(ji,jj+1,1,Kmm) & 360 & + vv(ji ,jj,1,Kmm) * vv(ji,jj ,1,Kmm) * e1e2v(ji ,jj) * e3v(ji ,jj,1,Kmm) & 361 & + vv(ji+1,jj,1,Kmm) * vv(ji+1,jj,1,Kmm) * e1e2v(ji+1,jj) * e3v(ji+1,jj,1,Kmm) ) & 362 & * r1_e1e2f(ji,jj) / e3f(ji,jj,1) * ssfmask(ji,jj) 363 END_2D 364 CALL lbc_lnk( 'diawri', z2d, 'F', 1. ) 365 CALL iom_put( "sKEf", z2d ) 366 ENDIF 367 ! 328 368 CALL iom_put( "hdiv", hdiv ) ! Horizontal divergence 329 369 … … 425 465 426 466 IF (ln_dia25h) CALL dia_25h( kt, Kmm ) ! 25h averaging 467 468 ! Output of vorticity terms 469 IF ( iom_use("relvor") .OR. iom_use("plavor") .OR. & 470 & iom_use("relpotvor") .OR. iom_use("abspotvor") .OR. & 471 & iom_use("Ens") ) THEN 472 ! 473 z2d(:,:) = 0._wp 474 ze3 = 0._wp 475 DO_2D( 1, 0, 1, 0 ) 476 z2d(ji,jj) = ( e2v(ji+1,jj ) * vv(ji+1,jj ,1,Kmm) - e2v(ji,jj) * vv(ji,jj,1,Kmm) & 477 & - e1u(ji ,jj+1) * uu(ji ,jj+1,1,Kmm) + e1u(ji,jj) * uu(ji,jj,1,Kmm) ) * r1_e1e2f(ji,jj) 478 END_2D 479 CALL lbc_lnk( 'diawri', z2d, 'F', 1. ) 480 CALL iom_put( "relvor", z2d ) ! relative vorticity ( zeta ) 481 ! 482 CALL iom_put( "plavor", ff_f ) ! planetary vorticity ( f ) 483 ! 484 DO_2D( 1, 0, 1, 0 ) 485 ze3 = ( e3t(ji,jj+1,1,Kmm) * e1e2t(ji,jj+1) + e3t(ji+1,jj+1,1,Kmm) * e1e2t(ji+1,jj+1) & 486 & + e3t(ji,jj ,1,Kmm) * e1e2t(ji,jj ) + e3t(ji+1,jj ,1,Kmm) * e1e2t(ji+1,jj ) ) * r1_e1e2f(ji,jj) 487 IF( ze3 /= 0._wp ) THEN ; ze3 = 4._wp / ze3 488 ELSE ; ze3 = 0._wp 489 ENDIF 490 z2d(ji,jj) = ze3 * z2d(ji,jj) 491 END_2D 492 CALL lbc_lnk( 'diawri', z2d, 'F', 1. ) 493 CALL iom_put( "relpotvor", z2d ) ! relative potential vorticity (zeta/h) 494 ! 495 DO_2D( 1, 0, 1, 0 ) 496 ze3 = ( e3t(ji,jj+1,1,Kmm) * e1e2t(ji,jj+1) + e3t(ji+1,jj+1,1,Kmm) * e1e2t(ji+1,jj+1) & 497 & + e3t(ji,jj ,1,Kmm) * e1e2t(ji,jj ) + e3t(ji+1,jj ,1,Kmm) * e1e2t(ji+1,jj ) ) * r1_e1e2f(ji,jj) 498 IF( ze3 /= 0._wp ) THEN ; ze3 = 4._wp / ze3 499 ELSE ; ze3 = 0._wp 500 ENDIF 501 z2d(ji,jj) = ze3 * ff_f(ji,jj) + z2d(ji,jj) 502 END_2D 503 CALL lbc_lnk( 'diawri', z2d, 'F', 1. ) 504 CALL iom_put( "abspotvor", z2d ) ! absolute potential vorticity ( q ) 505 ! 506 DO_2D( 1, 0, 1, 0 ) 507 z2d(ji,jj) = 0.5_wp * z2d(ji,jj) * z2d(ji,jj) 508 END_2D 509 CALL lbc_lnk( 'diawri', z2d, 'F', 1. ) 510 CALL iom_put( "Ens", z2d ) ! potential enstrophy ( 1/2*q2 ) 511 ! 512 ENDIF 427 513 428 514 IF( ln_timing ) CALL timing_stop('dia_wri') … … 998 1084 !! 999 1085 INTEGER :: inum, jk 1000 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, zgdept ! 3D workspace !!st patch to usesubstitution1086 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, zgdept ! 3D workspace for qco substitution 1001 1087 !!---------------------------------------------------------------------- 1002 1088 ! -
NEMO/trunk/src/OCE/DOM/dom_oce.F90
r13982 r14053 131 131 ! 132 132 REAL(wp), PUBLIC, ALLOCATABLE, SAVE , DIMENSION(:,:) :: e1e2t , r1_e1e2t !: associated metrics at t-point 133 REAL(wp), PUBLIC, ALLOCATABLE, SAVE , DIMENSION(:,:) :: e1e2u , e2_e1u, r1_e1e2u!: associated metrics at u-point134 REAL(wp), PUBLIC, ALLOCATABLE, SAVE , DIMENSION(:,:) :: e1e2v , e1_e2v, r1_e1e2v!: associated metrics at v-point133 REAL(wp), PUBLIC, ALLOCATABLE, SAVE , DIMENSION(:,:) :: e1e2u , r1_e1e2u , e2_e1u !: associated metrics at u-point 134 REAL(wp), PUBLIC, ALLOCATABLE, SAVE , DIMENSION(:,:) :: e1e2v , r1_e1e2v , e1_e2v !: associated metrics at v-point 135 135 REAL(wp), PUBLIC, ALLOCATABLE, SAVE , DIMENSION(:,:) :: e1e2f , r1_e1e2f !: associated metrics at f-point 136 136 ! … … 162 162 163 163 ! ! reference depths of cells 164 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdept_0 !: t- depth [m]165 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdepw_0 !: w- depth [m]166 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gde3w_0 !: w- depth (sum of e3w) [m]164 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdept_0 !: t- depth [m] 165 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdepw_0 !: w- depth [m] 166 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gde3w_0 !: w- depth (sum of e3w) [m] 167 167 ! ! time-dependent depths of cells 168 168 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: gdept, gdepw … … 205 205 206 206 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ssmask, ssumask, ssvmask, ssfmask !: surface mask at T-,U-, V- and F-pts 207 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: tmask, umask, vmask, wmask, fmask !: land/ocean mask at T-, U-, V-, W- and F-pts 208 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: wumask, wvmask !: land/ocean mask at WT-, WU- and WV-pts 209 207 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: tmask, umask, vmask, wmask, fmask !: land/ocean mask at T-, U-, V-, W- and F-pts 208 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: wumask, wvmask !: land/ocean mask at WU- and WV-pts 209 #if defined key_qco 210 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: fe3mask !: land/ocean mask at F-pts for qco 211 #endif 210 212 !!---------------------------------------------------------------------- 211 213 !! calendar variables … … 306 308 & e3w_0(jpi,jpj,jpk) , e3uw_0(jpi,jpj,jpk) , e3vw_0(jpi,jpj,jpk) , STAT=ierr(ii) ) 307 309 ! 308 #if ! defined key_qco 310 #if defined key_qco 311 ii = ii+1 312 ALLOCATE( r3t (jpi,jpj,jpt) , r3u (jpi,jpj,jpt) , r3v (jpi,jpj,jpt) , r3f (jpi,jpj) , & 313 & r3t_f(jpi,jpj) , r3u_f(jpi,jpj) , r3v_f(jpi,jpj) , STAT=ierr(ii) ) 314 #else 309 315 ii = ii+1 310 316 ALLOCATE( e3t(jpi,jpj,jpk,jpt) , e3u (jpi,jpj,jpk,jpt) , e3v (jpi,jpj,jpk,jpt) , e3f(jpi,jpj,jpk) , & … … 313 319 ! 314 320 ii = ii+1 315 ALLOCATE( r3t (jpi,jpj,jpt) , r3u (jpi,jpj,jpt) , r3v (jpi,jpj,jpt) , r3f (jpi,jpj) , &316 & r3t_f(jpi,jpj) , r3u_f(jpi,jpj) , r3v_f(jpi,jpj) , STAT=ierr(ii) )317 !318 ii = ii+1319 321 ALLOCATE( ht_0(jpi,jpj) , hu_0(jpi,jpj) , hv_0(jpi,jpj) , hf_0(jpi,jpj) , & 320 322 & r1_ht_0(jpi,jpj) , r1_hu_0(jpi,jpj) , r1_hv_0(jpi,jpj), r1_hf_0(jpi,jpj) , STAT=ierr(ii) ) … … 323 325 ii = ii+1 324 326 ALLOCATE( ht (jpi,jpj) , hu (jpi,jpj,jpt), hv (jpi,jpj,jpt) , & 325 & r1_hu (jpi,jpj,jpt), r1_hv (jpi,jpj,jpt) , STAT=ierr(ii) )326 #else327 ii = ii+1328 ALLOCATE( hu (jpi,jpj,jpt), hv (jpi,jpj,jpt) , &329 327 & r1_hu (jpi,jpj,jpt), r1_hv (jpi,jpj,jpt) , STAT=ierr(ii) ) 330 328 #endif … … 350 348 ii = ii+1 351 349 ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(ii) ) 350 #if defined key_qco 351 ! 352 ii = ii+1 353 ALLOCATE( fe3mask(jpi,jpj,jpk) , STAT=ierr(ii) ) 354 #endif 352 355 ! 353 356 dom_oce_alloc = MAXVAL(ierr) -
NEMO/trunk/src/OCE/DOM/domain.F90
r13982 r14053 15 15 !! 3.7 ! 2015-11 (G. Madec, A. Coward) time varying zgr by default 16 16 !! 4.0 ! 2016-10 (G. Madec, S. Flavoni) domain configuration / user defined interface 17 !! 4. x ! 2020-02 (G. Madec, S. Techene)introduce ssh to h0 ratio17 !! 4.1 ! 2020-02 (G. Madec, S. Techene) introduce ssh to h0 ratio 18 18 !!---------------------------------------------------------------------- 19 19 … … 28 28 USE oce ! ocean variables 29 29 USE dom_oce ! domain: ocean 30 #if defined key_qco 31 USE domqco ! quasi-eulerian 32 #else 33 USE domvvl ! variable volume 34 #endif 35 USE sshwzv , ONLY : ssh_init_rst ! set initial ssh 30 36 USE sbc_oce ! surface boundary condition: ocean 31 37 USE trc_oce ! shared ocean & passive tracers variab … … 35 41 USE dommsk ! domain: set the mask system 36 42 USE domwri ! domain: write the meshmask file 37 #if ! defined key_qco38 USE domvvl ! variable volume39 #else40 USE domqco ! variable volume41 #endif42 43 USE c1d ! 1D configuration 43 44 USE dyncor_c1d ! 1D configuration: Coriolis term (cor_c1d routine) 44 USE wet_dry , ONLY : ll_wd45 USE closea , ONLY : dom_clo ! closed seas45 USE wet_dry , ONLY : ll_wd ! wet & drying flag 46 USE closea , ONLY : dom_clo ! closed seas routine 46 47 ! 47 48 USE prtctl ! Print control (prt_ctl_info routine) … … 50 51 USE lbclnk ! ocean lateral boundary condition (or mpp link) 51 52 USE lib_mpp ! distributed memory computing library 53 USE restart ! only for lrst_oce 52 54 53 55 IMPLICIT NONE … … 58 60 PUBLIC dom_tile ! called by step.F90 59 61 62 !! * Substitutions 63 # include "do_loop_substitute.h90" 60 64 !!------------------------------------------------------------------------- 61 65 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 84 88 INTEGER :: ji, jj, jk, jt ! dummy loop indices 85 89 INTEGER :: iconf = 0 ! local integers 90 REAL(wp):: zrdt 86 91 CHARACTER (len=64) :: cform = "(A12, 3(A13, I7))" 87 92 INTEGER , DIMENSION(jpi,jpj) :: ik_top , ik_bot ! top and bottom ocean level … … 121 126 WRITE(numout,*) ' cn_cfg = ', TRIM( cn_cfg ), ' nn_cfg = ', nn_cfg 122 127 ENDIF 123 nn_wxios = 0 124 ln_xios_read = .FALSE. 128 125 129 ! 126 130 ! !== Reference coordinate system ==! … … 143 147 hv_0(:,:) = 0._wp 144 148 hf_0(:,:) = 0._wp 145 DO jk = 1, jpk 149 DO jk = 1, jpkm1 146 150 ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk) 147 151 hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk) 148 152 hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk) 149 hf_0(:,:) = hf_0(:,:) + e3f_0(:,:,jk) * fmask(:,:,jk)150 153 END DO 154 ! 155 DO jk = 1, jpkm1 156 hf_0(1:jpim1,:) = hf_0(1:jpim1,:) + e3f_0(1:jpim1,:,jk)*vmask(1:jpim1,:,jk)*vmask(2:jpi,:,jk) 157 END DO 158 CALL lbc_lnk('domain', hf_0, 'F', 1._wp) 159 ! 160 IF( lk_SWE ) THEN ! SWE case redefine hf_0 161 hf_0(:,:) = hf_0(:,:) + e3f_0(:,:,1) * ssfmask(:,:) 162 ENDIF 151 163 ! 152 164 r1_ht_0(:,:) = ssmask (:,:) / ( ht_0(:,:) + 1._wp - ssmask (:,:) ) … … 154 166 r1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) ) 155 167 r1_hf_0(:,:) = ssfmask(:,:) / ( hf_0(:,:) + 1._wp - ssfmask(:,:) ) 156 168 ! 169 IF( ll_wd ) THEN ! wet and drying (check ht_0 >= 0) 170 DO_2D( 1, 1, 1, 1 ) 171 IF( ht_0(ji,jj) < 0._wp .AND. ssmask(ji,jj) == 1._wp ) THEN 172 CALL ctl_stop( 'ssh_init_rst : ht_0 must be positive at potentially wet points' ) 173 ENDIF 174 END_2D 175 ENDIF 176 ! 177 ! !== initialisation of time varying coordinate ==! 178 ! 179 ! != ssh initialization 180 IF( .NOT.l_offline .AND. .NOT.l_SAS ) THEN 181 CALL ssh_init_rst( Kbb, Kmm, Kaa ) 182 ELSE 183 ssh(:,:,:) = 0._wp 184 ENDIF 157 185 ! 158 186 #if defined key_qco 159 ! !== initialisation of time varying coordinate ==!Quasi-Euerian coordinate case187 ! != Quasi-Euerian coordinate case 160 188 ! 161 189 IF( .NOT.l_offline ) CALL dom_qco_init( Kbb, Kmm, Kaa ) 162 !163 IF( ln_linssh ) CALL ctl_stop('STOP','domain: key_qco and ln_linssh = T are incompatible')164 !165 190 #else 166 ! !== time varying part of coordinate system ==! 167 ! 168 IF( ln_linssh ) THEN != Fix in time : set to the reference one for all 191 ! 192 IF( ln_linssh ) THEN != Fix in time : set to the reference one for all 169 193 ! 170 194 DO jt = 1, jpt ! depth of t- and w-grid-points … … 175 199 ! 176 200 DO jt = 1, jpt ! vertical scale factors 177 e3t (:,:,:,jt) = e3t_0(:,:,:)178 e3u (:,:,:,jt) = e3u_0(:,:,:)179 e3v (:,:,:,jt) = e3v_0(:,:,:)180 e3w (:,:,:,jt) = e3w_0(:,:,:)201 e3t (:,:,:,jt) = e3t_0(:,:,:) 202 e3u (:,:,:,jt) = e3u_0(:,:,:) 203 e3v (:,:,:,jt) = e3v_0(:,:,:) 204 e3w (:,:,:,jt) = e3w_0(:,:,:) 181 205 e3uw(:,:,:,jt) = e3uw_0(:,:,:) 182 206 e3vw(:,:,:,jt) = e3vw_0(:,:,:) 183 207 END DO 184 e3f (:,:,:) = e3f_0(:,:,:)208 e3f (:,:,:) = e3f_0(:,:,:) 185 209 ! 186 210 DO jt = 1, jpt ! water column thickness and its inverse 187 hu(:,:,jt)= hu_0(:,:)188 hv(:,:,jt)= hv_0(:,:)211 hu(:,:,jt) = hu_0(:,:) 212 hv(:,:,jt) = hv_0(:,:) 189 213 r1_hu(:,:,jt) = r1_hu_0(:,:) 190 214 r1_hv(:,:,jt) = r1_hv_0(:,:) 191 215 END DO 192 ht(:,:) = ht_0(:,:)193 ! 194 ELSE != time varying : initialize before/now/after variables216 ht (:,:) = ht_0(:,:) 217 ! 218 ELSE != Time varying : initialize before/now/after variables 195 219 ! 196 220 IF( .NOT.l_offline ) CALL dom_vvl_init( Kbb, Kmm, Kaa ) … … 373 397 USE ioipsl 374 398 !! 375 INTEGER :: ios ! Local integer 399 INTEGER :: ios ! Local integer 400 REAL(wp):: zrdt 401 !!---------------------------------------------------------------------- 376 402 ! 377 403 NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list, & … … 393 419 ENDIF 394 420 ! 421 ! !=======================! 422 ! !== namelist namdom ==! 423 ! !=======================! 424 ! 425 READ ( numnam_ref, namdom, IOSTAT = ios, ERR = 903) 426 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist' ) 427 READ ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 ) 428 904 IF( ios > 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist' ) 429 IF(lwm) WRITE( numond, namdom ) 430 ! 431 #if defined key_agrif 432 IF( .NOT. Agrif_Root() ) THEN ! AGRIF child, subdivide the Parent timestep 433 rn_Dt = Agrif_Parent (rn_Dt ) / Agrif_Rhot() 434 ENDIF 435 #endif 436 ! 437 IF(lwp) THEN 438 WRITE(numout,*) 439 WRITE(numout,*) ' Namelist : namdom --- space & time domain' 440 WRITE(numout,*) ' linear free surface (=T) ln_linssh = ', ln_linssh 441 WRITE(numout,*) ' create mesh/mask file ln_meshmask = ', ln_meshmask 442 WRITE(numout,*) ' ocean time step rn_Dt = ', rn_Dt 443 WRITE(numout,*) ' asselin time filter parameter rn_atfp = ', rn_atfp 444 WRITE(numout,*) ' online coarsening of dynamical fields ln_crs = ', ln_crs 445 ENDIF 446 ! 447 ! set current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3 448 rDt = 2._wp * rn_Dt 449 r1_Dt = 1._wp / rDt 450 ! 451 IF( l_SAS .AND. .NOT.ln_linssh ) THEN 452 CALL ctl_warn( 'SAS requires linear ssh : force ln_linssh = T' ) 453 ln_linssh = .TRUE. 454 ENDIF 455 ! 456 #if defined key_qco 457 IF( ln_linssh ) CALL ctl_stop( 'STOP','domain: key_qco and ln_linssh = T are incompatible' ) 458 #endif 459 ! 460 ! !=======================! 461 ! !== namelist namrun ==! 462 ! !=======================! 395 463 ! 396 464 READ ( numnam_ref, namrun, IOSTAT = ios, ERR = 901) … … 452 520 nleapy = nn_leapy 453 521 ninist = nn_istate 522 ! 523 ! !== Set parameters for restart reading using xIOS ==! 524 ! 525 IF( TRIM(Agrif_CFixed()) == '0' ) THEN 526 lrxios = ln_xios_read .AND. ln_rstart 527 IF( nn_wxios > 0 ) lwxios = .TRUE. !* set output file type for XIOS based on NEMO namelist 528 nxioso = nn_wxios 529 ENDIF 530 ! !== Check consistency between ln_rstart and ln_1st_euler ==! (i.e. set l_1st_euler) 454 531 l_1st_euler = ln_1st_euler 455 IF( .NOT. l_1st_euler .AND. .NOT. ln_rstart ) THEN 532 ! 533 IF( ln_rstart ) THEN !* Restart case 534 ! 535 IF(lwp) WRITE(numout,*) 536 IF(lwp) WRITE(numout,*) ' open the restart file' 537 CALL rst_read_open !- Open the restart file 538 ! 539 IF( iom_varid( numror, 'rdt', ldstop = .FALSE. ) > 0 ) THEN !- Check time-step consistency and force Euler restart if changed 540 CALL iom_get( numror, 'rdt', zrdt ) 541 IF( zrdt /= rn_Dt ) THEN 542 IF(lwp) WRITE( numout,*) 543 IF(lwp) WRITE( numout,*) ' rn_Dt = ', rn_Dt,' not equal to the READ one rdt = ', zrdt 544 IF(lwp) WRITE( numout,*) 545 IF(lwp) WRITE( numout,*) ' ==>>> forced euler first time-step' 546 l_1st_euler = .TRUE. 547 ENDIF 548 ENDIF 549 ! 550 IF( .NOT.l_SAS .AND. iom_varid( numror, 'sshb', ldstop = .FALSE. ) <= 0 ) THEN !- Check absence of one of the Kbb field (here sshb) 551 ! ! (any Kbb field is missing ==> all Kbb fields are missing) 552 IF( .NOT.l_1st_euler ) THEN 553 CALL ctl_warn('dom_nam : ssh at Kbb not found in restart files ', & 554 & 'l_1st_euler forced to .true. and ' , & 555 & 'ssh(Kbb) = ssh(Kmm) ' ) 556 l_1st_euler = .TRUE. 557 ENDIF 558 ENDIF 559 ELSEIF( .NOT.l_1st_euler ) THEN !* Initialization case 456 560 IF(lwp) WRITE(numout,*) 457 561 IF(lwp) WRITE(numout,*)' ==>>> Start from rest (ln_rstart=F)' 458 562 IF(lwp) WRITE(numout,*)' an Euler initial time step is used : l_1st_euler is forced to .true. ' 459 l_1st_euler = .true. 460 ENDIF 461 ! ! control of output frequency 462 IF( .NOT. ln_rst_list ) THEN ! we use nn_stock 563 l_1st_euler = .TRUE. 564 ENDIF 565 ! 566 ! !== control of output frequency ==! 567 ! 568 IF( .NOT. ln_rst_list ) THEN ! we use nn_stock 463 569 IF( nn_stock == -1 ) CALL ctl_warn( 'nn_stock = -1 --> no restart will be done' ) 464 570 IF( nn_stock == 0 .OR. nn_stock > nitend ) THEN … … 479 585 IF( Agrif_Root() ) THEN 480 586 IF(lwp) WRITE(numout,*) 481 SELECT CASE ( nleapy ) ! Choose calendar for IOIPSL587 SELECT CASE ( nleapy ) !== Choose calendar for IOIPSL ==! 482 588 CASE ( 1 ) 483 589 CALL ioconf_calendar('gregorian') … … 491 597 END SELECT 492 598 ENDIF 493 494 READ ( numnam_ref, namdom, IOSTAT = ios, ERR = 903) 495 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist' ) 496 READ ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 ) 497 904 IF( ios > 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist' ) 498 IF(lwm) WRITE( numond, namdom ) 499 ! 500 #if defined key_agrif 501 IF( .NOT. Agrif_Root() ) THEN 502 rn_Dt = Agrif_Parent(rn_Dt) / Agrif_Rhot() 503 ENDIF 504 #endif 505 ! 506 IF(lwp) THEN 507 WRITE(numout,*) 508 WRITE(numout,*) ' Namelist : namdom --- space & time domain' 509 WRITE(numout,*) ' linear free surface (=T) ln_linssh = ', ln_linssh 510 WRITE(numout,*) ' create mesh/mask file ln_meshmask = ', ln_meshmask 511 WRITE(numout,*) ' ocean time step rn_Dt = ', rn_Dt 512 WRITE(numout,*) ' asselin time filter parameter rn_atfp = ', rn_atfp 513 WRITE(numout,*) ' online coarsening of dynamical fields ln_crs = ', ln_crs 514 ENDIF 515 ! 516 !! Initialise current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3 517 rDt = 2._wp * rn_Dt 518 r1_Dt = 1._wp / rDt 519 599 ! 600 ! !========================! 601 ! !== namelist namtile ==! 602 ! !========================! 603 ! 520 604 READ ( numnam_ref, namtile, IOSTAT = ios, ERR = 905 ) 521 605 905 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtile in reference namelist' ) … … 537 621 ENDIF 538 622 ENDIF 539 540 IF( TRIM(Agrif_CFixed()) == '0' ) THEN 541 lrxios = ln_xios_read.AND.ln_rstart 542 !set output file type for XIOS based on NEMO namelist 543 IF (nn_wxios > 0) lwxios = .TRUE. 544 nxioso = nn_wxios 545 ENDIF 546 623 ! 547 624 #if defined key_netcdf4 548 ! ! NetCDF 4 case ("key_netcdf4" defined) 625 ! !=======================! 626 ! !== namelist namnc4 ==! NetCDF 4 case ("key_netcdf4" defined) 627 ! !=======================! 628 ! 549 629 READ ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907) 550 630 907 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist' ) … … 555 635 IF(lwp) THEN ! control print 556 636 WRITE(numout,*) 557 WRITE(numout,*) ' Namelist namnc4 - Netcdf4 chunking parameters '637 WRITE(numout,*) ' Namelist namnc4 - Netcdf4 chunking parameters ("key_netcdf4" defined)' 558 638 WRITE(numout,*) ' number of chunks in i-dimension nn_nchunks_i = ', nn_nchunks_i 559 639 WRITE(numout,*) ' number of chunks in j-dimension nn_nchunks_j = ', nn_nchunks_j … … 618 698 SUBROUTINE domain_cfg( cd_cfg, kk_cfg, kpi, kpj, kpk, kperio ) 619 699 !!---------------------------------------------------------------------- 620 !! *** ROUTINE dom _nam***700 !! *** ROUTINE domain_cfg *** 621 701 !! 622 702 !! ** Purpose : read the domain size in domain configuration file -
NEMO/trunk/src/OCE/DOM/dommsk.F90
r13461 r14053 181 181 ssvmask(:,:) = MAXVAL( vmask(:,:,:), DIM=3 ) 182 182 ssfmask(:,:) = MAXVAL( fmask(:,:,:), DIM=3 ) 183 IF( lk_SWE ) THEN ! Shallow Water Eq. case : redefine ssfmask 184 DO_2D( 0,0 , 0,0 ) 185 ssfmask(ji,jj) = MAX( ssmask(ji,jj+1), ssmask(ji+1,jj+1), & 186 & ssmask(ji,jj ), ssmask(ji+1,jj ) ) 187 END_2D 188 CALL lbc_lnk( 'dommsk', ssfmask, 'F', 1.0_wp ) 189 ENDIF 190 #if defined key_qco 191 fe3mask(:,:,:) = fmask(:,:,:) 192 #endif 183 193 184 194 ! Interior domain mask (used for global sum) -
NEMO/trunk/src/OCE/DOM/domqco.F90
r13970 r14053 8 8 !! 3.3 ! 2011-10 (M. Leclair) totally rewrote domvvl: vvl option includes z_star and z_tilde coordinates 9 9 !! 3.6 ! 2014-11 (P. Mathiot) add ice shelf capability 10 !! 4.1 ! 2019-08 (A. Coward, D. Storkey) rename dom_vvl_sf_swp -> dom_vvl_sf_update for new timestepping11 !! 4.x ! 2020-02 (G. Madec, S. Techene) pure z* (quasi-eulerian) coordinate12 !!---------------------------------------------------------------------- 13 14 !!---------------------------------------------------------------------- 15 !! dom_q e_init: define initial vertical scale factors, depths and column thickness16 !! dom_q e_r3c : Compute ssh/h_0 ratioat t-, u-, v-, and optionally f-points17 !! qe_rst_read : read/write restart file18 !! dom_qe_ctl: Check the vvl options10 !! 4.1 ! 2019-08 (A. Coward, D. Storkey) add time level indices for prognostic variables 11 !! - ! 2020-02 (S. Techene, G. Madec) quasi-eulerian coordinate (z* or s*) 12 !!---------------------------------------------------------------------- 13 14 !!---------------------------------------------------------------------- 15 !! dom_qco_init : define initial vertical scale factors, depths and column thickness 16 !! dom_qco_zgr : Set ssh/h_0 ratio at t 17 !! dom_qco_r3c : Compute ssh/h_0 ratio at t-, u-, v-, and optionally f-points 18 !! qco_ctl : Check the vvl options 19 19 !!---------------------------------------------------------------------- 20 20 USE oce ! ocean dynamics and tracers … … 55 55 LOGICAL , PUBLIC :: ln_vvl_dbg = .FALSE. ! debug control prints 56 56 57 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_td, vn_td ! thickness diffusion transport58 59 57 !! * Substitutions 60 58 # include "do_loop_substitute.h90" … … 79 77 !! 80 78 !!---------------------------------------------------------------------- 81 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa 79 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! time level indices 80 !!---------------------------------------------------------------------- 82 81 ! 83 82 IF(lwp) WRITE(numout,*) … … 85 84 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 86 85 ! 87 CALL dom_qco_ctl ! choose vertical coordinate (z_star, z_tilde or layer) 88 ! 89 ! ! Read or initialize e3t_(b/n), tilde_e3t_(b/n) and hdiv_lf 90 CALL qe_rst_read( nit000, Kbb, Kmm ) 91 ! 92 CALL dom_qco_zgr(Kbb, Kmm, Kaa) ! interpolation scale factor, depth and water column 86 CALL qco_ctl ! choose vertical coordinate (z_star, z_tilde or layer) 87 ! 88 CALL dom_qco_zgr( Kbb, Kmm ) ! interpolation scale factor, depth and water column 89 ! 90 #if defined key_agrif 91 ! We need to define r3[tuv](Kaa) for AGRIF initialisation (should not be a 92 ! problem for the restartability...) 93 r3t(:,:,Kaa) = r3t(:,:,Kmm) 94 r3u(:,:,Kaa) = r3u(:,:,Kmm) 95 r3v(:,:,Kaa) = r3v(:,:,Kmm) 96 #endif 93 97 ! 94 98 END SUBROUTINE dom_qco_init 95 99 96 100 97 SUBROUTINE dom_qco_zgr( Kbb, Kmm, Kaa)101 SUBROUTINE dom_qco_zgr( Kbb, Kmm ) 98 102 !!---------------------------------------------------------------------- 99 103 !! *** ROUTINE dom_qco_init *** 100 104 !! 101 !! ** Purpose : Initialization of all ssh. to h._0 ratio 102 !! 103 !! ** Method : - interpolate scale factors 104 !! 105 !! ** Action : - r3(t/u/v)_b 106 !! - r3(t/u/v/f)_n 107 !! 108 !! Reference : Leclair, M., and G. Madec, 2011, Ocean Modelling. 109 !!---------------------------------------------------------------------- 110 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa 105 !! ** Purpose : Initialization of all r3. = ssh./h._0 ratios 106 !! 107 !! ** Method : Call domqco using Kbb and Kmm 108 !! NB: dom_qco_zgr is called by dom_qco_init it uses ssh from ssh_init 109 !! 110 !! ** Action : - r3(t/u/v)(Kbb) 111 !! - r3(t/u/v/f)(Kmm) 112 !!---------------------------------------------------------------------- 113 INTEGER, INTENT(in) :: Kbb, Kmm ! time level indices 111 114 !!---------------------------------------------------------------------- 112 115 ! 113 116 ! !== Set of all other vertical scale factors ==! (now and before) 114 117 ! ! Horizontal interpolation of e3t 115 CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) )118 CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) ) 116 119 CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm), r3f(:,:) ) 117 120 ! … … 143 146 ! !== ratio at u-,v-point ==! 144 147 ! 145 IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging) 148 !!st IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging) 149 #if ! defined key_qcoTest_FluxForm 150 ! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 146 151 DO_2D( 0, 0, 0, 0 ) 147 152 pr3u(ji,jj) = 0.5_wp * ( e1e2t(ji ,jj) * pssh(ji ,jj) & … … 150 155 & + e1e2t(ji,jj+1) * pssh(ji,jj+1) ) * r1_hv_0(ji,jj) * r1_e1e2v(ji,jj) 151 156 END_2D 152 ELSE !- Flux Form (simple averaging) 157 !!st ELSE !- Flux Form (simple averaging) 158 #else 153 159 DO_2D( 0, 0, 0, 0 ) 154 pr3u(ji,jj) = 0.5_wp * ( pssh(ji ,jj) + pssh(ji+1,jj) ) * r1_hu_0(ji,jj)155 pr3v(ji,jj) = 0.5_wp * ( pssh(ji,jj ) + pssh(ji,jj+1) ) * r1_hv_0(ji,jj)160 pr3u(ji,jj) = 0.5_wp * ( pssh(ji,jj) + pssh(ji+1,jj ) ) * r1_hu_0(ji,jj) 161 pr3v(ji,jj) = 0.5_wp * ( pssh(ji,jj) + pssh(ji ,jj+1) ) * r1_hv_0(ji,jj) 156 162 END_2D 157 ENDIF 163 !!st ENDIF 164 #endif 158 165 ! 159 166 IF( .NOT.PRESENT( pr3f ) ) THEN !- lbc on ratio at u-, v-points only … … 163 170 ELSE !== ratio at f-point ==! 164 171 ! 165 IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging) 166 DO_2D( 1, 0, 1, 0 ) ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 172 !!st IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging) 173 #if ! defined key_qcoTest_FluxForm 174 ! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 175 176 DO_2D( 0, 0, 0, 0 ) ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 167 177 pr3f(ji,jj) = 0.25_wp * ( e1e2t(ji ,jj ) * pssh(ji ,jj ) & 168 178 & + e1e2t(ji+1,jj ) * pssh(ji+1,jj ) & … … 170 180 & + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1) ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) 171 181 END_2D 172 ELSE !- Flux Form (simple averaging) 173 DO_2D( 1, 0, 1, 0 ) ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 174 pr3f(ji,jj) = 0.25_wp * ( pssh(ji ,jj ) + pssh(ji+1,jj ) & 175 & + pssh(ji ,jj+1) + pssh(ji+1,jj+1) ) * r1_hf_0(ji,jj) 182 !!st ELSE !- Flux Form (simple averaging) 183 #else 184 DO_2D( 0, 0, 0, 0 ) ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 185 pr3f(ji,jj) = 0.25_wp * ( pssh(ji,jj ) + pssh(ji+1,jj ) & 186 & + pssh(ji,jj+1) + pssh(ji+1,jj+1) ) * r1_hf_0(ji,jj) 176 187 END_2D 177 ENDIF 188 !!st ENDIF 189 #endif 178 190 ! ! lbc on ratio at u-,v-,f-points 179 191 CALL lbc_lnk_multi( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp ) … … 184 196 185 197 186 SUBROUTINE q e_rst_read( kt, Kbb, Kmm )198 SUBROUTINE qco_ctl 187 199 !!--------------------------------------------------------------------- 188 !! *** ROUTINE qe_rst_read *** 189 !! 190 !! ** Purpose : Read ssh in restart file 191 !! 192 !! ** Method : use of IOM library 193 !! if the restart does not contain ssh, 194 !! it is set to the _0 values. 195 !!---------------------------------------------------------------------- 196 INTEGER , INTENT(in) :: kt ! ocean time-step 197 INTEGER , INTENT(in) :: Kbb, Kmm ! ocean time level indices 198 ! 199 INTEGER :: ji, jj, jk 200 INTEGER :: id1, id2 ! local integers 201 !!---------------------------------------------------------------------- 202 ! 203 IF( ln_rstart ) THEN !* Read the restart file 204 CALL rst_read_open ! open the restart file if necessary 205 ! 206 id1 = iom_varid( numror, 'sshb', ldstop = .FALSE. ) 207 id2 = iom_varid( numror, 'sshn', ldstop = .FALSE. ) 208 ! 209 ! ! --------- ! 210 ! ! all cases ! 211 ! ! --------- ! 212 ! 213 IF( MIN( id1, id2 ) > 0 ) THEN ! all required arrays exist 214 CALL iom_get( numror, jpdom_auto, 'sshb' , ssh(:,:,Kbb) ) 215 CALL iom_get( numror, jpdom_auto, 'sshn' , ssh(:,:,Kmm) ) 216 ! needed to restart if land processor not computed 217 IF(lwp) write(numout,*) 'qe_rst_read : ssh(:,:,Kbb) and ssh(:,:,Kmm) found in restart files' 218 WHERE ( ssmask(:,:) == 0.0_wp ) !!gm/st ==> sm should not be necessary on ssh when it was required on e3 219 ssh(:,:,Kmm) = 0._wp 220 ssh(:,:,Kbb) = 0._wp 221 END WHERE 222 IF( l_1st_euler ) THEN 223 ssh(:,:,Kbb) = ssh(:,:,Kmm) 224 ENDIF 225 ELSE IF( id1 > 0 ) THEN 226 IF(lwp) write(numout,*) 'qe_rst_read WARNING : ssh(:,:,Kmm) not found in restart files' 227 IF(lwp) write(numout,*) 'sshn set equal to sshb.' 228 IF(lwp) write(numout,*) 'neuler is forced to 0' 229 CALL iom_get( numror, jpdom_auto, 'sshb', ssh(:,:,Kbb) ) 230 ssh(:,:,Kmm) = ssh(:,:,Kbb) 231 l_1st_euler = .TRUE. 232 ELSE IF( id2 > 0 ) THEN 233 IF(lwp) write(numout,*) 'qe_rst_read WARNING : ssh(:,:,Kbb) not found in restart files' 234 IF(lwp) write(numout,*) 'sshb set equal to sshn.' 235 IF(lwp) write(numout,*) 'neuler is forced to 0' 236 CALL iom_get( numror, jpdom_auto, 'sshn', ssh(:,:,Kmm) ) 237 ssh(:,:,Kbb) = ssh(:,:,Kmm) 238 l_1st_euler = .TRUE. 239 ELSE 240 IF(lwp) write(numout,*) 'qe_rst_read WARNING : ssh(:,:,Kmm) not found in restart file' 241 IF(lwp) write(numout,*) 'ssh_b and ssh_n set to zero' 242 IF(lwp) write(numout,*) 'neuler is forced to 0' 243 ssh(:,:,:) = 0._wp 244 l_1st_euler = .TRUE. 245 ENDIF 246 ! 247 ELSE !* Initialize at "rest" 248 ! 249 IF( ll_wd ) THEN ! MJB ll_wd edits start here - these are essential 250 ! 251 IF( cn_cfg == 'wad' ) THEN ! Wetting and drying test case 252 CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb) ) 253 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones 254 ssh(:,: ,Kmm) = ssh(:,: ,Kbb) 255 uu (:,:,: ,Kmm) = uu (:,:,: ,Kbb) 256 vv (:,:,: ,Kmm) = vv (:,:,: ,Kbb) 257 ELSE ! if not test case 258 ssh(:,:,Kmm) = -ssh_ref 259 ssh(:,:,Kbb) = -ssh_ref 260 ! 261 DO_2D( 1, 1, 1, 1 ) 262 IF( ht_0(ji,jj)-ssh_ref < rn_wdmin1 ) THEN ! if total depth is less than min depth 263 ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 264 ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) ) 265 ENDIF 266 END_2D 267 ENDIF 268 269 DO ji = 1, jpi 270 DO jj = 1, jpj 271 IF ( ht_0(ji,jj) .LE. 0.0 .AND. NINT( ssmask(ji,jj) ) .EQ. 1) THEN 272 CALL ctl_stop( 'qe_rst_read: ht_0 must be positive at potentially wet points' ) 273 ENDIF 274 END DO 275 END DO 276 ! 277 ELSE 278 ! 279 ! Just to read set ssh in fact, called latter once vertical grid 280 ! is set up: 281 ! CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb) ) 282 ! ! 283 ssh(:,:,:) = 0._wp 284 ! 285 ENDIF ! end of ll_wd edits 286 ! 287 ENDIF 288 ! 289 END SUBROUTINE qe_rst_read 290 291 292 SUBROUTINE dom_qco_ctl 293 !!--------------------------------------------------------------------- 294 !! *** ROUTINE dom_qco_ctl *** 200 !! *** ROUTINE qco_ctl *** 295 201 !! 296 202 !! ** Purpose : Control the consistency between namelist options … … 312 218 IF(lwp) THEN ! Namelist print 313 219 WRITE(numout,*) 314 WRITE(numout,*) ' dom_qco_ctl : choice/control of the variable vertical coordinate'315 WRITE(numout,*) '~~~~~~~~ ~~~'220 WRITE(numout,*) 'qco_ctl : choice/control of the variable vertical coordinate' 221 WRITE(numout,*) '~~~~~~~~' 316 222 WRITE(numout,*) ' Namelist nam_vvl : chose a vertical coordinate' 317 223 WRITE(numout,*) ' zstar ln_vvl_zstar = ', ln_vvl_zstar … … 357 263 #endif 358 264 ! 359 END SUBROUTINE dom_qco_ctl265 END SUBROUTINE qco_ctl 360 266 361 267 !!====================================================================== -
NEMO/trunk/src/OCE/DOM/domvvl.F90
r13982 r14053 9 9 !! 3.6 ! 2014-11 (P. Mathiot) add ice shelf capability 10 10 !! 4.1 ! 2019-08 (A. Coward, D. Storkey) rename dom_vvl_sf_swp -> dom_vvl_sf_update for new timestepping 11 !! 4.x !2020-02 (G. Madec, S. Techene) introduce ssh to h0 ratio11 !! - ! 2020-02 (G. Madec, S. Techene) introduce ssh to h0 ratio 12 12 !!---------------------------------------------------------------------- 13 13 … … 766 766 !! ** Purpose : Read or write VVL file in restart file 767 767 !! 768 !! ** Method : use of IOM library 769 !! if the restart does not contain vertical scale factors, 770 !! they are set to the _0 values 771 !! if the restart does not contain vertical scale factors increments (z_tilde), 772 !! they are set to 0. 768 !! ** Method : * restart comes from a linear ssh simulation : 769 !! an attempt to read e3t_n stops simulation 770 !! * restart comes from a z-star, z-tilde, or layer : 771 !! read e3t_n and e3t_b 772 !! * restart comes from a z-star : 773 !! set tilde_e3t_n, tilde_e3t_n, and hdiv_lf to 0 774 !! * restart comes from layer : 775 !! read tilde_e3t_n and tilde_e3t_b 776 !! set hdiv_lf to 0 777 !! * restart comes from a z-tilde: 778 !! read tilde_e3t_n, tilde_e3t_b, and hdiv_lf 779 !! 780 !! NB: if l_1st_euler = T (ln_1st_euler or ssh_b not found) 781 !! Kbb fields set to Kmm ones 773 782 !!---------------------------------------------------------------------- 774 783 INTEGER , INTENT(in) :: kt ! ocean time-step … … 776 785 CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag 777 786 ! 778 INTEGER :: ji, jj, jk 779 INTEGER :: id 1, id2, id3, id4, id5! local integers780 !!---------------------------------------------------------------------- 781 ! 782 IF( TRIM(cdrw) == 'READ' ) THEN ! Read/initialise783 ! ! ===============784 IF( ln_rstart ) THEN !* Read the restart file785 CALL rst_read_open ! open the restart file if necessary786 CALL iom_get( numror, jpdom_auto, 'sshn' , ssh(:,:,Kmm) )787 INTEGER :: ji, jj, jk ! dummy loop indices 788 INTEGER :: id3, id4, id5 ! local integers 789 !!---------------------------------------------------------------------- 790 ! 791 ! !=====================! 792 IF( TRIM(cdrw) == 'READ' ) THEN ! Read / initialise ! 793 ! !=====================! 794 ! 795 IF( ln_rstart ) THEN !== Read the restart file ==! 787 796 ! 788 id1 = iom_varid( numror, 'e3t_b', ldstop = .FALSE. ) 789 id2 = iom_varid( numror, 'e3t_n', ldstop = .FALSE. ) 790 id3 = iom_varid( numror, 'tilde_e3t_b', ldstop = .FALSE. ) 797 CALL rst_read_open !* open the restart file if necessary 798 ! ! --------- ! 799 ! ! all cases ! 800 ! ! --------- ! 801 ! 802 id3 = iom_varid( numror, 'tilde_e3t_b', ldstop = .FALSE. ) !* check presence 791 803 id4 = iom_varid( numror, 'tilde_e3t_n', ldstop = .FALSE. ) 792 id5 = iom_varid( numror, 'hdiv_lf' , ldstop = .FALSE. )804 id5 = iom_varid( numror, 'hdiv_lf' , ldstop = .FALSE. ) 793 805 ! 794 ! ! --------- ! 795 ! ! all cases ! 796 ! ! --------- ! 797 ! 798 IF( MIN( id1, id2 ) > 0 ) THEN ! all required arrays exist 806 ! !* scale factors 807 IF(lwp) WRITE(numout,*) ' Kmm scale factor read in the restart file' 808 CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm) ) 809 WHERE ( tmask(:,:,:) == 0.0_wp ) 810 e3t(:,:,:,Kmm) = e3t_0(:,:,:) 811 END WHERE 812 IF( l_1st_euler ) THEN ! euler 813 IF(lwp) WRITE(numout,*) ' Euler first time step : e3t(Kbb) = e3t(Kmm)' 814 e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 815 ELSE ! leap frog 816 IF(lwp) WRITE(numout,*) ' Kbb scale factor read in the restart file' 799 817 CALL iom_get( numror, jpdom_auto, 'e3t_b', e3t(:,:,:,Kbb) ) 800 CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm) )801 ! needed to restart if land processor not computed802 IF(lwp) write(numout,*) 'dom_vvl_rst : e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) found in restart files'803 818 WHERE ( tmask(:,:,:) == 0.0_wp ) 804 e3t(:,:,:,Kmm) = e3t_0(:,:,:)805 819 e3t(:,:,:,Kbb) = e3t_0(:,:,:) 806 820 END WHERE 807 IF( l_1st_euler ) THEN808 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)809 ENDIF810 ELSE IF( id1 > 0 ) THEN811 IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart files'812 IF(lwp) write(numout,*) 'e3t_n set equal to e3t_b.'813 IF(lwp) write(numout,*) 'l_1st_euler is forced to true'814 CALL iom_get( numror, jpdom_auto, 'e3t_b', e3t(:,:,:,Kbb) )815 e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb)816 l_1st_euler = .true.817 ELSE IF( id2 > 0 ) THEN818 IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kbb) not found in restart files'819 IF(lwp) write(numout,*) 'e3t_b set equal to e3t_n.'820 IF(lwp) write(numout,*) 'l_1st_euler is forced to true'821 CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm) )822 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)823 l_1st_euler = .true.824 ELSE825 IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart file'826 IF(lwp) write(numout,*) 'Compute scale factor from sshn'827 IF(lwp) write(numout,*) 'l_1st_euler is forced to true'828 DO jk = 1, jpk829 e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) &830 & / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk) &831 & + e3t_0(:,:,jk) * (1._wp -tmask(:,:,jk))832 END DO833 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)834 l_1st_euler = .true.835 821 ENDIF 836 ! !----------- !837 IF( ln_vvl_zstar ) THEN !z_star case !838 ! !----------- !822 ! ! ------------ ! 823 IF( ln_vvl_zstar ) THEN ! z_star case ! 824 ! ! ------------ ! 839 825 IF( MIN( id3, id4 ) > 0 ) THEN 840 826 CALL ctl_stop( 'dom_vvl_rst: z_star cannot restart from a z_tilde or layer run' ) 841 827 ENDIF 842 ! ! ----------------------- ! 843 ELSE ! z_tilde and layer cases ! 844 ! ! ----------------------- ! 845 IF( MIN( id3, id4 ) > 0 ) THEN ! all required arrays exist 846 CALL iom_get( numror, jpdom_auto, 'tilde_e3t_b', tilde_e3t_b(:,:,:) ) 828 ! ! ------------------------ ! 829 ELSE ! z_tilde and layer cases ! 830 ! ! ------------------------ ! 831 ! 832 IF( id4 > 0 ) THEN !* scale factor increments 833 IF(lwp) WRITE(numout,*) ' Kmm scale factor increments read in the restart file' 847 834 CALL iom_get( numror, jpdom_auto, 'tilde_e3t_n', tilde_e3t_n(:,:,:) ) 848 ELSE ! one at least array is missing 835 IF( l_1st_euler ) THEN ! euler 836 IF(lwp) WRITE(numout,*) ' Euler first time step : tilde_e3t(Kbb) = tilde_e3t(Kmm)' 837 tilde_e3t_b(:,:,:) = tilde_e3t_n(:,:,:) 838 ELSE ! leap frog 839 IF(lwp) WRITE(numout,*) ' Kbb scale factor increments read in the restart file' 840 CALL iom_get( numror, jpdom_auto, 'tilde_e3t_b', tilde_e3t_b(:,:,:) ) 841 ENDIF 842 ELSE 849 843 tilde_e3t_b(:,:,:) = 0.0_wp 850 844 tilde_e3t_n(:,:,:) = 0.0_wp 851 845 ENDIF 852 ! ! ------------ !853 IF( ln_vvl_ztilde ) THEN ! z_tilde case !854 ! ! ------------ !846 ! ! ------------ ! 847 IF( ln_vvl_ztilde ) THEN ! z_tilde case ! 848 ! ! ------------ ! 855 849 IF( id5 > 0 ) THEN ! required array exists 856 850 CALL iom_get( numror, jpdom_auto, 'hdiv_lf', hdiv_lf(:,:,:) ) 857 851 ELSE ! array is missing 858 hdiv_lf(:,:,:) = 0.0_wp 852 hdiv_lf(:,:,:) = 0.0_wp 859 853 ENDIF 860 854 ENDIF 861 855 ENDIF 862 856 ! 863 ELSE ! * Initialize at "rest"857 ELSE !== Initialize at "rest" with ssh ==! 864 858 ! 865 866 IF( ll_wd ) THEN ! MJB ll_wd edits start here - these are essential 867 ! 868 IF( cn_cfg == 'wad' ) THEN 869 ! Wetting and drying test case 870 CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb) ) 871 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones 872 ssh (:,:,Kmm) = ssh(:,:,Kbb) 873 uu (:,:,:,Kmm) = uu (:,:,:,Kbb) 874 vv (:,:,:,Kmm) = vv (:,:,:,Kbb) 875 ELSE 876 ! if not test case 877 ssh(:,:,Kmm) = -ssh_ref 878 ssh(:,:,Kbb) = -ssh_ref 879 880 DO_2D( 1, 1, 1, 1 ) 881 IF( ht_0(ji,jj)-ssh_ref < rn_wdmin1 ) THEN ! if total depth is less than min depth 882 ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 883 ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) ) 884 ENDIF 885 END_2D 886 ENDIF !If test case else 887 888 ! Adjust vertical metrics for all wad 889 DO jk = 1, jpk 890 e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) & 891 & / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk) & 892 & + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) ) 893 END DO 894 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 895 896 DO_2D( 1, 1, 1, 1 ) 897 IF ( ht_0(ji,jj) .LE. 0.0 .AND. NINT( ssmask(ji,jj) ) .EQ. 1) THEN 898 CALL ctl_stop( 'dom_vvl_rst: ht_0 must be positive at potentially wet points' ) 899 ENDIF 900 END_2D 901 ! 902 ELSE 903 ! 904 ! Just to read set ssh in fact, called latter once vertical grid 905 ! is set up: 906 ! CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb) ) 907 ! ! 908 ! DO jk=1,jpk 909 ! e3t(:,:,jk,Kbb) = e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb) ) & 910 ! & / ( ht_0(:,:) + 1._wp -ssmask(:,:) ) * tmask(:,:,jk) 911 ! END DO 912 ! e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 913 ssh(:,:,Kmm)=0._wp 914 e3t(:,:,:,Kmm)=e3t_0(:,:,:) 915 e3t(:,:,:,Kbb)=e3t_0(:,:,:) 916 ! 917 END IF ! end of ll_wd edits 918 859 DO jk = 1, jpk 860 e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( 1._wp + ssh(:,:,Kmm) * r1_ht_0(:,:) * tmask(:,:,jk) ) 861 END DO 862 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 863 ! 919 864 IF( ln_vvl_ztilde .OR. ln_vvl_layer) THEN 920 865 tilde_e3t_b(:,:,:) = 0._wp 921 866 tilde_e3t_n(:,:,:) = 0._wp 922 867 IF( ln_vvl_ztilde ) hdiv_lf(:,:,:) = 0._wp 923 END 868 ENDIF 924 869 ENDIF 925 ! 926 ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN ! Create restart file 927 ! ! =================== 870 ! !=======================! 871 ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN ! Create restart file ! 872 ! !=======================! 873 ! 928 874 IF(lwp) WRITE(numout,*) '---- dom_vvl_rst ----' 929 875 ! ! --------- ! -
NEMO/trunk/src/OCE/DOM/domzgr_substitute.h90
r13237 r14053 15 15 # define e3u(i,j,k,t) (e3u_0(i,j,k)*(1._wp+r3u(i,j,t)*umask(i,j,k))) 16 16 # define e3v(i,j,k,t) (e3v_0(i,j,k)*(1._wp+r3v(i,j,t)*vmask(i,j,k))) 17 # define e3f(i,j,k) (e3f_0(i,j,k)*(1._wp+r3f(i,j)*fmask(i,j,k))) 17 # define e3f(i,j,k) (e3f_0(i,j,k)*(1._wp+r3f(i,j)*fe3mask(i,j,k))) 18 # define e3f_vor(i,j,k) (e3f_0vor(i,j,k)*(1._wp+r3f(i,j)*fe3mask(i,j,k))) 18 19 # define e3w(i,j,k,t) (e3w_0(i,j,k)*(1._wp+r3t(i,j,t))) 19 20 # define e3uw(i,j,k,t) (e3uw_0(i,j,k)*(1._wp+r3u(i,j,t))) 20 21 # define e3vw(i,j,k,t) (e3vw_0(i,j,k)*(1._wp+r3v(i,j,t))) 21 # define ht(i,j) (ht_0(i,j) +ssh(i,j,Kmm))22 # define ht(i,j) (ht_0(i,j)*(1._wp+r3t(i,j,Kmm))) 22 23 # define hu(i,j,t) (hu_0(i,j)*(1._wp+r3u(i,j,t))) 23 24 # define hv(i,j,t) (hv_0(i,j)*(1._wp+r3v(i,j,t))) … … 29 30 #endif 30 31 !!---------------------------------------------------------------------- 32 !!# define e3t_f(i,j,k) (e3t_0(i,j,k)*(1._wp+r3t_f(i,j)*tmask(i,j,k))) 33 !!# define e3u_f(i,j,k) (e3u_0(i,j,k)*(1._wp+r3u_f(i,j)*umask(i,j,k))) 34 !!# define e3v_f(i,j,k) (e3v_0(i,j,k)*(1._wp+r3v_f(i,j)*vmask(i,j,k))) -
NEMO/trunk/src/OCE/DOM/istate.F90
r13295 r14053 42 42 PRIVATE 43 43 44 PUBLIC istate_init ! routine called by step.F9044 PUBLIC istate_init ! routine called by nemogcm.F90 45 45 46 46 !! * Substitutions … … 59 59 !! 60 60 !! ** Purpose : Initialization of the dynamics and tracer fields. 61 !! 62 !! ** Method : 61 63 !!---------------------------------------------------------------------- 62 64 INTEGER, INTENT( in ) :: Kbb, Kmm, Kaa ! ocean time level indices 63 65 ! 64 66 INTEGER :: ji, jj, jk ! dummy loop indices 65 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgdept ! 3D table !!st patch to use gdept subtitute67 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgdept ! 3D table for qco substitute 66 68 !!gm see comment further down 67 69 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: zuvd ! U & V data workspace … … 73 75 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 74 76 75 !!gm Why not include in the first call of dta_tsd ?76 !!gm probably associated with the use of internal damping...77 77 CALL dta_tsd_init ! Initialisation of T & S input data 78 !!gm to be moved in usrdef of C1D case 78 79 79 ! IF( lk_c1d ) CALL dta_uvd_init ! Initialization of U & V input data 80 !!gm81 80 82 rhd (:,:,: ) = 0._wp ; rhop (:,:,: ) = 0._wp ! set one for all to 0 at level jpk83 rn2b (:,:,: ) = 0._wp ; rn2 (:,:,: ) = 0._wp ! set one for all to 0 at levels 1 and jpk84 ts (:,:,:,:,Kaa) = 0._wp ! set one for all to 0 at level jpk85 rab_b(:,:,:,: ) = 0._wp ; rab_n(:,:,:,:) = 0._wp ! set one for all to 0 at level jpk81 rhd (:,:,: ) = 0._wp ; rhop (:,:,: ) = 0._wp ! set one for all to 0 at level jpk 82 rn2b (:,:,: ) = 0._wp ; rn2 (:,:,: ) = 0._wp ! set one for all to 0 at levels 1 and jpk 83 ts (:,:,:,:,Kaa) = 0._wp ! set one for all to 0 at level jpk 84 rab_b(:,:,:,: ) = 0._wp ; rab_n(:,:,:,:) = 0._wp ! set one for all to 0 at level jpk 86 85 #if defined key_agrif 87 86 uu (:,:,: ,Kaa) = 0._wp ! used in agrif_oce_sponge at initialization … … 96 95 CALL agrif_istate( Kbb, Kmm, Kaa ) ! Interp from parent 97 96 ! 98 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) 99 ssh (:,:,Kmm) = ssh(:,:,Kbb) 100 uu (:,:,:,Kmm) = uu (:,:,:,Kbb) 101 vv (:,:,:,Kmm) = vv (:,:,:,Kbb) 97 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) 98 !!st 99 !!st need for a recent agrif version to be displaced toward ssh_init_rst with agrif_istate_ssh 100 ssh(:,:, Kmm) = ssh(:,: ,Kbb) 101 !!st end 102 uu (:,:,: ,Kmm) = uu (:,:,: ,Kbb) 103 vv (:,:,: ,Kmm) = vv (:,:,: ,Kbb) 102 104 ELSE 103 105 #endif … … 117 119 CALL dta_tsd( nit000, ts(:,:,:,:,Kbb) ) ! read 3D T and S data at nit000 118 120 ! 119 ssh(:,:,Kbb) = 0._wp ! set the ocean at rest 120 uu (:,:,:,Kbb) = 0._wp 121 vv (:,:,:,Kbb) = 0._wp 121 uu (:,:,:,Kbb) = 0._wp 122 vv (:,:,:,Kbb) = 0._wp 122 123 ! 123 IF( ll_wd ) THEN124 ssh(:,:,Kbb) = -ssh_ref ! Added in 30 here for bathy that adds 30 as Iterative test CEOD125 !126 ! Apply minimum wetdepth criterion127 !128 DO_2D( 1, 1, 1, 1 )129 IF( ht_0(ji,jj) + ssh(ji,jj,Kbb) < rn_wdmin1 ) THEN130 ssh(ji,jj,Kbb) = tmask(ji,jj,1)*( rn_wdmin1 - (ht_0(ji,jj)) )131 ENDIF132 END_2D133 ENDIF134 !135 124 ELSE ! user defined initial T and S 136 125 DO jk = 1, jpk 137 126 zgdept(:,:,jk) = gdept(:,:,jk,Kbb) 138 127 END DO 139 CALL usr_def_istate( zgdept, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb) , ssh(:,:,Kbb) )128 CALL usr_def_istate( zgdept, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb) ) 140 129 ENDIF 141 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones 142 ssh (:,:,Kmm) = ssh(:,:,Kbb) 143 uu (:,:,:,Kmm) = uu (:,:,:,Kbb) 144 vv (:,:,:,Kmm) = vv (:,:,:,Kbb) 145 146 !!gm POTENTIAL BUG : 147 !!gm ISSUE : if ssh(:,:,Kbb) /= 0 then, in non linear free surface, the e3._n, e3._b should be recomputed 148 !! as well as gdept_ and gdepw_.... !!!!! 149 !! ===>>>> probably a call to domvvl initialisation here.... 150 130 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones 131 uu (:,:,: ,Kmm) = uu (:,:,: ,Kbb) 132 vv (:,:,: ,Kmm) = vv (:,:,: ,Kbb) 151 133 152 134 ! 153 !!gm to be moved in usrdef of C1D case154 !IF ( ln_uvd_init .AND. lk_c1d ) THEN ! read 3D U and V data at nit000155 !ALLOCATE( zuvd(jpi,jpj,jpk,2) )156 ! CALL dta_uvd( nit000, zuvd )157 ! uu(:,:,:,Kbb) = zuvd(:,:,:,1); uu(:,:,:,Kmm) = uu(:,:,:,Kbb)158 ! vv(:,:,:,Kbb) = zuvd(:,:,:,2); vv(:,:,:,Kmm) = vv(:,:,:,Kbb)159 !DEALLOCATE( zuvd )160 !ENDIF135 !!gm ==>>> to be moved in usrdef_istate of C1D case 136 IF ( ln_uvd_init .AND. lk_c1d ) THEN ! read 3D U and V data at nit000 137 ALLOCATE( zuvd(jpi,jpj,jpk,2) ) 138 CALL dta_uvd( nit000, Kbb, zuvd ) 139 uu(:,:,:,Kbb) = zuvd(:,:,:,1) ; uu(:,:,:,Kmm) = uu(:,:,:,Kbb) 140 vv(:,:,:,Kbb) = zuvd(:,:,:,2) ; vv(:,:,:,Kmm) = vv(:,:,:,Kbb) 141 DEALLOCATE( zuvd ) 142 ENDIF 161 143 ! 162 !!gm This is to be changed !!!!163 ! ! - ML - ssh(:,:,Kmm) could be modified by istate_eel, so that initialization of e3t(:,:,:,Kbb) is done here164 ! IF( .NOT.ln_linssh ) THEN165 ! DO jk = 1, jpk166 ! e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm)167 ! END DO168 ! ENDIF169 !!gm170 144 ! 171 145 ENDIF -
NEMO/trunk/src/OCE/DOM/phycst.F90
r12489 r14053 66 66 REAL(wp), PUBLIC :: r1_rhos !: 1 / rhos 67 67 REAL(wp), PUBLIC :: r1_rcpi !: 1 / rcpi 68 68 69 !!---------------------------------------------------------------------- 69 70 !! NEMO/OCE 4.0 , NEMO Consortium (2018) -
NEMO/trunk/src/OCE/DYN/dynadv.F90
r12377 r14053 127 127 IF( ioptio /= 1 ) CALL ctl_stop( 'choose ONE and only ONE advection scheme' ) 128 128 IF( nn_dynkeg /= nkeg_C2 .AND. nn_dynkeg /= nkeg_HW ) CALL ctl_stop( 'KEG scheme wrong value of nn_dynkeg' ) 129 129 #if defined key_qcoTest_FluxForm 130 IF( ln_dynadv_vec ) THEN CALL ctl_stop( 'STOP', 'key_qcoTest_FluxForm requires flux form advection' ) 131 #endif 130 132 131 133 IF(lwp) THEN ! Print the choice -
NEMO/trunk/src/OCE/DYN/dynatf_qco.F90
r13295 r14053 1 MODULE dynatf qco1 MODULE dynatf_qco 2 2 !!========================================================================= 3 !! *** MODULE dynatf qco ***3 !! *** MODULE dynatf_qco *** 4 4 !! Ocean dynamics: time filtering 5 5 !!========================================================================= … … 50 50 USE prtctl ! Print control 51 51 USE timing ! Timing 52 #if defined key_agrif53 USE agrif_oce_interp54 #endif55 52 56 53 IMPLICIT NONE … … 199 196 ! JC: Would be more clever to swap variables than to make a full vertical 200 197 ! integration 201 ! 198 ! CAUTION : calculation need to be done in the same way than see GM 202 199 uu_b(:,:,Kaa) = e3u(:,:,1,Kaa) * puu(:,:,1,Kaa) * umask(:,:,1) 203 uu_b(:,:,Kmm) = e3u(:,:,1,Kmm) * puu(:,:,1,Kmm) * umask(:,:,1)200 uu_b(:,:,Kmm) = (e3u_0(:,:,1) * ( 1._wp + r3u_f(:,:) * umask(:,:,1) )) * puu(:,:,1,Kmm) * umask(:,:,1) 204 201 vv_b(:,:,Kaa) = e3v(:,:,1,Kaa) * pvv(:,:,1,Kaa) * vmask(:,:,1) 205 vv_b(:,:,Kmm) = e3v(:,:,1,Kmm) * pvv(:,:,1,Kmm) * vmask(:,:,1)202 vv_b(:,:,Kmm) = (e3v_0(:,:,1) * ( 1._wp + r3v_f(:,:) * vmask(:,:,1))) * pvv(:,:,1,Kmm) * vmask(:,:,1) 206 203 DO jk = 2, jpkm1 207 204 uu_b(:,:,Kaa) = uu_b(:,:,Kaa) + e3u(:,:,jk,Kaa) * puu(:,:,jk,Kaa) * umask(:,:,jk) 208 uu_b(:,:,Kmm) = uu_b(:,:,Kmm) + e3u(:,:,jk,Kmm) * puu(:,:,jk,Kmm) * umask(:,:,jk)205 uu_b(:,:,Kmm) = uu_b(:,:,Kmm) + (e3u_0(:,:,jk) * ( 1._wp + r3u_f(:,:) * umask(:,:,jk) )) * puu(:,:,jk,Kmm) * umask(:,:,jk) 209 206 vv_b(:,:,Kaa) = vv_b(:,:,Kaa) + e3v(:,:,jk,Kaa) * pvv(:,:,jk,Kaa) * vmask(:,:,jk) 210 vv_b(:,:,Kmm) = vv_b(:,:,Kmm) + e3v(:,:,jk,Kmm) * pvv(:,:,jk,Kmm) * vmask(:,:,jk)207 vv_b(:,:,Kmm) = vv_b(:,:,Kmm) + (e3v_0(:,:,jk) * ( 1._wp + r3v_f(:,:) * vmask(:,:,jk) )) * pvv(:,:,jk,Kmm) * vmask(:,:,jk) 211 208 END DO 212 209 uu_b(:,:,Kaa) = uu_b(:,:,Kaa) * r1_hu(:,:,Kaa) 213 210 vv_b(:,:,Kaa) = vv_b(:,:,Kaa) * r1_hv(:,:,Kaa) 214 uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * r1_hu(:,:,Kmm)215 vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * r1_hv(:,:,Kmm)211 uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * (r1_hu_0(:,:)/( 1._wp + r3u_f(:,:) )) 212 vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * (r1_hv_0(:,:)/( 1._wp + r3v_f(:,:) )) 216 213 ! 217 214 IF( .NOT.ln_dynspg_ts ) THEN ! output the barotropic currents … … 235 232 236 233 !!========================================================================= 237 END MODULE dynatf qco234 END MODULE dynatf_qco -
NEMO/trunk/src/OCE/DYN/dynldf_lap_blp.F90
r13497 r14053 5 5 !!====================================================================== 6 6 !! History : 3.7 ! 2014-01 (G. Madec, S. Masson) Original code, re-entrant laplacian 7 !! 4.0 ! 2020-04 (A. Nasser, G. Madec) Add symmetric mixing tensor 7 8 !!---------------------------------------------------------------------- 8 9 … … 19 20 USE in_out_manager ! I/O manager 20 21 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 21 22 USE lib_mpp 23 22 24 IMPLICIT NONE 23 25 PRIVATE … … 47 49 !! 48 50 !! ** Action : - pu_rhs, pv_rhs increased by the harmonic operator applied on pu, pv. 51 !! 52 !! Reference : S.Griffies, R.Hallberg 2000 Mon.Wea.Rev., DOI:/ 49 53 !!---------------------------------------------------------------------- 50 54 INTEGER , INTENT(in ) :: kt ! ocean time-step index … … 57 61 REAL(wp) :: zsign ! local scalars 58 62 REAL(wp) :: zua, zva ! local scalars 59 REAL(wp), DIMENSION(jpi,jpj) :: zcur, zdiv 63 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zcur, zdiv 64 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zten, zshe ! tension (diagonal) and shearing (anti-diagonal) terms 60 65 !!---------------------------------------------------------------------- 61 66 ! … … 70 75 ENDIF 71 76 ! 72 ! ! =============== 73 DO jk = 1, jpkm1 ! Horizontal slab 74 ! ! =============== 75 DO_2D( 0, 1, 0, 1 ) 76 ! ! ahm * e3 * curl (computed from 1 to jpim1/jpjm1) 77 zcur(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) * e3f(ji-1,jj-1,jk) * r1_e1e2f(ji-1,jj-1) & ! ahmf already * by fmask 78 & * ( e2v(ji ,jj-1) * pv(ji ,jj-1,jk) - e2v(ji-1,jj-1) * pv(ji-1,jj-1,jk) & 79 & - e1u(ji-1,jj ) * pu(ji-1,jj ,jk) + e1u(ji-1,jj-1) * pu(ji-1,jj-1,jk) ) 80 ! ! ahm * div (computed from 2 to jpi/jpj) 81 zdiv(ji,jj) = ahmt(ji,jj,jk) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kbb) & ! ahmt already * by tmask 82 & * ( e2u(ji,jj)*e3u(ji,jj,jk,Kbb) * pu(ji,jj,jk) - e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kbb) * pu(ji-1,jj,jk) & 83 & + e1v(ji,jj)*e3v(ji,jj,jk,Kbb) * pv(ji,jj,jk) - e1v(ji,jj-1)*e3v(ji,jj-1,jk,Kbb) * pv(ji,jj-1,jk) ) 84 END_2D 77 SELECT CASE( nn_dynldf_typ ) 78 ! 79 CASE ( np_typ_rot ) !== Vorticity-Divergence operator ==! 85 80 ! 86 DO_2D( 0, 0, 0, 0 ) ! - curl( curl) + grad( div ) 87 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * umask(ji,jj,jk) * ( & ! * by umask is mandatory for dyn_ldf_blp use 88 & - ( zcur(ji ,jj) - zcur(ji,jj-1) ) * r1_e2u(ji,jj) / e3u(ji,jj,jk,Kmm) & 89 & + ( zdiv(ji+1,jj) - zdiv(ji,jj ) ) * r1_e1u(ji,jj) ) 81 ALLOCATE( zcur(jpi,jpj) , zdiv(jpi,jpj) ) 82 ! 83 DO jk = 1, jpkm1 ! Horizontal slab 84 ! 85 DO_2D( 0, 1, 0, 1 ) 86 ! ! ahm * e3 * curl (computed from 1 to jpim1/jpjm1) 87 zcur(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) * e3f(ji-1,jj-1,jk) * r1_e1e2f(ji-1,jj-1) & ! ahmf already * by fmask 88 & * ( e2v(ji ,jj-1) * pv(ji ,jj-1,jk) - e2v(ji-1,jj-1) * pv(ji-1,jj-1,jk) & 89 & - e1u(ji-1,jj ) * pu(ji-1,jj ,jk) + e1u(ji-1,jj-1) * pu(ji-1,jj-1,jk) ) 90 ! ! ahm * div (computed from 2 to jpi/jpj) 91 zdiv(ji,jj) = ahmt(ji,jj,jk) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kbb) & ! ahmt already * by tmask 92 & * ( e2u(ji,jj)*e3u(ji,jj,jk,Kbb) * pu(ji,jj,jk) - e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kbb) * pu(ji-1,jj,jk) & 93 & + e1v(ji,jj)*e3v(ji,jj,jk,Kbb) * pv(ji,jj,jk) - e1v(ji,jj-1)*e3v(ji,jj-1,jk,Kbb) * pv(ji,jj-1,jk) ) 94 END_2D 95 ! 96 DO_2D( 0, 0, 0, 0 ) ! - curl( curl) + grad( div ) 97 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * umask(ji,jj,jk) * ( & ! * by umask is mandatory for dyn_ldf_blp use 98 & - ( zcur(ji ,jj) - zcur(ji,jj-1) ) * r1_e2u(ji,jj) / e3u(ji,jj,jk,Kmm) & 99 & + ( zdiv(ji+1,jj) - zdiv(ji,jj ) ) * r1_e1u(ji,jj) ) 90 100 ! 91 pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zsign * vmask(ji,jj,jk) * ( & ! * by vmask is mandatory for dyn_ldf_blp use 92 & ( zcur(ji,jj ) - zcur(ji-1,jj) ) * r1_e1v(ji,jj) / e3v(ji,jj,jk,Kmm) & 93 & + ( zdiv(ji,jj+1) - zdiv(ji ,jj) ) * r1_e2v(ji,jj) ) 94 END_2D 95 ! ! =============== 96 END DO ! End of slab 97 ! ! =============== 101 pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zsign * vmask(ji,jj,jk) * ( & ! * by vmask is mandatory for dyn_ldf_blp use 102 & ( zcur(ji,jj ) - zcur(ji-1,jj) ) * r1_e1v(ji,jj) / e3v(ji,jj,jk,Kmm) & 103 & + ( zdiv(ji,jj+1) - zdiv(ji ,jj) ) * r1_e2v(ji,jj) ) 104 END_2D 105 ! 106 END DO ! End of slab 107 ! 108 DEALLOCATE( zcur , zdiv ) 109 ! 110 CASE ( np_typ_sym ) !== Symmetric operator ==! 111 ! 112 ALLOCATE( zten(jpi,jpj) , zshe(jpi,jpj) ) 113 ! 114 DO jk = 1, jpkm1 ! Horizontal slab 115 ! 116 DO_2D( 0, 1, 0, 1 ) 117 ! ! shearing stress component (F-point) NB : ahmf has already been multiplied by fmask 118 zshe(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) & 119 & * ( e1f(ji-1,jj-1) * r1_e2f(ji-1,jj-1) & 120 & * ( pu(ji-1,jj ,jk) * r1_e1u(ji-1,jj ) - pu(ji-1,jj-1,jk) * r1_e1u(ji-1,jj-1) ) & 121 & + e2f(ji-1,jj-1) * r1_e1f(ji-1,jj-1) & 122 & * ( pv(ji ,jj-1,jk) * r1_e2v(ji ,jj-1) - pv(ji-1,jj-1,jk) * r1_e2v(ji-1,jj-1) ) ) 123 ! ! tension stress component (T-point) NB : ahmt has already been multiplied by tmask 124 zten(ji,jj) = ahmt(ji,jj,jk) & 125 & * ( e2t(ji,jj) * r1_e1t(ji,jj) & 126 & * ( pu(ji,jj,jk) * r1_e2u(ji,jj) - pu(ji-1,jj,jk) * r1_e2u(ji-1,jj) ) & 127 & - e1t(ji,jj) * r1_e2t(ji,jj) & 128 & * ( pv(ji,jj,jk) * r1_e1v(ji,jj) - pv(ji,jj-1,jk) * r1_e1v(ji,jj-1) ) ) 129 END_2D 130 ! 131 DO_2D( 0, 0, 0, 0 ) 132 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) & 133 & * ( ( zten(ji+1,jj ) * e2t(ji+1,jj )*e2t(ji+1,jj ) * e3t(ji+1,jj ,jk,Kmm) & 134 & - zten(ji ,jj ) * e2t(ji ,jj )*e2t(ji ,jj ) * e3t(ji ,jj ,jk,Kmm) ) * r1_e2u(ji,jj) & 135 & + ( zshe(ji ,jj ) * e1f(ji ,jj )*e1f(ji ,jj ) * e3f(ji ,jj ,jk) & 136 & - zshe(ji ,jj-1) * e1f(ji ,jj-1)*e1f(ji ,jj-1) * e3f(ji ,jj-1,jk) ) * r1_e1u(ji,jj) ) 137 ! 138 pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zsign * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) & 139 & * ( ( zshe(ji ,jj ) * e2f(ji ,jj )*e2f(ji ,jj ) * e3f(ji ,jj ,jk) & 140 & - zshe(ji-1,jj ) * e2f(ji-1,jj )*e2f(ji-1,jj ) * e3f(ji-1,jj ,jk) ) * r1_e2v(ji,jj) & 141 & - ( zten(ji ,jj+1) * e1t(ji ,jj+1)*e1t(ji ,jj+1) * e3t(ji ,jj+1,jk,Kmm) & 142 & - zten(ji ,jj ) * e1t(ji ,jj )*e1t(ji ,jj ) * e3t(ji ,jj ,jk,Kmm) ) * r1_e1v(ji,jj) ) 143 ! 144 END_2D 145 ! 146 END DO 147 ! 148 DEALLOCATE( zten , zshe ) 149 ! 150 END SELECT 98 151 ! 99 152 END SUBROUTINE dyn_ldf_lap -
NEMO/trunk/src/OCE/DYN/dynspg_ts.F90
r13970 r14053 306 306 ENDIF 307 307 ! 308 ! != Add atmospheric pressureforcing =!309 ! ! ------------------ ----------------!310 IF( ln_bt_fw ) THEN ! Add wind forcing308 ! != Add wind forcing =! 309 ! ! ------------------ ! 310 IF( ln_bt_fw ) THEN 311 311 DO_2D( 0, 0, 0, 0 ) 312 312 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) … … 386 386 ! 387 387 IF( ln_linssh ) THEN ! mid-step ocean depth is fixed (hup2_e=hu_n=hu_0) 388 zhup2_e(:,:) = hu (:,:,Kmm)389 zhvp2_e(:,:) = hv (:,:,Kmm)390 zhtp2_e(:,:) = ht (:,:)391 ENDIF 392 ! 393 IF (ln_bt_fw) THEN! FORWARD integration: start from NOW fields394 sshn_e(:,:) = pssh (:,:,Kmm)388 zhup2_e(:,:) = hu_0(:,:) 389 zhvp2_e(:,:) = hv_0(:,:) 390 zhtp2_e(:,:) = ht_0(:,:) 391 ENDIF 392 ! 393 IF( ln_bt_fw ) THEN ! FORWARD integration: start from NOW fields 394 sshn_e(:,:) = pssh (:,:,Kmm) 395 395 un_e (:,:) = puu_b(:,:,Kmm) 396 396 vn_e (:,:) = pvv_b(:,:,Kmm) … … 401 401 hvr_e (:,:) = r1_hv(:,:,Kmm) 402 402 ELSE ! CENTRED integration: start from BEFORE fields 403 sshn_e(:,:) = pssh (:,:,Kbb)403 sshn_e(:,:) = pssh (:,:,Kbb) 404 404 un_e (:,:) = puu_b(:,:,Kbb) 405 405 vn_e (:,:) = pvv_b(:,:,Kbb) … … 412 412 ! 413 413 ! Initialize sums: 414 puu_b 415 pvv_b 414 puu_b (:,:,Kaa) = 0._wp ! After barotropic velocities (or transport if flux form) 415 pvv_b (:,:,Kaa) = 0._wp 416 416 pssh (:,:,Kaa) = 0._wp ! Sum for after averaged sea level 417 un_adv(:,:) = 0._wp ! Sum for now transport issued from ts loop418 vn_adv(:,:) = 0._wp417 un_adv(:,:) = 0._wp ! Sum for now transport issued from ts loop 418 vn_adv(:,:) = 0._wp 419 419 ! 420 420 IF( ln_wd_dl ) THEN … … 475 475 ! 476 476 ! ! ocean u- and v-depth at mid-step (separate DO-loops remove the need of a lbc_lnk) 477 #if defined key_qcoTest_FluxForm 478 ! ! 'key_qcoTest_FluxForm' : simple ssh average 479 DO_2D( 1, 1, 1, 0 ) ! not jpi-column 480 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * ( zsshp2_e(ji,jj) + zsshp2_e(ji+1,jj ) ) * ssumask(ji,jj) 481 END_2D 482 DO_2D( 1, 0, 1, 1 ) 483 zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * ( zsshp2_e(ji,jj) + zsshp2_e(ji ,jj+1) ) * ssvmask(ji,jj) 484 END_2D 485 #else 486 ! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 477 487 DO_2D( 1, 1, 1, 0 ) ! not jpi-column 478 488 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj) & … … 485 495 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 486 496 END_2D 497 #endif 487 498 ! 488 499 ENDIF … … 540 551 ! 541 552 ! Sea Surface Height at u-,v-points (vvl case only) 542 IF( .NOT.ln_linssh ) THEN 553 IF( .NOT.ln_linssh ) THEN 554 #if defined key_qcoTest_FluxForm 555 ! ! 'key_qcoTest_FluxForm' : simple ssh average 556 DO_2D( 1, 1, 1, 0 ) 557 zsshu_a(ji,jj) = r1_2 * ( ssha_e(ji,jj) + ssha_e(ji+1,jj ) ) * ssumask(ji,jj) 558 END_2D 559 DO_2D( 1, 0, 1, 1 ) 560 zsshv_a(ji,jj) = r1_2 * ( ssha_e(ji,jj) + ssha_e(ji ,jj+1) ) * ssvmask(ji,jj) 561 END_2D 562 #else 543 563 DO_2D( 0, 0, 0, 0 ) 544 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 545 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 546 & + e1e2t(ji+1,jj ) * ssha_e(ji+1,jj ) ) 547 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 548 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 549 & + e1e2t(ji ,jj+1) * ssha_e(ji ,jj+1) ) 550 END_2D 551 ENDIF 564 zsshu_a(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 565 & + e1e2t(ji+1,jj ) * ssha_e(ji+1,jj ) ) * ssumask(ji,jj) 566 zsshv_a(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 567 & + e1e2t(ji ,jj+1) * ssha_e(ji ,jj+1) ) * ssvmask(ji,jj) 568 END_2D 569 #endif 570 ENDIF 552 571 ! 553 572 ! Half-step back interpolation of SSH for surface pressure computation at step jit+1/2 … … 624 643 ! ! hu_e, hv_e hold depth at jn, zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 625 644 ! ! backward interpolated depth used in spg terms at jn+1/2 645 #if defined key_qcoTest_FluxForm 646 ! ! 'key_qcoTest_FluxForm' : simple ssh average 647 zhu_bck = hu_0(ji,jj) + r1_2 * ( zsshp2_e(ji,jj) + zsshp2_e(ji+1,jj ) ) * ssumask(ji,jj) 648 zhv_bck = hv_0(ji,jj) + r1_2 * ( zsshp2_e(ji,jj) + zsshp2_e(ji ,jj+1) ) * ssvmask(ji,jj) 649 #else 626 650 zhu_bck = hu_0(ji,jj) + r1_2*r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 627 651 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 628 652 zhv_bck = hv_0(ji,jj) + r1_2*r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 629 653 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 654 #endif 630 655 ! ! inverse depth at jn+1 631 656 z1_hu = ssumask(ji,jj) / ( hu_0(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj) ) … … 646 671 IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 647 672 DO_2D( 0, 0, 0, 0 ) 648 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rDt_e * zCdU_u(ji,jj) * hur_e(ji,jj))649 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rDt_e * zCdU_v(ji,jj) * hvr_e(ji,jj))673 ua_e(ji,jj) = ua_e(ji,jj) / ( 1._wp - rDt_e * zCdU_u(ji,jj) * hur_e(ji,jj) ) 674 va_e(ji,jj) = va_e(ji,jj) / ( 1._wp - rDt_e * zCdU_v(ji,jj) * hvr_e(ji,jj) ) 650 675 END_2D 651 676 ENDIF 652 677 653 678 IF( .NOT.ln_linssh ) THEN !* Update ocean depth (variable volume case only) 654 hu_e (2:jpim1,2:jpjm1) = hu_0(2:jpim1,2:jpjm1) + zsshu_a(2:jpim1,2:jpjm1)655 hur_e(2:jpim1,2:jpjm1) = ssumask(2:jpim1,2:jpjm1) / ( hu_e(2:jpim1,2:jpjm1) + 1._wp - ssumask(2:jpim1,2:jpjm1))656 hv_e (2:jpim1,2:jpjm1) = hv_0(2:jpim1,2:jpjm1) + zsshv_a(2:jpim1,2:jpjm1)657 hvr_e(2:jpim1,2:jpjm1) = ssvmask(2:jpim1,2:jpjm1) / ( hv_e(2:jpim1,2:jpjm1) + 1._wp - ssvmask(2:jpim1,2:jpjm1))679 hu_e (2:jpim1,2:jpjm1) = hu_0(2:jpim1,2:jpjm1) + zsshu_a(2:jpim1,2:jpjm1) 680 hur_e(2:jpim1,2:jpjm1) = ssumask(2:jpim1,2:jpjm1) / ( hu_e(2:jpim1,2:jpjm1) + 1._wp - ssumask(2:jpim1,2:jpjm1) ) 681 hv_e (2:jpim1,2:jpjm1) = hv_0(2:jpim1,2:jpjm1) + zsshv_a(2:jpim1,2:jpjm1) 682 hvr_e(2:jpim1,2:jpjm1) = ssvmask(2:jpim1,2:jpjm1) / ( hv_e(2:jpim1,2:jpjm1) + 1._wp - ssvmask(2:jpim1,2:jpjm1) ) 658 683 ENDIF 659 684 ! … … 743 768 ELSE 744 769 ! At this stage, pssh(:,:,:,Krhs) has been corrected: compute new depths at velocity points 770 #if defined key_qcoTest_FluxForm 771 ! ! 'key_qcoTest_FluxForm' : simple ssh average 745 772 DO_2D( 1, 0, 1, 0 ) 746 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 747 & * ( e1e2t(ji ,jj) * pssh(ji ,jj,Kaa) & 748 & + e1e2t(ji+1,jj) * pssh(ji+1,jj,Kaa) ) 749 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 750 & * ( e1e2t(ji,jj ) * pssh(ji,jj ,Kaa) & 751 & + e1e2t(ji,jj+1) * pssh(ji,jj+1,Kaa) ) 752 END_2D 773 zsshu_a(ji,jj) = r1_2 * ( pssh(ji,jj,Kaa) + pssh(ji+1,jj ,Kaa) ) * ssumask(ji,jj) 774 zsshv_a(ji,jj) = r1_2 * ( pssh(ji,jj,Kaa) + pssh(ji ,jj+1,Kaa) ) * ssvmask(ji,jj) 775 END_2D 776 #else 777 DO_2D( 1, 0, 1, 0 ) 778 zsshu_a(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj) * pssh(ji ,jj,Kaa) & 779 & + e1e2t(ji+1,jj) * pssh(ji+1,jj,Kaa) ) * ssumask(ji,jj) 780 zsshv_a(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * pssh(ji,jj ,Kaa) & 781 & + e1e2t(ji,jj+1) * pssh(ji,jj+1,Kaa) ) * ssvmask(ji,jj) 782 END_2D 783 #endif 753 784 CALL lbc_lnk_multi( 'dynspg_ts', zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 754 785 ! 755 786 DO jk=1,jpkm1 756 puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_Dt_b 757 pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_Dt_b 787 puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) & 788 & * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_Dt_b 789 pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) & 790 & * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_Dt_b 758 791 END DO 759 792 ! Save barotropic velocities not transport: … … 899 932 IF( TRIM(cdrw) == 'READ' ) THEN ! Read/initialise 900 933 ! ! --------------- 901 IF( ln_rstart .AND. ln_bt_fw .AND. (.NOT.l_1st_euler)) THEN !* Read the restart file934 IF( ln_rstart .AND. ln_bt_fw .AND. .NOT.l_1st_euler ) THEN !* Read the restart file 902 935 CALL iom_get( numror, jpdom_auto, 'ub2_b' , ub2_b (:,:), cd_type = 'U', psgn = -1._wp ) 903 936 CALL iom_get( numror, jpdom_auto, 'vb2_b' , vb2_b (:,:), cd_type = 'V', psgn = -1._wp ) … … 1060 1093 !! although they should be updated in the variable volume case. Not a big approximation. 1061 1094 !! To remove this approximation, copy lines below inside barotropic loop 1062 !! and update depths at T- F points (ht and zhf resp.) at each barotropic time step1095 !! and update depths at T- points (ht) at each barotropic time step 1063 1096 !! 1064 1097 !! Compute zwz = f / ( height of the water colomn ) … … 1067 1100 INTEGER :: ji ,jj, jk ! dummy loop indices 1068 1101 REAL(wp) :: z1_ht 1069 REAL(wp), DIMENSION(jpi,jpj) :: zhf1070 1102 !!---------------------------------------------------------------------- 1071 1103 ! 1072 1104 SELECT CASE( nvor_scheme ) 1073 CASE( np_EEN ) != EEN scheme using e3f energy & enstrophy scheme1074 SELECT CASE( nn_e en_e3f )!* ff_f/e3 at F-point1105 CASE( np_EEN, np_ENE, np_ENS , np_MIX ) != schemes using the same e3f definition 1106 SELECT CASE( nn_e3f_typ ) !* ff_f/e3 at F-point 1075 1107 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 1076 DO_2D( 1, 0, 1, 0 )1077 zwz(ji,jj) = ( ht(ji ,jj+1) + ht(ji+1,jj+1) +&1078 & ht(ji ,jj ) + ht(ji+1,jj )) * 0.25_wp1108 DO_2D( 0, 0, 0, 0 ) 1109 zwz(ji,jj) = ( ht(ji,jj+1) + ht(ji+1,jj+1) & 1110 & + ht(ji,jj ) + ht(ji+1,jj ) ) * 0.25_wp 1079 1111 IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 1080 1112 END_2D 1081 1113 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 1082 DO_2D( 1, 0, 1, 0 )1083 zwz(ji,jj) = ( ht (ji ,jj+1) + ht(ji+1,jj+1) &1084 & + ht (ji ,jj ) + ht(ji+1,jj ) ) &1085 & / ( MAX( 1._wp, ssmask(ji,jj+1) + ssmask(ji+1,jj+1) &1086 & + ssmask(ji ,jj ) + ssmask(ji+1,jj )) )1114 DO_2D( 0, 0, 0, 0 ) 1115 zwz(ji,jj) = ( ht(ji,jj+1) + ht(ji+1,jj+1) & 1116 & + ht(ji,jj ) + ht(ji+1,jj ) ) & 1117 & / ( MAX(ssmask(ji,jj+1) + ssmask(ji+1,jj+1) & 1118 & + ssmask(ji,jj ) + ssmask(ji+1,jj ) , 1._wp ) ) 1087 1119 IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 1088 1120 END_2D 1089 1121 END SELECT 1090 1122 CALL lbc_lnk( 'dynspg_ts', zwz, 'F', 1._wp ) 1091 ! 1092 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1123 END SELECT 1124 ! 1125 SELECT CASE( nvor_scheme ) 1126 CASE( np_EEN ) 1127 ! 1128 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1093 1129 DO_2D( 0, 1, 0, 1 ) 1094 1130 ftne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1) … … 1098 1134 END_2D 1099 1135 ! 1100 CASE( np_EET ) != EEN scheme using e3t energy conserving scheme1101 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ;ftsw(1,:) = 0._wp1136 CASE( np_EET ) != EEN scheme using e3t energy conserving scheme 1137 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1102 1138 DO_2D( 0, 1, 0, 1 ) 1103 1139 z1_ht = ssmask(ji,jj) / ( ht(ji,jj) + 1._wp - ssmask(ji,jj) ) … … 1108 1144 END_2D 1109 1145 ! 1110 CASE( np_ENE, np_ENS , np_MIX ) != all other schemes (ENE, ENS, MIX) except ENT !1111 !1112 zwz(:,:) = 0._wp1113 zhf(:,:) = 0._wp1114 1115 !!gm assume 0 in both cases (which is almost surely WRONG ! ) as hvatf has been removed1116 !!gm A priori a better value should be something like :1117 !!gm zhf(i,j) = masked sum of ht(i,j) , ht(i+1,j) , ht(i,j+1) , (i+1,j+1)1118 !!gm divided by the sum of the corresponding mask1119 !!gm1120 !!1121 IF( .NOT.ln_sco ) THEN1122 1123 !!gm agree the JC comment : this should be done in a much clear way1124 1125 ! JC: It not clear yet what should be the depth at f-points over land in z-coordinate case1126 ! Set it to zero for the time being1127 ! IF( rn_hmin < 0._wp ) THEN ; jk = - INT( rn_hmin ) ! from a nb of level1128 ! ELSE ; jk = MINLOC( gdepw_0, mask = gdepw_0 > rn_hmin, dim = 1 ) ! from a depth1129 ! ENDIF1130 ! zhf(:,:) = gdepw_0(:,:,jk+1)1131 !1132 ELSE1133 !1134 !zhf(:,:) = hbatf(:,:)1135 DO_2D( 1, 0, 1, 0 )1136 zhf(ji,jj) = ( ht_0 (ji,jj ) + ht_0 (ji+1,jj ) &1137 & + ht_0 (ji,jj+1) + ht_0 (ji+1,jj+1) ) &1138 & / MAX( ssmask(ji,jj ) + ssmask(ji+1,jj ) &1139 & + ssmask(ji,jj+1) + ssmask(ji+1,jj+1) , 1._wp )1140 END_2D1141 ENDIF1142 !1143 DO jj = 1, jpjm11144 zhf(:,jj) = zhf(:,jj) * (1._wp- umask(:,jj,1) * umask(:,jj+1,1))1145 END DO1146 !1147 DO jk = 1, jpkm11148 DO jj = 1, jpjm11149 zhf(:,jj) = zhf(:,jj) + e3f(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk)1150 END DO1151 END DO1152 CALL lbc_lnk( 'dynspg_ts', zhf, 'F', 1._wp )1153 ! JC: TBC. hf should be greater than 01154 DO_2D( 1, 1, 1, 1 )1155 IF( zhf(ji,jj) /= 0._wp ) zwz(ji,jj) = 1._wp / zhf(ji,jj)1156 END_2D1157 zwz(:,:) = ff_f(:,:) * zwz(:,:)1158 1146 END SELECT 1159 1147 1160 1148 END SUBROUTINE dyn_cor_2d_init 1161 1162 1149 1163 1150 … … 1353 1340 END SUBROUTINE wad_spg 1354 1341 1355 1356 1342 1357 1343 SUBROUTINE dyn_drg_init( Kbb, Kmm, puu, pvv, puu_b ,pvv_b, pu_RHSi, pv_RHSi, pCdU_u, pCdU_v ) -
NEMO/trunk/src/OCE/DYN/dynvor.F90
r14007 r14053 21 21 !! - ! 2018-03 (G. Madec) add two new schemes (ln_dynvor_enT and ln_dynvor_eet) 22 22 !! - ! 2018-04 (G. Madec) add pre-computed gradient for metric term calculation 23 !! 4.x ! 2020-03 (G. Madec, A. Nasser) make ln_dynvor_msk truly efficient on relative vorticity 23 24 !! 4.2 ! 2020-12 (G. Madec, E. Clementi) add vortex force trends (ln_vortex_force=T) 24 25 !!---------------------------------------------------------------------- … … 26 27 !!---------------------------------------------------------------------- 27 28 !! dyn_vor : Update the momentum trend with the vorticity trend 29 !! vor_enT : energy conserving scheme at T-pt (ln_dynvor_enT=T) 30 !! vor_ene : energy conserving scheme (ln_dynvor_ene=T) 28 31 !! vor_ens : enstrophy conserving scheme (ln_dynvor_ens=T) 29 !! vor_ene : energy conserving scheme (ln_dynvor_ene=T)30 32 !! vor_een : energy and enstrophy conserving (ln_dynvor_een=T) 33 !! vor_eeT : energy conserving at T-pt (ln_dynvor_eeT=T) 31 34 !! dyn_vor_init : set and control of the different vorticity option 32 35 !!---------------------------------------------------------------------- … … 58 61 LOGICAL, PUBLIC :: ln_dynvor_eeT !: t-point energy conserving scheme (EET) 59 62 LOGICAL, PUBLIC :: ln_dynvor_een !: energy & enstrophy conserving scheme (EEN) 60 INTEGER, PUBLIC :: nn_een_e3f !: e3f=masked averaging of e3t divided by 4 (=0) or by the sum of mask (=1)61 63 LOGICAL, PUBLIC :: ln_dynvor_mix !: mixed scheme (MIX) 62 64 LOGICAL, PUBLIC :: ln_dynvor_msk !: vorticity multiplied by fmask (=T) or not (=F) (all vorticity schemes) 65 INTEGER, PUBLIC :: nn_e3f_typ !: e3f=masked averaging of e3t divided by 4 (=0) or by the sum of mask (=1) 63 66 64 67 INTEGER, PUBLIC :: nvor_scheme !: choice of the type of advection scheme … … 81 84 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: di_e2u_2 ! = di(e2u)/2 used in T-point metric term calculation 82 85 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1v_2 ! = dj(e1v)/2 - - - - 83 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: di_e2v_2e1e2f ! = di(e2v)/(2*e1e2f) used in F-point metric term calculation 84 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1u_2e1e2f ! = dj(e1u)/(2*e1e2f) - - - - 86 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: di_e2v_2e1e2f ! = di(e2u)/(2*e1e2f) used in F-point metric term calculation 87 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1u_2e1e2f ! = dj(e1v)/(2*e1e2f) - - - - 88 ! 89 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: e3f_0vor ! e3f used in EEN, ENE and ENS cases (key_qco only) 85 90 86 91 REAL(wp) :: r1_4 = 0.250_wp ! =1/4 … … 235 240 INTEGER :: ji, jj, jk ! dummy loop indices 236 241 REAL(wp) :: zx1, zy1, zx2, zy2 ! local scalars 237 REAL(wp), DIMENSION(jpi,jpj) 238 REAL(wp), DIMENSION(jpi,jpj,jpkm1) :: zwz! 3D workspace, jpkm1 -> avoid lbc_lnk on jpk that is not defined242 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwt ! 2D workspace 243 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zwz ! 3D workspace, jpkm1 -> avoid lbc_lnk on jpk that is not defined 239 244 !!---------------------------------------------------------------------- 240 245 ! … … 246 251 ! 247 252 ! 248 SELECT CASE( kvor ) !== volume weighted vorticity considered ==! 249 CASE ( np_RVO ) !* relative vorticity 250 DO jk = 1, jpkm1 ! Horizontal slab 253 SELECT CASE( kvor ) !== relative vorticity considered ==! 254 ! 255 CASE ( np_RVO , np_CRV ) !* relative vorticity at f-point is used 256 ALLOCATE( zwz(jpi,jpj,jpk) ) 257 DO jk = 1, jpkm1 ! Horizontal slab 251 258 DO_2D( 1, 0, 1, 0 ) 252 259 zwz(ji,jj,jk) = ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) & 253 260 & - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj) 254 261 END_2D 255 IF( ln_dynvor_msk ) THEN ! mask /unmaskrelative vorticity262 IF( ln_dynvor_msk ) THEN ! mask relative vorticity 256 263 DO_2D( 1, 0, 1, 0 ) 257 264 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk) … … 259 266 ENDIF 260 267 END DO 261 262 268 CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp ) 263 264 CASE ( np_CRV ) !* Coriolis + relative vorticity 265 DO jk = 1, jpkm1 ! Horizontal slab 266 DO_2D( 1, 0, 1, 0 ) ! relative vorticity 267 zwz(ji,jj,jk) = ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) & 268 & - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj) 269 END_2D 270 IF( ln_dynvor_msk ) THEN ! mask/unmask relative vorticity 271 DO_2D( 1, 0, 1, 0 ) 272 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk) 273 END_2D 274 ENDIF 275 END DO 276 277 CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp ) 278 269 ! 279 270 END SELECT 280 271 281 272 ! ! =============== 282 273 DO jk = 1, jpkm1 ! Horizontal slab 283 !! ===============284 274 ! ! =============== 275 ! 285 276 SELECT CASE( kvor ) !== volume weighted vorticity considered ==! 277 ! 286 278 CASE ( np_COR ) !* Coriolis (planetary vorticity) 287 279 zwt(:,:) = ff_t(:,:) * e1e2t(:,:)*e3t(:,:,jk,Kmm) 288 280 CASE ( np_RVO ) !* relative vorticity 289 281 DO_2D( 0, 1, 0, 1 ) 290 zwt(ji,jj) = r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) &291 & + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) )&292 & 282 zwt(ji,jj) = r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) & 283 & + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) & 284 & * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm) 293 285 END_2D 294 286 CASE ( np_MET ) !* metric term 295 287 DO_2D( 0, 1, 0, 1 ) 296 zwt(ji,jj) = ( ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj) &297 & - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj) ) &298 & 288 zwt(ji,jj) = ( ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj) & 289 & - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj) ) & 290 & * e3t(ji,jj,jk,Kmm) 299 291 END_2D 300 292 CASE ( np_CRV ) !* Coriolis + relative vorticity 301 293 DO_2D( 0, 1, 0, 1 ) 302 zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) &303 & + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) ) &304 & 294 zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) & 295 & + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) ) & 296 & * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm) 305 297 END_2D 306 298 CASE ( np_CME ) !* Coriolis + metric 307 299 DO_2D( 0, 1, 0, 1 ) 308 zwt(ji,jj) = ( ff_t(ji,jj) * e1e2t(ji,jj) &309 & + ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj) &310 & - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj) ) &311 & 300 zwt(ji,jj) = ( ff_t(ji,jj) * e1e2t(ji,jj) & 301 & + ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj) & 302 & - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj) ) & 303 & * e3t(ji,jj,jk,Kmm) 312 304 END_2D 313 305 CASE DEFAULT ! error 314 CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' 306 CALL ctl_stop('STOP','dyn_vor: wrong value for kvor') 315 307 END SELECT 316 308 ! … … 328 320 END DO ! End of slab 329 321 ! ! =============== 322 ! 323 SELECT CASE( kvor ) ! deallocate zwz if necessary 324 CASE ( np_RVO , np_CRV ) ; DEALLOCATE( zwz ) 325 END SELECT 326 ! 330 327 END SUBROUTINE vor_enT 331 328 … … 358 355 ! 359 356 INTEGER :: ji, jj, jk ! dummy loop indices 360 REAL(wp) :: zx1, zy1, zx2, zy2 ! local scalars357 REAL(wp) :: zx1, zy1, zx2, zy2, ze3f, zmsk ! local scalars 361 358 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz ! 2D workspace 362 359 !!---------------------------------------------------------------------- … … 380 377 & - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj) 381 378 END_2D 379 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity 380 DO_2D( 1, 0, 1, 0 ) 381 zwz(ji,jj) = zwz(ji,jj) * fmask(ji,jj,jk) 382 END_2D 383 ENDIF 382 384 CASE ( np_MET ) !* metric term 383 385 DO_2D( 1, 0, 1, 0 ) … … 390 392 & - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj) 391 393 END_2D 394 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity (NOT the Coriolis term) 395 DO_2D( 1, 0, 1, 0 ) 396 zwz(ji,jj) = ( zwz(ji,jj) - ff_f(ji,jj) ) * fmask(ji,jj,jk) + ff_f(ji,jj) 397 END_2D 398 ENDIF 392 399 CASE ( np_CME ) !* Coriolis + metric 393 400 DO_2D( 1, 0, 1, 0 ) … … 399 406 END SELECT 400 407 ! 401 IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==! 402 DO_2D( 1, 0, 1, 0 ) 403 zwz(ji,jj) = zwz(ji,jj) * fmask(ji,jj,jk) 404 END_2D 405 ENDIF 406 407 IF( ln_sco ) THEN 408 zwz(:,:) = zwz(:,:) / e3f(:,:,jk) 409 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 410 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 411 ELSE 412 zwx(:,:) = e2u(:,:) * pu(:,:,jk) 413 zwy(:,:) = e1v(:,:) * pv(:,:,jk) 414 ENDIF 408 #if defined key_qco 409 DO_2D( 1, 0, 1, 0 ) !== potential vorticity ==! (key_qco) 410 zwz(ji,jj) = zwz(ji,jj) / e3f_vor(ji,jj,jk) 411 END_2D 412 #else 413 SELECT CASE( nn_e3f_typ ) !== potential vorticity ==! 414 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 415 DO_2D( 1, 0, 1, 0 ) 416 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 417 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 418 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 419 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 420 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * 4._wp / ze3f 421 ELSE ; zwz(ji,jj) = 0._wp 422 ENDIF 423 END_2D 424 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 425 DO_2D( 1, 0, 1, 0 ) 426 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 427 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 428 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 429 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 430 zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) & 431 & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) ) 432 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * zmsk / ze3f 433 ELSE ; zwz(ji,jj) = 0._wp 434 ENDIF 435 END_2D 436 END SELECT 437 #endif 438 ! !== horizontal fluxes ==! 439 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 440 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 441 ! 415 442 ! !== compute and add the vorticity term trend =! 416 443 DO_2D( 0, 0, 0, 0 ) … … 455 482 ! 456 483 INTEGER :: ji, jj, jk ! dummy loop indices 457 REAL(wp) :: zuav, zvau ! local scalars484 REAL(wp) :: zuav, zvau, ze3f, zmsk ! local scalars 458 485 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz, zww ! 2D workspace 459 486 !!---------------------------------------------------------------------- … … 476 503 & - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj) 477 504 END_2D 505 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity 506 DO_2D( 1, 0, 1, 0 ) 507 zwz(ji,jj) = zwz(ji,jj) * fmask(ji,jj,jk) 508 END_2D 509 ENDIF 478 510 CASE ( np_MET ) !* metric term 479 511 DO_2D( 1, 0, 1, 0 ) … … 486 518 & - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj) 487 519 END_2D 520 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity (NOT the Coriolis term) 521 DO_2D( 1, 0, 1, 0 ) 522 zwz(ji,jj) = ( zwz(ji,jj) - ff_f(ji,jj) ) * fmask(ji,jj,jk) + ff_f(ji,jj) 523 END_2D 524 ENDIF 488 525 CASE ( np_CME ) !* Coriolis + metric 489 526 DO_2D( 1, 0, 1, 0 ) … … 495 532 END SELECT 496 533 ! 497 IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==! 498 DO_2D( 1, 0, 1, 0 ) 499 zwz(ji,jj) = zwz(ji,jj) * fmask(ji,jj,jk) 500 END_2D 501 ENDIF 502 ! 503 IF( ln_sco ) THEN !== horizontal fluxes ==! 504 zwz(:,:) = zwz(:,:) / e3f(:,:,jk) 505 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 506 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 507 ELSE 508 zwx(:,:) = e2u(:,:) * pu(:,:,jk) 509 zwy(:,:) = e1v(:,:) * pv(:,:,jk) 510 ENDIF 534 ! 535 #if defined key_qco 536 DO_2D( 1, 0, 1, 0 ) !== potential vorticity ==! (key_qco) 537 zwz(ji,jj) = zwz(ji,jj) / e3f_vor(ji,jj,jk) 538 END_2D 539 #else 540 SELECT CASE( nn_e3f_typ ) !== potential vorticity ==! 541 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 542 DO_2D( 1, 0, 1, 0 ) 543 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 544 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 545 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 546 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 547 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * 4._wp / ze3f 548 ELSE ; zwz(ji,jj) = 0._wp 549 ENDIF 550 END_2D 551 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 552 DO_2D( 1, 0, 1, 0 ) 553 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 554 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 555 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 556 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 557 zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) & 558 & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) ) 559 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * zmsk / ze3f 560 ELSE ; zwz(ji,jj) = 0._wp 561 ENDIF 562 END_2D 563 END SELECT 564 #endif 565 ! !== horizontal fluxes ==! 566 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 567 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 568 ! 511 569 ! !== compute and add the vorticity term trend =! 512 570 DO_2D( 0, 0, 0, 0 ) … … 566 624 ! ! =============== 567 625 ! 568 SELECT CASE( nn_een_e3f ) ! == reciprocal of e3 at F-point 626 #if defined key_qco 627 DO_2D( 1, 0, 1, 0 ) ! == reciprocal of e3 at F-point (key_qco) 628 z1_e3f(ji,jj) = 1._wp / e3f_vor(ji,jj,jk) 629 END_2D 630 #else 631 SELECT CASE( nn_e3f_typ ) ! == reciprocal of e3 at F-point 569 632 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 570 633 DO_2D( 1, 0, 1, 0 ) … … 590 653 END_2D 591 654 END SELECT 655 #endif 592 656 ! 593 657 SELECT CASE( kvor ) !== vorticity considered ==! 658 ! 594 659 CASE ( np_COR ) !* Coriolis (planetary vorticity) 595 660 DO_2D( 1, 0, 1, 0 ) … … 601 666 & - e1u(ji ,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)*z1_e3f(ji,jj) 602 667 END_2D 668 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity 669 DO_2D( 1, 0, 1, 0 ) 670 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk) 671 END_2D 672 ENDIF 603 673 CASE ( np_MET ) !* metric term 604 674 DO_2D( 1, 0, 1, 0 ) … … 612 682 & * r1_e1e2f(ji,jj) ) * z1_e3f(ji,jj) 613 683 END_2D 684 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity 685 DO_2D( 1, 0, 1, 0 ) 686 zwz(ji,jj,jk) = ( zwz(ji,jj,jk) - ff_f(ji,jj) ) * fmask(ji,jj,jk) + ff_f(ji,jj) 687 END_2D 688 ENDIF 614 689 CASE ( np_CME ) !* Coriolis + metric 615 690 DO_2D( 1, 0, 1, 0 ) … … 620 695 CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' ) 621 696 END SELECT 622 ! 623 IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==! 624 DO_2D( 1, 0, 1, 0 ) 625 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk) 626 END_2D 627 ENDIF 697 ! ! =============== 628 698 END DO ! End of slab 629 ! 699 ! ! =============== 700 ! 630 701 CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp ) 631 702 ! 703 ! ! =============== 632 704 DO jk = 1, jpkm1 ! Horizontal slab 705 ! ! =============== 633 706 ! 634 707 ! !== horizontal fluxes ==! 635 708 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 636 709 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 637 710 ! 638 711 ! !== compute and add the vorticity term trend =! 639 jj = 2 640 ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0 641 DO ji = 2, jpi ! split in 2 parts due to vector opt. 642 ztne(ji,jj) = zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) 643 ztnw(ji,jj) = zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) 644 ztse(ji,jj) = zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) 645 ztsw(ji,jj) = zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) 646 END DO 647 DO jj = 3, jpj 648 DO ji = 2, jpi ! vector opt. ok because we start at jj = 3 649 ztne(ji,jj) = zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) 650 ztnw(ji,jj) = zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) 651 ztse(ji,jj) = zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) 652 ztsw(ji,jj) = zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) 653 END DO 654 END DO 712 DO_2D( 0, 1, 0, 1 ) 713 ztne(ji,jj) = zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) 714 ztnw(ji,jj) = zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) 715 ztse(ji,jj) = zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) 716 ztsw(ji,jj) = zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) 717 END_2D 718 ! 655 719 DO_2D( 0, 0, 0, 0 ) 656 720 zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) & … … 667 731 668 732 669 670 733 SUBROUTINE vor_eeT( kt, Kmm, kvor, pu, pv, pu_rhs, pv_rhs ) 671 734 !!---------------------------------------------------------------------- … … 685 748 !! References : Arakawa and Lamb 1980, Mon. Wea. Rev., 109, 18-36 686 749 !!---------------------------------------------------------------------- 687 INTEGER , INTENT(in ) :: kt ! ocean time-step index750 INTEGER , INTENT(in ) :: kt ! ocean time-step index 688 751 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 689 INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric690 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities691 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs 752 INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric 753 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities 754 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs ! total v-trend 692 755 ! 693 756 INTEGER :: ji, jj, jk ! dummy loop indices … … 702 765 IF( kt == nit000 ) THEN 703 766 IF(lwp) WRITE(numout,*) 704 IF(lwp) WRITE(numout,*) 'dyn:vor_ee n: vorticity term: energy and enstrophy conserving scheme'767 IF(lwp) WRITE(numout,*) 'dyn:vor_eeT : vorticity term: energy and enstrophy conserving scheme' 705 768 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 706 769 ENDIF … … 722 785 & * r1_e1e2f(ji,jj) 723 786 END_2D 787 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity 788 DO_2D( 1, 0, 1, 0 ) 789 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk) 790 END_2D 791 ENDIF 724 792 CASE ( np_MET ) !* metric term 725 793 DO_2D( 1, 0, 1, 0 ) … … 733 801 & * r1_e1e2f(ji,jj) ) 734 802 END_2D 803 IF( ln_dynvor_msk ) THEN ! mask the relative vorticity 804 DO_2D( 1, 0, 1, 0 ) 805 zwz(ji,jj,jk) = ( zwz(ji,jj,jk) - ff_f(ji,jj) ) * fmask(ji,jj,jk) + ff_f(ji,jj) 806 END_2D 807 ENDIF 735 808 CASE ( np_CME ) !* Coriolis + metric 736 809 DO_2D( 1, 0, 1, 0 ) … … 742 815 END SELECT 743 816 ! 744 IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==! 745 DO_2D( 1, 0, 1, 0 ) 746 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk) 747 END_2D 748 ENDIF 749 END DO 817 ! ! =============== 818 END DO ! End of slab 819 ! ! =============== 750 820 ! 751 821 CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp ) 752 822 ! 823 ! ! =============== 753 824 DO jk = 1, jpkm1 ! Horizontal slab 754 755 ! !== horizontal fluxes ==! 825 ! ! =============== 826 ! 827 ! !== horizontal fluxes ==! 756 828 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 757 829 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 758 830 ! 759 831 ! !== compute and add the vorticity term trend =! 760 jj = 2 761 ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0 762 DO ji = 2, jpi ! split in 2 parts due to vector opt. 763 z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 764 ztne(ji,jj) = ( zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) ) * z1_e3t 765 ztnw(ji,jj) = ( zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) ) * z1_e3t 766 ztse(ji,jj) = ( zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) ) * z1_e3t 767 ztsw(ji,jj) = ( zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) ) * z1_e3t 768 END DO 769 DO jj = 3, jpj 770 DO ji = 2, jpi ! vector opt. ok because we start at jj = 3 771 z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 772 ztne(ji,jj) = ( zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) ) * z1_e3t 773 ztnw(ji,jj) = ( zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) ) * z1_e3t 774 ztse(ji,jj) = ( zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) ) * z1_e3t 775 ztsw(ji,jj) = ( zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) ) * z1_e3t 776 END DO 777 END DO 832 DO_2D( 0, 1, 0, 1 ) 833 z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 834 ztne(ji,jj) = ( zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) ) * z1_e3t 835 ztnw(ji,jj) = ( zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) ) * z1_e3t 836 ztse(ji,jj) = ( zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) ) * z1_e3t 837 ztsw(ji,jj) = ( zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) ) * z1_e3t 838 END_2D 839 ! 778 840 DO_2D( 0, 0, 0, 0 ) 779 841 zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) & … … 799 861 INTEGER :: ji, jj, jk ! dummy loop indices 800 862 INTEGER :: ioptio, ios ! local integer 863 REAL(wp) :: zmsk ! local scalars 801 864 !! 802 865 NAMELIST/namdyn_vor/ ln_dynvor_ens, ln_dynvor_ene, ln_dynvor_enT, ln_dynvor_eeT, & 803 & ln_dynvor_een, nn_e en_e3f, ln_dynvor_mix, ln_dynvor_msk866 & ln_dynvor_een, nn_e3f_typ , ln_dynvor_mix, ln_dynvor_msk 804 867 !!---------------------------------------------------------------------- 805 868 ! … … 823 886 WRITE(numout,*) ' energy conserving scheme (een using e3t) ln_dynvor_eeT = ', ln_dynvor_eeT 824 887 WRITE(numout,*) ' enstrophy and energy conserving scheme ln_dynvor_een = ', ln_dynvor_een 825 WRITE(numout,*) ' e3f = averaging /4 (=0) or /sum(tmask) (=1) nn_e en_e3f = ', nn_een_e3f888 WRITE(numout,*) ' e3f = averaging /4 (=0) or /sum(tmask) (=1) nn_e3f_typ = ', nn_e3f_typ 826 889 WRITE(numout,*) ' mixed enstrophy/energy conserving scheme ln_dynvor_mix = ', ln_dynvor_mix 827 890 WRITE(numout,*) ' masked (=T) or unmasked(=F) vorticity ln_dynvor_msk = ', ln_dynvor_msk 828 891 ENDIF 829 830 IF( ln_dynvor_msk ) CALL ctl_stop( 'dyn_vor_init: masked vorticity is not currently not available')831 892 832 893 !!gm this should be removed when choosing a unique strategy for fmask at the coast … … 891 952 ! 892 953 END SELECT 893 954 #if defined key_qco 955 SELECT CASE( nvor_scheme ) ! qco case: pre-computed a specific e3f_0 for some vorticity schemes 956 CASE( np_ENS , np_ENE , np_EEN , np_MIX ) 957 ! 958 ALLOCATE( e3f_0vor(jpi,jpj,jpk) ) 959 ! 960 SELECT CASE( nn_e3f_typ ) 961 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 962 DO_3D( 0, 0, 0, 0, 1, jpk ) 963 e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) & 964 & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) & 965 & + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) & 966 & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) * 0.25_wp 967 END_3D 968 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 969 DO_3D( 0, 0, 0, 0, 1, jpk ) 970 zmsk = (tmask(ji,jj+1,jk) +tmask(ji+1,jj+1,jk) & 971 & + tmask(ji,jj ,jk) +tmask(ji+1,jj ,jk) ) 972 ! 973 IF( zmsk /= 0._wp ) THEN 974 e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) & 975 & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) & 976 & + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) & 977 & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) / zmsk 978 ENDIF 979 END_3D 980 END SELECT 981 ! 982 CALL lbc_lnk( 'dynvor', e3f_0vor, 'F', 1._wp ) 983 ! ! insure e3f_0vor /= 0 984 WHERE( e3f_0vor(:,:,:) == 0._wp ) e3f_0vor(:,:,:) = e3f_0(:,:,:) 985 ! 986 END SELECT 987 ! 988 #endif 894 989 IF(lwp) THEN ! Print the choice 895 990 WRITE(numout,*) … … 898 993 CASE( np_ENE ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (Coriolis at F-points) (ENE)' 899 994 CASE( np_ENT ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (Coriolis at T-points) (ENT)' 995 IF( ln_dynadv_vec ) CALL ctl_warn('dyn_vor_init: ENT scheme may not work in vector form') 900 996 CASE( np_EET ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (EEN scheme using e3t) (EET)' 901 997 CASE( np_EEN ) ; WRITE(numout,*) ' ==>>> energy and enstrophy conserving scheme (EEN)' -
NEMO/trunk/src/OCE/DYN/sshwzv.F90
r13497 r14053 6 6 !! History : 3.1 ! 2009-02 (G. Madec, M. Leclair) Original code 7 7 !! 3.3 ! 2010-04 (M. Leclair, G. Madec) modified LF-RA 8 !! - ! 2010-05 (K. Mogensen, A. Weaver, M. Martin, D. Lea) Assimilation interface 9 !! - ! 2010-09 (D.Storkey and E.O'Dea) bug fixes for BDY module 10 !! 3.3 ! 2011-10 (M. Leclair) split former ssh_wzv routine and remove all vvl related work 11 !! 4.0 ! 2018-12 (A. Coward) add mixed implicit/explicit advection 12 !! 4.1 ! 2019-08 (A. Coward, D. Storkey) Rename ssh_nxt -> ssh_atf. Now only does time filtering. 8 !! - ! 2010-05 (K. Mogensen, A. Weaver, M. Martin, D. Lea) Assimilation interface 9 !! - ! 2010-09 (D.Storkey and E.O'Dea) bug fixes for BDY module 10 !! 3.3 ! 2011-10 (M. Leclair) split former ssh_wzv routine and remove all vvl related work 11 !! 4.0 ! 2018-12 (A. Coward) add mixed implicit/explicit advection 12 !! 4.1 ! 2019-08 (A. Coward, D. Storkey) Rename ssh_nxt -> ssh_atf. Now only does time filtering. 13 !! - ! 2020-08 (S. Techene, G. Madec) add here ssh initiatlisation 13 14 !!---------------------------------------------------------------------- 14 15 … … 17 18 !! ssh_atf : time filter the ssh arrays 18 19 !! wzv : compute now vertical velocity 20 !! ssh_init_rst : ssh set from restart or domcfg.nc file or usr_def_istat_ssh 19 21 !!---------------------------------------------------------------------- 20 22 USE oce ! ocean dynamics and tracers variables … … 40 42 USE timing ! Timing 41 43 USE wet_dry ! Wetting/Drying flux limiting 42 44 USE usrdef_istate, ONLY : usr_def_istate_ssh ! user defined ssh initial state 45 43 46 IMPLICIT NONE 44 47 PRIVATE 45 48 46 PUBLIC ssh_nxt ! called by step.F90 47 PUBLIC wzv ! called by step.F90 48 PUBLIC wAimp ! called by step.F90 49 PUBLIC ssh_atf ! called by step.F90 49 PUBLIC ssh_nxt ! called by step.F90 50 PUBLIC wzv ! called by step.F90 51 PUBLIC wAimp ! called by step.F90 52 PUBLIC ssh_atf ! called by step.F90 53 PUBLIC ssh_init_rst ! called by domain.F90 50 54 51 55 !! * Substitutions 52 56 # include "do_loop_substitute.h90" 53 57 # include "domzgr_substitute.h90" 54 55 58 !!---------------------------------------------------------------------- 56 59 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 299 302 ! ! filtered "now" field 300 303 pssh(:,:,Kmm) = pssh(:,:,Kmm) + rn_atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 304 ! 301 305 IF( .NOT.ln_linssh ) THEN ! "now" <-- with forcing removed 302 306 zcoef = rn_atfp * rn_Dt * r1_rho0 … … 307 311 308 312 ! ice sheet coupling 309 IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) pssh(:,:,Kbb) = pssh(:,:,Kbb) - rn_atfp * rn_Dt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 313 IF( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1 ) & 314 & pssh(:,:,Kbb) = pssh(:,:,Kbb) - rn_atfp * rn_Dt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 310 315 311 316 ENDIF 312 317 ENDIF 313 318 ! 314 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=pssh(:,:,Kmm), clinfo1=' pssh(:,:,Kmm) -: ', mask1=tmask )319 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=pssh(:,:,Kmm), clinfo1=' atf - pssh(:,:,Kmm): ', mask1=tmask ) 315 320 ! 316 321 IF( ln_timing ) CALL timing_stop('ssh_atf') … … 431 436 ! 432 437 END SUBROUTINE wAimp 438 439 440 SUBROUTINE ssh_init_rst( Kbb, Kmm, Kaa ) 441 !!--------------------------------------------------------------------- 442 !! *** ROUTINE ssh_init_rst *** 443 !! 444 !! ** Purpose : ssh initialization of the sea surface height (ssh) 445 !! 446 !! ** Method : set ssh from restart or read configuration, or user_def 447 !! * ln_rstart = T 448 !! USE of IOM library to read ssh in the restart file 449 !! Leap-Frog: Kbb and Kmm are read except for l_1st_euler=T 450 !! 451 !! * otherwise 452 !! call user defined ssh or 453 !! set to -ssh_ref in wet and drying case with domcfg.nc 454 !! 455 !! NB: ssh_b/n are written by restart.F90 456 !!---------------------------------------------------------------------- 457 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 458 ! 459 INTEGER :: ji, jj, jk 460 !!---------------------------------------------------------------------- 461 ! 462 IF(lwp) THEN 463 WRITE(numout,*) 464 WRITE(numout,*) 'ssh_init_rst : ssh initialization' 465 WRITE(numout,*) '~~~~~~~~~~~~ ' 466 ENDIF 467 ! 468 ! !=============================! 469 IF( ln_rstart ) THEN !== Read the restart file ==! 470 ! !=============================! 471 ! 472 ! !* Read ssh at Kmm 473 IF(lwp) WRITE(numout,*) 474 IF(lwp) WRITE(numout,*) ' Kmm sea surface height read in the restart file' 475 CALL iom_get( numror, jpdom_auto, 'sshn' , ssh(:,:,Kmm) ) 476 ! 477 IF( l_1st_euler ) THEN !* Euler at first time-step 478 IF(lwp) WRITE(numout,*) 479 IF(lwp) WRITE(numout,*) ' Euler first time step : ssh(Kbb) = ssh(Kmm)' 480 ssh(:,:,Kbb) = ssh(:,:,Kmm) 481 ! 482 ELSE !* read ssh at Kbb 483 IF(lwp) WRITE(numout,*) 484 IF(lwp) WRITE(numout,*) ' Kbb sea surface height read in the restart file' 485 CALL iom_get( numror, jpdom_auto, 'sshb', ssh(:,:,Kbb) ) 486 ENDIF 487 ! !============================! 488 ELSE !== Initialize at "rest" ==! 489 ! !============================! 490 ! 491 IF(lwp) WRITE(numout,*) 492 IF(lwp) WRITE(numout,*) ' initialization at rest' 493 ! 494 IF( ll_wd ) THEN !* wet and dry 495 ! 496 IF( ln_read_cfg ) THEN ! read configuration : ssh_ref is read in domain_cfg file 497 !!st why ssh is not masked : i.e. ssh(:,:,Kmm) = -ssh_ref*ssmask(:,:), 498 !!st since at the 1st time step lbclnk will be applied on ssh at Kaa but not initially at Kbb and Kmm 499 ssh(:,:,Kbb) = -ssh_ref 500 ! 501 DO_2D( 1, 1, 1, 1 ) 502 IF( ht_0(ji,jj)-ssh_ref < rn_wdmin1 ) THEN ! if total depth is less than min depth 503 ssh(ji,jj,Kbb) = rn_wdmin1 - ht_0(ji,jj) 504 ENDIF 505 END_2D 506 ELSE ! user define configuration case 507 CALL usr_def_istate_ssh( tmask, ssh(:,:,Kbb) ) 508 ENDIF 509 ! 510 ELSE !* user defined configuration 511 CALL usr_def_istate_ssh( tmask, ssh(:,:,Kbb) ) 512 ! 513 ENDIF 514 ! 515 ssh(:,:,Kmm) = ssh(:,:,Kbb) !* set now values from to before ones 516 ssh(:,:,Kaa) = 0._wp 517 ENDIF 518 ! 519 END SUBROUTINE ssh_init_rst 520 433 521 !!====================================================================== 434 522 END MODULE sshwzv -
NEMO/trunk/src/OCE/IOM/iom.F90
r14039 r14053 174 174 CALL set_grid( "V", glamv, gphiv, .FALSE., .FALSE. ) 175 175 CALL set_grid( "W", glamt, gphit, .FALSE., .FALSE. ) 176 CALL set_grid( "F", glamf, gphif, .FALSE., .FALSE. ) 176 177 CALL set_grid_znl( gphit ) 177 178 ! … … 180 181 CALL iom_set_domain_attr("grid_U", area = real( e1e2u(Nis0:Nie0, Njs0:Nje0), dp)) 181 182 CALL iom_set_domain_attr("grid_V", area = real( e1e2v(Nis0:Nie0, Njs0:Nje0), dp)) 182 CALL iom_set_domain_attr("grid_W", area = real( e1e2t(Nis0:Nie0, Njs0:Nje0), dp)) 183 CALL iom_set_domain_attr("grid_W", area = REAL( e1e2t(Nis0:Nie0, Njs0:Nje0), dp)) 184 CALL iom_set_domain_attr("grid_F", area = real( e1e2f(Nis0:Nie0, Njs0:Nje0), dp)) 183 185 CALL set_grid_bounds( "T", glamf, gphif, glamt, gphit ) 184 186 CALL set_grid_bounds( "U", glamv, gphiv, glamu, gphiu ) 185 187 CALL set_grid_bounds( "V", glamu, gphiu, glamv, gphiv ) 186 188 CALL set_grid_bounds( "W", glamf, gphif, glamt, gphit ) 189 CALL set_grid_bounds( "F", glamt, gphit, glamf, gphif ) 187 190 ENDIF 188 191 ENDIF … … 191 194 CALL dom_grid_crs ! Save the parent grid information & Switch to coarse grid domain 192 195 ! 193 CALL set_grid( "T", glamt_crs, gphit_crs, .FALSE., .FALSE. ) 194 CALL set_grid( "U", glamu_crs, gphiu_crs, .FALSE., .FALSE. ) 195 CALL set_grid( "V", glamv_crs, gphiv_crs, .FALSE., .FALSE. ) 196 CALL set_grid( "W", glamt_crs, gphit_crs, .FALSE., .FALSE. ) 196 CALL set_grid( "T", glamt_crs, gphit_crs, .FALSE., .FALSE. ) 197 CALL set_grid( "U", glamu_crs, gphiu_crs, .FALSE., .FALSE. ) 198 CALL set_grid( "V", glamv_crs, gphiv_crs, .FALSE., .FALSE. ) 199 CALL set_grid( "W", glamt_crs, gphit_crs, .FALSE., .FALSE. ) 197 200 CALL set_grid_znl( gphit_crs ) 198 201 ! … … 217 220 CALL iom_set_axis_attr( "depthv", paxis = gdept_1d ) 218 221 CALL iom_set_axis_attr( "depthw", paxis = gdepw_1d ) 222 CALL iom_set_axis_attr( "depthf", paxis = gdept_1d ) 219 223 220 224 ! ABL … … 238 242 CALL iom_set_axis_attr( "depthv", bounds=zw_bnds ) 239 243 CALL iom_set_axis_attr( "depthw", bounds=zt_bnds ) 244 CALL iom_set_axis_attr( "depthf", bounds=zw_bnds ) 240 245 241 246 ! ABL -
NEMO/trunk/src/OCE/IOM/restart.F90
r13970 r14053 11 11 !! 3.7 ! 2014-01 (G. Madec) suppression of curl and hdiv from the restart 12 12 !! - ! 2014-12 (G. Madec) remove KPP scheme 13 !! 4.1 ! 2020-11 (S. Techene, G. Madec) move ssh initiatlisation in DYN/sshwzv:ssh_init_rst 13 14 !!---------------------------------------------------------------------- 14 15 … … 139 140 !! ** Method : Write in numrow when kt == nitrst in NetCDF 140 141 !! file, save fields which are necessary for restart 142 !! 143 !! NB: ssh is written here (rst_write) 144 !! but is read or set in DYN/sshwzv:shh_init_rst 141 145 !!---------------------------------------------------------------------- 142 146 INTEGER, INTENT(in) :: kt ! ocean time-step … … 233 237 !! *** ROUTINE rst_read *** 234 238 !! 235 !! ** Purpose : Read files for NetCDF restart 236 !! 237 !! ** Method : Read in restart.nc file fields which are necessary for restart 239 !! ** Purpose : Read velocity and T-S fields in the restart file 240 !! 241 !! ** Method : Read in restart.nc fields which are necessary for restart 242 !! 243 !! NB: restart file openned in DOM/domain.F90:dom_init 244 !! before field in restart tested in DOM/domain.F90:dom_init 245 !! (sshb) 246 !! 247 !! NB: ssh is read or set in DYN/sshwzv:shh_init_rst 248 !! but is written in IOM/restart:rst_write 238 249 !!---------------------------------------------------------------------- 239 250 INTEGER, INTENT(in) :: Kbb, Kmm ! ocean time level indices 240 REAL(wp) :: zrdt241 251 INTEGER :: jk 242 252 REAL(wp), DIMENSION(jpi, jpj, jpk) :: w3d 243 253 !!---------------------------------------------------------------------- 244 245 CALL rst_read_open ! open restart for reading (if not already opened) 246 247 ! Check dynamics and tracer time-step consistency and force Euler restart if changed 248 IF( iom_varid( numror, 'rdt', ldstop = .FALSE. ) > 0 ) THEN 249 CALL iom_get( numror, 'rdt', zrdt ) 250 IF( zrdt /= rn_Dt ) THEN 251 IF(lwp) WRITE( numout,*) 252 IF(lwp) WRITE( numout,*) 'rst_read: rdt not equal to the read one' 253 IF(lwp) WRITE( numout,*) 254 IF(lwp) WRITE( numout,*) ' ==>>> forced euler first time-step' 255 l_1st_euler = .TRUE. 256 ENDIF 257 ENDIF 258 254 ! 259 255 IF(.NOT.lrxios ) CALL iom_delay_rst( 'READ', 'OCE', numror ) ! read only ocean delayed global communication variables 260 261 ! Diurnal DSST262 IF( ln_diurnal ) CALL iom_get( numror, jpdom_auto, 'Dsst' , x_dsst )256 ! 257 ! !* Diurnal DSST 258 IF( ln_diurnal ) CALL iom_get( numror, jpdom_auto, 'Dsst' , x_dsst ) 263 259 IF ( ln_diurnal_only ) THEN 264 260 IF(lwp) WRITE( numout, * ) & … … 269 265 RETURN 270 266 ENDIF 271 272 IF( iom_varid( numror, 'ub', ldstop = .FALSE. ) > 0 ) THEN 273 ! before fields 267 ! 268 ! !* Read Kmm fields 269 IF(lwp) WRITE(numout,*) ' Kmm u, v and T-S fields read in the restart file' 270 CALL iom_get( numror, jpdom_auto, 'un' , uu(:,:,: ,Kmm), cd_type = 'U', psgn = -1._wp ) 271 CALL iom_get( numror, jpdom_auto, 'vn' , vv(:,:,: ,Kmm), cd_type = 'V', psgn = -1._wp ) 272 CALL iom_get( numror, jpdom_auto, 'tn' , ts(:,:,:,jp_tem,Kmm) ) 273 CALL iom_get( numror, jpdom_auto, 'sn' , ts(:,:,:,jp_sal,Kmm) ) 274 ! 275 IF( l_1st_euler ) THEN !* Euler restart 276 IF(lwp) WRITE(numout,*) ' Kbb u, v and T-S fields set to Kmm values' 277 ts(:,:,:,:,Kbb) = ts(:,:,:,:,Kmm) ! all before fields set to now values 278 uu(:,:,: ,Kbb) = uu(:,:,: ,Kmm) 279 vv(:,:,: ,Kbb) = vv(:,:,: ,Kmm) 280 ELSE !* Leap frog restart 281 IF(lwp) WRITE(numout,*) ' Kbb u, v and T-S fields read in the restart file' 274 282 CALL iom_get( numror, jpdom_auto, 'ub' , uu(:,:,: ,Kbb), cd_type = 'U', psgn = -1._wp ) 275 283 CALL iom_get( numror, jpdom_auto, 'vb' , vv(:,:,: ,Kbb), cd_type = 'V', psgn = -1._wp ) 276 284 CALL iom_get( numror, jpdom_auto, 'tb' , ts(:,:,:,jp_tem,Kbb) ) 277 285 CALL iom_get( numror, jpdom_auto, 'sb' , ts(:,:,:,jp_sal,Kbb) ) 278 CALL iom_get( numror, jpdom_auto, 'sshb' ,ssh(:,: ,Kbb) ) 279 ELSE 280 l_1st_euler = .TRUE. ! before field not found, forced euler 1st time-step 281 ENDIF 282 ! 283 ! now fields 284 CALL iom_get( numror, jpdom_auto, 'un' , uu(:,:,: ,Kmm), cd_type = 'U', psgn = -1._wp ) 285 CALL iom_get( numror, jpdom_auto, 'vn' , vv(:,:,: ,Kmm), cd_type = 'V', psgn = -1._wp ) 286 CALL iom_get( numror, jpdom_auto, 'tn' , ts(:,:,:,jp_tem,Kmm) ) 287 CALL iom_get( numror, jpdom_auto, 'sn' , ts(:,:,:,jp_sal,Kmm) ) 288 CALL iom_get( numror, jpdom_auto, 'sshn' ,ssh(:,: ,Kmm) ) 286 ENDIF 287 ! 289 288 IF( iom_varid( numror, 'rhop', ldstop = .FALSE. ) > 0 ) THEN 290 289 CALL iom_get( numror, jpdom_auto, 'rhop' , rhop ) ! now potential density … … 293 292 ENDIF 294 293 ! 295 IF( l_1st_euler ) THEN ! Euler restart296 ts (:,:,:,:,Kbb) = ts (:,:,:,:,Kmm) ! all before fields set to now values297 uu (:,:,: ,Kbb) = uu (:,:,: ,Kmm)298 vv (:,:,: ,Kbb) = vv (:,:,: ,Kmm)299 ssh (:,: ,Kbb) = ssh (:,: ,Kmm)300 ENDIF301 !302 294 END SUBROUTINE rst_read 303 295 -
NEMO/trunk/src/OCE/ISF/isfcpl.F90
r13970 r14053 10 10 11 11 !!---------------------------------------------------------------------- 12 !! isfrst : read/write iceshelf variables in/from restart12 !! isfrst : read/write iceshelf variables in/from restart 13 13 !!---------------------------------------------------------------------- 14 USE isf_oce ! ice shelf variable 14 USE oce ! ocean dynamics and tracers 15 #if defined key_qco 16 USE domqco , ONLY : dom_qco_zgr ! vertical scale factor interpolation 17 #else 18 USE domvvl , ONLY : dom_vvl_zgr ! vertical scale factor interpolation 19 #endif 20 USE domutl , ONLY : dom_ngb ! find the closest grid point from a given lon/lat position 21 USE isf_oce ! ice shelf variable 15 22 USE isfutils, ONLY : debug 16 USE lib_mpp , ONLY: mpp_sum, mpp_max ! mpp routine17 #if ! defined key_qco18 USE domvvl , ONLY: dom_vvl_zgr ! vertical scale factor interpolation19 #else20 USE domqco , ONLY: dom_qco_zgr ! vertical scale factor interpolation21 #endif22 USE domutl , ONLY: dom_ngb ! find the closest grid point from a given lon/lat position23 23 ! 24 USE oce ! ocean dynamics and tracers25 24 USE in_out_manager ! I/O manager 26 25 USE iom ! I/O library 26 USE lib_mpp , ONLY : mpp_sum, mpp_max ! mpp routine 27 27 ! 28 28 IMPLICIT NONE … … 34 34 35 35 TYPE isfcons 36 INTEGER :: ii ! i global37 INTEGER :: jj ! j global38 INTEGER :: kk ! k level39 REAL(wp):: dvol ! volume increment40 REAL(wp):: dsal ! salt increment41 REAL(wp):: dtem ! heat increment42 REAL(wp):: lon ! lon43 REAL(wp):: lat ! lat44 INTEGER :: ngb ! 0/1 (valid location or not (ie on halo or no neigbourg))36 INTEGER :: ii ! i global 37 INTEGER :: jj ! j global 38 INTEGER :: kk ! k level 39 REAL(wp):: dvol ! volume increment 40 REAL(wp):: dsal ! salt increment 41 REAL(wp):: dtem ! heat increment 42 REAL(wp):: lon ! lon 43 REAL(wp):: lat ! lat 44 INTEGER :: ngb ! 0/1 (valid location or not (ie on halo or no neigbourg)) 45 45 END TYPE 46 46 ! … … 121 121 #endif 122 122 END SUBROUTINE isfcpl_init 123 ! 124 SUBROUTINE isfcpl_rst_write(kt, Kmm) 123 124 125 SUBROUTINE isfcpl_rst_write( kt, Kmm ) 125 126 !!--------------------------------------------------------------------- 126 127 !! *** ROUTINE iscpl_rst_write *** … … 133 134 !!---------------------------------------------------------------------- 134 135 INTEGER :: jk ! loop index 135 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, ze3u, ze3v, zgdepw ! e3t , e3u, e3v !!st patch to usesubstitution136 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, ze3u, ze3v, zgdepw ! for qco substitution 136 137 !!---------------------------------------------------------------------- 137 138 ! … … 153 154 END SUBROUTINE isfcpl_rst_write 154 155 156 155 157 SUBROUTINE isfcpl_ssh(Kbb, Kmm, Kaa) 156 158 !!---------------------------------------------------------------------- … … 184 186 zdssmask(:,:) = ssmask(:,:) - zssmask0(:,:) 185 187 DO_2D( 0, 0, 0, 0 ) 186 jip1=ji+1 ; jim1=ji-1;187 jjp1=jj+1 ; jjm1=jj-1;188 jip1=ji+1 ; jim1=ji-1 189 jjp1=jj+1 ; jjm1=jj-1 188 190 ! 189 191 zsummsk = zssmask0(jip1,jj) + zssmask0(jim1,jj) + zssmask0(ji,jjp1) + zssmask0(ji,jjm1) … … 191 193 IF (zdssmask(ji,jj) == 1._wp .AND. zsummsk /= 0._wp) THEN 192 194 ssh(ji,jj,Kmm)=( zssh(jip1,jj)*zssmask0(jip1,jj) & 193 & + zssh(jim1,jj)*zssmask0(jim1,jj) &194 & + zssh(ji,jjp1)*zssmask0(ji,jjp1) &195 & + zssh(ji,jjm1)*zssmask0(ji,jjm1)) / zsummsk195 & + zssh(jim1,jj)*zssmask0(jim1,jj) & 196 & + zssh(ji,jjp1)*zssmask0(ji,jjp1) & 197 & + zssh(ji,jjm1)*zssmask0(ji,jjm1)) / zsummsk 196 198 zssmask_b(ji,jj) = 1._wp 197 199 ENDIF … … 222 224 CALL dom_vvl_zgr(Kbb, Kmm, Kaa) 223 225 #else 224 CALL dom_qco_zgr(Kbb, Kmm , Kaa)226 CALL dom_qco_zgr(Kbb, Kmm) 225 227 #endif 226 228 ! 227 229 END SUBROUTINE isfcpl_ssh 228 230 231 229 232 SUBROUTINE isfcpl_tra(Kmm) 230 233 !!---------------------------------------------------------------------- … … 375 378 ! 376 379 END SUBROUTINE isfcpl_tra 380 377 381 378 382 SUBROUTINE isfcpl_vol(Kmm) … … 466 470 risfcpl_ssh(:,:) = risfcpl_ssh(:,:) + risfcpl_vol(:,:,jk) * r1_e1e2t(:,:) 467 471 END DO 468 472 ! 469 473 END SUBROUTINE isfcpl_vol 470 474 475 471 476 SUBROUTINE isfcpl_cons(Kmm) 472 477 !!---------------------------------------------------------------------- -
NEMO/trunk/src/OCE/ISF/isfdynatf.F90
r13237 r14053 15 15 USE phycst , ONLY: r1_rho0 ! physical constant 16 16 USE dom_oce ! time and space domain 17 USE oce, ONLY : ssh ! sea-surface height !!st needed forsubstitution17 USE oce, ONLY : ssh ! sea-surface height for qco substitution 18 18 19 19 USE in_out_manager -
NEMO/trunk/src/OCE/ISF/isfrst.F90
r13970 r14053 28 28 !!---------------------------------------------------------------------- 29 29 CONTAINS 30 !31 SUBROUTINE isfrst_read( cdisf, ptsc, pfwf, ptsc_b, pfwf_b )30 31 SUBROUTINE isfrst_read( cdisf, ptsc, pfwf, ptsc_b, pfwf_b ) 32 32 !!--------------------------------------------------------------------- 33 33 !! … … 51 51 ! 52 52 ! read restart 53 IF( iom_varid( numror, cfwf_b, ldstop = .FALSE. ) > 0) THEN53 IF( .NOT.l_1st_euler ) THEN 54 54 IF(lwp) WRITE(numout,*) ' nit000-1 isf tracer content forcing fields read in the restart file' 55 55 CALL iom_get( numror, jpdom_auto, cfwf_b, pfwf_b(:,:) ) ! before ice shelf melt … … 62 62 ! 63 63 END SUBROUTINE isfrst_read 64 ! 65 SUBROUTINE isfrst_write(kt, cdisf, ptsc, pfwf ) 64 65 66 SUBROUTINE isfrst_write( kt, cdisf, ptsc, pfwf ) 66 67 !!--------------------------------------------------------------------- 67 68 !! … … 94 95 ! 95 96 END SUBROUTINE isfrst_write 96 ! 97 98 !!====================================================================== 97 99 END MODULE isfrst -
NEMO/trunk/src/OCE/LBC/mppini.F90
r13982 r14053 217 217 ! then we calculate them here now that we have our communicator size 218 218 IF(lwp) THEN 219 WRITE(numout,*) 219 220 WRITE(numout,*) 'mpp_init:' 220 221 WRITE(numout,*) '~~~~~~~~ ' 221 WRITE(numout,*)222 222 ENDIF 223 223 IF( jpni < 1 .OR. jpnj < 1 ) THEN -
NEMO/trunk/src/OCE/LDF/ldfdyn.F90
r13497 r14053 34 34 ! !!* Namelist namdyn_ldf : lateral mixing on momentum * 35 35 LOGICAL , PUBLIC :: ln_dynldf_OFF !: No operator (i.e. no explicit diffusion) 36 INTEGER , PUBLIC :: nn_dynldf_typ !: operator type (0: div-rot ; 1: symmetric) 36 37 LOGICAL , PUBLIC :: ln_dynldf_lap !: laplacian operator 37 38 LOGICAL , PUBLIC :: ln_dynldf_blp !: bilaplacian operator … … 52 53 53 54 ! !!* Parameter to control the type of lateral viscous operator 54 INTEGER, PARAMETER, PUBLIC :: np_ERROR =-10 !: error in setting the operator 55 INTEGER, PARAMETER, PUBLIC :: np_no_ldf = 00 !: without operator (i.e. no lateral viscous trend) 55 INTEGER, PARAMETER, PUBLIC :: np_ERROR =-10 !: error in setting the operator 56 INTEGER, PARAMETER, PUBLIC :: np_no_ldf = 00 !: without operator (i.e. no lateral viscous trend) 57 ! 58 INTEGER, PARAMETER, PUBLIC :: np_typ_rot = 0 !: div-rot operator 59 INTEGER, PARAMETER, PUBLIC :: np_typ_sym = 1 !: symmetric operator 60 ! 56 61 ! !! laplacian ! bilaplacian ! 57 62 INTEGER, PARAMETER, PUBLIC :: np_lap = 10 , np_blp = 20 !: iso-level operator … … 109 114 CHARACTER(len=5) :: cl_Units ! units (m2/s or m4/s) 110 115 !! 111 NAMELIST/namdyn_ldf/ ln_dynldf_OFF, ln_dynldf_lap, ln_dynldf_blp, & ! type of operator112 & ln_dynldf_lev, ln_dynldf_hor, ln_dynldf_iso, & ! acting direction of the operator113 & nn_ahm_ijk_t , rn_Uv , rn_Lv, rn_ahm_b,& ! lateral eddy coefficient114 & rn_csmc , rn_minfac , rn_maxfac ! Smagorinsky settings116 NAMELIST/namdyn_ldf/ ln_dynldf_OFF, nn_dynldf_typ, ln_dynldf_lap, ln_dynldf_blp, & ! type of operator 117 & ln_dynldf_lev, ln_dynldf_hor, ln_dynldf_iso, & ! acting direction of the operator 118 & nn_ahm_ijk_t , rn_Uv , rn_Lv , rn_ahm_b, & ! lateral eddy coefficient 119 & rn_csmc , rn_minfac , rn_maxfac ! Smagorinsky settings 115 120 !!---------------------------------------------------------------------- 116 121 ! … … 130 135 WRITE(numout,*) ' type :' 131 136 WRITE(numout,*) ' no explicit diffusion ln_dynldf_OFF = ', ln_dynldf_OFF 137 WRITE(numout,*) ' type of operator (div-rot or sym) nn_dynldf_typ = ', nn_dynldf_typ 132 138 WRITE(numout,*) ' laplacian operator ln_dynldf_lap = ', ln_dynldf_lap 133 139 WRITE(numout,*) ' bilaplacian operator ln_dynldf_blp = ', ln_dynldf_blp … … 147 153 WRITE(numout,*) ' Smagorinsky coefficient rn_csmc = ', rn_csmc 148 154 WRITE(numout,*) ' factor multiplier for eddy visc.' 149 WRITE(numout,*) ' lower limit (default 1.0) rn_minfac = ', rn_minfac150 WRITE(numout,*) ' upper limit (default 1.0) rn_maxfac = ', rn_maxfac155 WRITE(numout,*) ' lower limit (default 1.0) rn_minfac = ', rn_minfac 156 WRITE(numout,*) ' upper limit (default 1.0) rn_maxfac = ', rn_maxfac 151 157 ENDIF 152 158 … … 160 166 IF( ln_dynldf_lap ) THEN ; ioptio = ioptio + 1 ; ENDIF 161 167 IF( ln_dynldf_blp ) THEN ; ioptio = ioptio + 1 ; ENDIF 162 IF( ioptio /= 1 ) CALL ctl_stop( ' dyn_ldf_init: use ONE of the 3 operator options (NONE/lap/blp)' )168 IF( ioptio /= 1 ) CALL ctl_stop( 'ldf_dyn_init: use ONE of the 3 operator options (NONE/lap/blp)' ) 163 169 ! 164 170 IF(.NOT.ln_dynldf_OFF ) THEN !== direction ==>> type of operator ==! 171 ! 172 SELECT CASE( nn_dynldf_typ ) ! div-rot or symmetric 173 CASE( np_typ_rot ) ; IF(lwp) WRITE(numout,*) ' ==>>> use div-rot operator ' 174 CASE( np_typ_sym ) ; IF(lwp) WRITE(numout,*) ' ==>>> use symmetric operator ' 175 CASE DEFAULT ! error 176 CALL ctl_stop('ldf_dyn_init: wrong value for nn_dynldf_typ (0 or 1)' ) 177 END SELECT 178 ! 165 179 ioptio = 0 166 180 IF( ln_dynldf_lev ) ioptio = ioptio + 1 167 181 IF( ln_dynldf_hor ) ioptio = ioptio + 1 168 182 IF( ln_dynldf_iso ) ioptio = ioptio + 1 169 IF( ioptio /= 1 ) CALL ctl_stop( ' dyn_ldf_init: use ONE of the 3 direction options (level/hor/iso)' )183 IF( ioptio /= 1 ) CALL ctl_stop( 'ldf_dyn_init: use ONE of the 3 direction options (level/hor/iso)' ) 170 184 ! 171 185 ! ! Set nldf_dyn, the type of lateral diffusion, from ln_dynldf_... logicals -
NEMO/trunk/src/OCE/SBC/sbcapr.F90
r13970 r14053 148 148 ! ! ---------------------------------------- ! 149 149 ! !* Restart: read in restart file 150 IF( ln_rstart .AND. iom_varid( numror, 'ssh_ibb', ldstop = .FALSE. ) > 0) THEN150 IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN 151 151 IF(lwp) WRITE(numout,*) 'sbc_apr: ssh_ibb read in the restart file' 152 152 CALL iom_get( numror, jpdom_auto, 'ssh_ibb', ssh_ibb ) ! before inv. barometer ssh -
NEMO/trunk/src/OCE/SBC/sbcice_cice.F90
r13295 r14053 12 12 USE oce ! ocean dynamics and tracers 13 13 USE dom_oce ! ocean space and time domain 14 # if !defined key_qco15 USE dom vvl14 # if defined key_qco 15 USE domqco ! Variable volume 16 16 # else 17 USE dom qco17 USE domvvl ! Variable volume 18 18 # endif 19 19 USE phycst, only : rcp, rho0, r1_rho0, rhos, rhoi … … 238 238 !!gm especially here it is assumed zstar coordinate, but it can be ztilde.... 239 239 #if defined key_qco 240 IF( .NOT.ln_linssh ) CALL dom_qco_zgr( Kbb, Kmm , Kaa) ! interpolation scale factor, depth and water column240 IF( .NOT.ln_linssh ) CALL dom_qco_zgr( Kbb, Kmm ) ! interpolation scale factor, depth and water column 241 241 #else 242 242 IF( .NOT.ln_linssh ) THEN -
NEMO/trunk/src/OCE/SBC/sbcmod.F90
r14030 r14053 523 523 IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! 524 524 ! ! ---------------------------------------- ! 525 IF( ln_rstart .AND. & !* Restart: read in restart file 526 & iom_varid( numror, 'utau_b', ldstop = .FALSE. ) > 0 ) THEN 527 IF(lwp) WRITE(numout,*) ' nit000-1 surface forcing fields red in the restart file' 528 CALL iom_get( numror, jpdom_auto, 'utau_b', utau_b ) ! before i-stress (U-point) 529 CALL iom_get( numror, jpdom_auto, 'vtau_b', vtau_b ) ! before j-stress (V-point) 530 CALL iom_get( numror, jpdom_auto, 'qns_b', qns_b ) ! before non solar heat flux (T-point) 531 ! The 3D heat content due to qsr forcing is treated in traqsr 532 ! CALL iom_get( numror, jpdom_auto, 'qsr_b' , qsr_b ) ! before solar heat flux (T-point) 533 CALL iom_get( numror, jpdom_auto, 'emp_b', emp_b ) ! before freshwater flux (T-point) 525 IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN !* Restart: read in restart file 526 IF(lwp) WRITE(numout,*) ' nit000-1 surface forcing fields read in the restart file' 527 CALL iom_get( numror, jpdom_auto, 'utau_b', utau_b ) ! i-stress 528 CALL iom_get( numror, jpdom_auto, 'vtau_b', vtau_b ) ! j-stress 529 CALL iom_get( numror, jpdom_auto, 'qns_b', qns_b ) ! non solar heat flux 530 CALL iom_get( numror, jpdom_auto, 'emp_b', emp_b ) ! freshwater flux 531 ! NB: The 3D heat content due to qsr forcing (qsr_hc_b) is treated in traqsr 534 532 ! To ensure restart capability with 3.3x/3.4 restart files !! to be removed in v3.6 535 533 IF( iom_varid( numror, 'sfx_b', ldstop = .FALSE. ) > 0 ) THEN … … 566 564 ! ! ---------------------------------------- ! 567 565 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 568 CALL iom_put( "empmr" , emp 569 CALL iom_put( "empbmr" , emp_b 570 CALL iom_put( "saltflx", sfx )! downward salt flux (includes virtual salt flux beneath ice in linear free surface case)571 CALL iom_put( "fmmflx" , fmmflx )! Freezing-melting water flux572 CALL iom_put( "qt" , qns + qsr )! total heat flux573 CALL iom_put( "qns" , qns )! solar heat flux574 CALL iom_put( "qsr" , qsr )! solar heat flux566 CALL iom_put( "empmr" , emp - rnf ) ! upward water flux 567 CALL iom_put( "empbmr" , emp_b - rnf ) ! before upward water flux ( needed to recalculate the time evolution of ssh in offline ) 568 CALL iom_put( "saltflx", sfx ) ! downward salt flux (includes virtual salt flux beneath ice in linear free surface case) 569 CALL iom_put( "fmmflx" , fmmflx ) ! Freezing-melting water flux 570 CALL iom_put( "qt" , qns + qsr ) ! total heat flux 571 CALL iom_put( "qns" , qns ) ! solar heat flux 572 CALL iom_put( "qsr" , qsr ) ! solar heat flux 575 573 IF( nn_ice > 0 .OR. ll_opa ) CALL iom_put( "ice_cover", fr_i ) ! ice fraction 576 CALL iom_put( "taum" , taum )! wind stress module577 CALL iom_put( "wspd" , wndm )! wind speed module over free ocean or leads in presence of sea-ice578 CALL iom_put( "qrp" , qrp )! heat flux damping579 CALL iom_put( "erp" , erp )! freshwater flux damping574 CALL iom_put( "taum" , taum ) ! wind stress module 575 CALL iom_put( "wspd" , wndm ) ! wind speed module over free ocean or leads in presence of sea-ice 576 CALL iom_put( "qrp" , qrp ) ! heat flux damping 577 CALL iom_put( "erp" , erp ) ! freshwater flux damping 580 578 ENDIF 581 579 ! 582 580 IF(sn_cfctl%l_prtctl) THEN ! print mean trends (used for debugging) 583 CALL prt_ctl(tab2d_1=fr_i , clinfo1=' fr_i - : ', mask1=tmask )584 CALL prt_ctl(tab2d_1=(emp-rnf) , clinfo1=' emp-rnf - : ', mask1=tmask )585 CALL prt_ctl(tab2d_1=(sfx-rnf) , clinfo1=' sfx-rnf - : ', mask1=tmask )586 CALL prt_ctl(tab2d_1=qns , clinfo1=' qns - : ', mask1=tmask )587 CALL prt_ctl(tab2d_1=qsr , clinfo1=' qsr - : ', mask1=tmask )588 CALL prt_ctl(tab3d_1=tmask , clinfo1=' tmask - : ', mask1=tmask, kdim=jpk )581 CALL prt_ctl(tab2d_1=fr_i , clinfo1=' fr_i - : ', mask1=tmask ) 582 CALL prt_ctl(tab2d_1=(emp-rnf) , clinfo1=' emp-rnf - : ', mask1=tmask ) 583 CALL prt_ctl(tab2d_1=(sfx-rnf) , clinfo1=' sfx-rnf - : ', mask1=tmask ) 584 CALL prt_ctl(tab2d_1=qns , clinfo1=' qns - : ', mask1=tmask ) 585 CALL prt_ctl(tab2d_1=qsr , clinfo1=' qsr - : ', mask1=tmask ) 586 CALL prt_ctl(tab3d_1=tmask , clinfo1=' tmask - : ', mask1=tmask, kdim=jpk ) 589 587 CALL prt_ctl(tab3d_1=ts(:,:,:,jp_tem,Kmm), clinfo1=' sst - : ', mask1=tmask, kdim=1 ) 590 588 CALL prt_ctl(tab3d_1=ts(:,:,:,jp_sal,Kmm), clinfo1=' sss - : ', mask1=tmask, kdim=1 ) -
NEMO/trunk/src/OCE/SBC/sbcrnf.F90
r14032 r14053 157 157 IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! 158 158 ! ! ---------------------------------------- ! 159 IF( ln_rstart .AND. & !* Restart: read in restart file 160 & iom_varid( numror, 'rnf_b', ldstop = .FALSE. ) > 0 ) THEN 159 IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN !* Restart: read in restart file 161 160 IF(lwp) WRITE(numout,*) ' nit000-1 runoff forcing fields red in the restart file', lrxios 162 CALL iom_get( numror, jpdom_auto, 'rnf_b' , rnf_b )! before runoff161 CALL iom_get( numror, jpdom_auto, 'rnf_b' , rnf_b ) ! before runoff 163 162 CALL iom_get( numror, jpdom_auto, 'rnf_hc_b', rnf_tsc_b(:,:,jp_tem) ) ! before heat content of runoff 164 163 CALL iom_get( numror, jpdom_auto, 'rnf_sc_b', rnf_tsc_b(:,:,jp_sal) ) ! before salinity content of runoff 165 ELSE 164 ELSE !* no restart: set from nit000 values 166 165 IF(lwp) WRITE(numout,*) ' nit000-1 runoff forcing fields set to nit000' 167 166 rnf_b (:,: ) = rnf (:,: ) … … 176 175 & 'at it= ', kt,' date= ', ndastp 177 176 IF(lwp) WRITE(numout,*) '~~~~' 178 CALL iom_rstput( kt, nitrst, numrow, 'rnf_b' , rnf)177 CALL iom_rstput( kt, nitrst, numrow, 'rnf_b' , rnf ) 179 178 CALL iom_rstput( kt, nitrst, numrow, 'rnf_hc_b', rnf_tsc(:,:,jp_tem) ) 180 179 CALL iom_rstput( kt, nitrst, numrow, 'rnf_sc_b', rnf_tsc(:,:,jp_sal) ) -
NEMO/trunk/src/OCE/TRA/traatf_qco.F90
r13982 r14053 1 MODULE traatf qco1 MODULE traatf_qco 2 2 !!====================================================================== 3 !! *** MODULE traatf qco ***3 !! *** MODULE traatf_qco *** 4 4 !! Ocean active tracers: Asselin time filtering for temperature and salinity 5 5 !!====================================================================== … … 45 45 USE prtctl ! Print control 46 46 USE timing ! Timing 47 #if defined key_agrif48 USE agrif_oce_interp49 #endif50 47 51 48 IMPLICIT NONE … … 149 146 ENDIF 150 147 ! 151 CALL lbc_lnk_multi( 'traatfqco', pts(:,:,:,jp_tem,Kmm) , 'T', 1. , pts(:,:,:,jp_sal,Kmm) , 'T', 1.)152 148 CALL lbc_lnk_multi( 'traatfqco', pts(:,:,:,jp_tem,Kmm) , 'T', 1._wp, pts(:,:,:,jp_sal,Kmm) , 'T', 1._wp ) 149 ! 153 150 ENDIF 154 151 ! … … 370 367 371 368 !!====================================================================== 372 END MODULE traatf qco369 END MODULE traatf_qco -
NEMO/trunk/src/OCE/TRA/traqsr.F90
r13982 r14053 144 144 145 145 IF( kt == nit000 ) THEN !== 1st time step ==! 146 IF( ln_rstart .AND. iom_varid( numror, 'qsr_hc_b', ldstop = .FALSE. ) > 0 .AND..NOT.l_1st_euler ) THEN ! read in restart146 IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN ! read in restart 147 147 z1_2 = 0.5_wp 148 148 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile … … 150 150 CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b ) ! before heat content trend due to Qsr flux 151 151 ENDIF 152 ELSE ! No restart or restart not found: Euler forward time stepping152 ELSE ! No restart or Euler forward at 1st time step 153 153 z1_2 = 1._wp 154 154 DO_3D( isj, iej, isi, iei, 1, jpk ) -
NEMO/trunk/src/OCE/TRA/trasbc.F90
r13982 r14053 72 72 !! - send trends to trdtra module for further diagnostics(l_trdtra=T) 73 73 !!---------------------------------------------------------------------- 74 INTEGER, INTENT(in ) :: kt ! ocean time-step index75 INTEGER, INTENT(in ) :: Kmm, Krhs ! time level indices76 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation74 INTEGER, INTENT(in ) :: kt ! ocean time-step index 75 INTEGER, INTENT(in ) :: Kmm, Krhs ! time level indices 76 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer Eq. 77 77 ! 78 78 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 117 117 ! !== Set before sbc tracer content fields ==! 118 118 IF( kt == nit000 ) THEN !* 1st time-step 119 IF( ln_rstart .AND. & ! Restart: read in restart file 120 & iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN 119 IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN ! Restart: read in restart file 121 120 zfact = 0.5_wp 122 121 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile … … 126 125 CALL iom_get( numror, jpdom_auto, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal) ) ! before salt content sbc trend 127 126 ENDIF 128 ELSE ! No restart or restart not found: Euler forward time stepping127 ELSE ! No restart or restart not found: Euler forward time stepping 129 128 zfact = 1._wp 130 129 DO_2D( isj, iej, isi, iei ) -
NEMO/trunk/src/OCE/USR/usrdef_istate.F90
r13497 r14053 7 7 !! User defined : set the initial state of a user configuration 8 8 !!====================================================================== 9 !! History : 4.0 ! 2016-03 (S. Flavoni) Original code 9 !! History : 4.0 ! 2016-03 (S. Flavoni) Original code 10 !! ! 2020-11 (S. Techene, G. Madec) separate tsuv from ssh 10 11 !!---------------------------------------------------------------------- 11 12 … … 22 23 PRIVATE 23 24 24 PUBLIC usr_def_istate ! called in istate.F90 25 PUBLIC usr_def_istate ! called in istate.F90 26 PUBLIC usr_def_istate_ssh ! called by domqco.F90 25 27 26 28 !! * Substitutions … … 33 35 CONTAINS 34 36 35 SUBROUTINE usr_def_istate( pdept, ptmask, pts, pu, pv , pssh)37 SUBROUTINE usr_def_istate( pdept, ptmask, pts, pu, pv ) 36 38 !!---------------------------------------------------------------------- 37 39 !! *** ROUTINE usr_def_istate *** … … 48 50 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT( out) :: pu ! i-component of the velocity [m/s] 49 51 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT( out) :: pv ! j-component of the velocity [m/s] 50 REAL(wp), DIMENSION(jpi,jpj) , INTENT( out) :: pssh ! sea-surface height51 52 ! 52 53 INTEGER :: ji, jj, jk ! dummy loop indices … … 59 60 pu (:,:,:) = 0._wp ! ocean at rest 60 61 pv (:,:,:) = 0._wp 61 pssh(:,:) = 0._wp62 62 ! 63 63 DO_3D( 1, 1, 1, 1, 1, jpk ) ! horizontally uniform T & S profiles … … 80 80 END SUBROUTINE usr_def_istate 81 81 82 83 SUBROUTINE usr_def_istate_ssh( ptmask, pssh ) 84 !!---------------------------------------------------------------------- 85 !! *** ROUTINE usr_def_istate_ssh *** 86 !! 87 !! ** Purpose : Initialization of ssh 88 !! 89 !! ** Method : Set ssh as null, ptmask is required for test cases 90 !!---------------------------------------------------------------------- 91 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: ptmask ! t-point ocean mask [m] 92 REAL(wp), DIMENSION(jpi,jpj) , INTENT( out) :: pssh ! sea-surface height [m] 93 !!---------------------------------------------------------------------- 94 ! 95 IF(lwp) WRITE(numout,*) 96 IF(lwp) WRITE(numout,*) 'usr_def_istate_ssh : GYRE configuration, analytical definition of initial state' 97 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~~ Ocean at rest, ssh is zero' 98 ! 99 ! Sea level: 100 pssh(:,:) = 0._wp 101 ! 102 END SUBROUTINE usr_def_istate_ssh 103 82 104 !!====================================================================== 83 105 END MODULE usrdef_istate -
NEMO/trunk/src/OCE/ZDF/zdfddm.F90
r13497 r14053 31 31 !! * Substitutions 32 32 # include "do_loop_substitute.h90" 33 # include "domzgr_substitute.h90" 33 34 !!---------------------------------------------------------------------- 34 35 !! NEMO/OCE 4.0 , NEMO Consortium (2018) -
NEMO/trunk/src/OCE/nemogcm.F90
r14032 r14053 42 42 !!---------------------------------------------------------------------- 43 43 USE step_oce ! module used in the ocean time stepping module (step.F90) 44 ! 44 45 USE phycst ! physical constant (par_cst routine) 45 46 USE domain ! domain initialization (dom_init & dom_cfg routines) 46 USE closea ! treatment of closed seas (for ln_closea)47 USE usrdef_nam ! user defined configuration 48 USE tide_mod, ONLY : tide_init ! tidal components initialization (tide_init routine)49 USE bdyini 47 USE wet_dry ! Wetting and drying setting (wad_init routine) 48 USE usrdef_nam ! user defined configuration namelist 49 USE tide_mod, ONLY : tide_init ! tidal components initialization (tide_init routine) 50 USE bdyini , ONLY : bdy_init ! open boundary cond. setting (bdy_init routine) 50 51 USE istate ! initial state setting (istate_init routine) 51 USE ldfdyn ! lateral viscosity setting (ldfdyn_init routine)52 USE ldftra ! lateral diffusivity setting (ldftra_init routine)53 52 USE trdini ! dyn/tra trends initialization (trd_init routine) 54 USE asminc ! assimilation increments 55 USE asmbkg ! writing out state trajectory 56 USE diadct ! sections transports (dia_dct_init routine) 57 USE diaobs ! Observation diagnostics (dia_obs_init routine) 58 USE diacfl ! CFL diagnostics (dia_cfl_init routine) 59 USE diamlr ! IOM context management for multiple-linear-regression analysis 53 USE icbini ! handle bergs, initialisation 54 USE icbstp , ONLY : icb_end ! handle bergs, close iceberg files 55 USE cpl_oasis3 ! OASIS3 coupling 56 USE dyndmp ! Momentum damping (C1D only) 57 USE step_diu ! diurnal bulk SST timestepping (called from here if run offline) 58 USE crsini ! initialise grid coarsening utility 59 USE dia25h , ONLY : dia_25h_init ! 25h mean output (initialisation) 60 USE c1d ! 1D configuration 61 USE step_c1d ! Time stepping loop for the 1D configuration 62 #if defined key_top 63 USE trcini ! passive tracer initialisation 64 #endif 65 #if defined key_nemocice_decomp 66 USE ice_domain_size, only: nx_global, ny_global 67 #endif 60 68 #if defined key_qco 61 USE st epMLF! NEMO time-stepping (stp_MLF routine)69 USE stpmlf ! NEMO time-stepping (stp_MLF routine) 62 70 #else 63 71 USE step ! NEMO time-stepping (stp routine) 64 72 #endif 65 USE isfstp ! ice shelf (isf_stp_init routine)66 USE icbini ! handle bergs, initialisation67 USE icbstp ! handle bergs, calving, themodynamics and transport68 USE cpl_oasis3 ! OASIS3 coupling69 USE c1d ! 1D configuration70 USE step_c1d ! Time stepping loop for the 1D configuration71 USE dyndmp ! Momentum damping72 USE stopar ! Stochastic param.: ???73 USE stopts ! Stochastic param.: ???74 USE diu_layers ! diurnal bulk SST and coolskin75 USE step_diu ! diurnal bulk SST timestepping (called from here if run offline)76 USE crsini ! initialise grid coarsening utility77 USE dia25h ! 25h mean output78 USE diadetide ! Weights computation for daily detiding of model diagnostics79 USE sbc_oce , ONLY : lk_oasis80 USE wet_dry ! Wetting and drying setting (wad_init routine)81 #if defined key_top82 USE trcini ! passive tracer initialisation83 #endif84 #if defined key_nemocice_decomp85 USE ice_domain_size, only: nx_global, ny_global86 #endif87 73 ! 88 USE prtctl ! Print control89 USE in_out_manager ! I/O manager90 74 USE lib_mpp ! distributed memory computing 91 75 USE mppini ! shared/distributed memory setting (mpp_init routine) 92 76 USE lbcnfd , ONLY : isendto, nsndto ! Setup of north fold exchanges 93 77 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 94 #if defined key_iomput 95 USE xios ! xIOserver 96 #endif 97 #if defined key_agrif 98 USE agrif_all_update ! Master Agrif update 99 #endif 100 USE halo_mng 78 USE halo_mng ! Halo manager 101 79 102 80 IMPLICIT NONE … … 182 160 ! 183 161 DO WHILE( istp <= nitend .AND. nstop == 0 ) 184 # if defined key_qco162 # if defined key_qco 185 163 CALL stp_MLF 186 # else164 # else 187 165 CALL stp 188 # endif166 # endif 189 167 istp = istp + 1 190 168 END DO … … 195 173 ! 196 174 DO WHILE( istp <= nitend .AND. nstop == 0 ) 197 175 ! 198 176 ncom_stp = istp 199 177 IF( ln_timing ) THEN … … 202 180 IF ( istp == nitend ) elapsed_time = zstptiming - elapsed_time 203 181 ENDIF 204 205 # if defined key_qco182 ! 183 # if defined key_qco 206 184 CALL stp_MLF ( istp ) 207 # else185 # else 208 186 CALL stp ( istp ) 209 # endif187 # endif 210 188 istp = istp + 1 211 189 ! 212 190 IF( lwp .AND. ln_timing ) WRITE(numtime,*) 'timing step ', istp-1, ' : ', MPI_Wtime() - zstptiming 213 191 ! 214 192 END DO 215 193 ! … … 279 257 INTEGER :: ios, ilocal_comm ! local integers 280 258 !! 281 NAMELIST/namctl/ sn_cfctl, ln_timing, ln_diacfl, 282 & nn_isplt, nn_jsplt, nn_ictls, nn_ictle, nn_jctls, nn_jctle259 NAMELIST/namctl/ sn_cfctl, ln_timing, ln_diacfl, nn_isplt, nn_jsplt , nn_ictls, & 260 & nn_ictle, nn_jctls , nn_jctle 283 261 NAMELIST/namcfg/ ln_read_cfg, cn_domcfg, ln_closea, ln_write_cfg, cn_domcfg_out, ln_use_jattr 284 262 !!---------------------------------------------------------------------- … … 350 328 IF(lwp) THEN ! open listing units 351 329 ! 352 IF( .NOT. 330 IF( .NOT.lwm ) & ! alreay opened for narea == 1 353 331 & CALL ctl_opn( numout, 'ocean.output', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, -1, .FALSE., narea ) 354 332 ! … … 357 335 WRITE(numout,*) ' NEMO team' 358 336 WRITE(numout,*) ' Ocean General Circulation Model' 359 WRITE(numout,*) ' NEMO version 4.0 (20 19) '337 WRITE(numout,*) ' NEMO version 4.0 (2020) ' 360 338 WRITE(numout,*) 361 339 WRITE(numout,*) " ._ ._ ._ ._ ._ " … … 373 351 WRITE(numout,*) " ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ " 374 352 WRITE(numout,*) 375 376 ! Print the working precision to ocean.output 377 IF (wp == dp) THEN 378 WRITE(numout,*) "Working precision = double-precision" 379 ELSE 380 WRITE(numout,*) "Working precision = single-precision" 353 ! 354 WRITE(numout,cform_aaa) ! Flag AAAAAAA 355 ! 356 ! ! Control print of the working precision 357 WRITE(numout,*) 358 IF( wp == dp ) THEN ; WRITE(numout,*) "par_kind : wp = Working precision = dp = double-precision" 359 ELSE ; WRITE(numout,*) "par_kind : wp = Working precision = sp = single-precision" 381 360 ENDIF 382 WRITE(numout,*) 383 ! 384 WRITE(numout,cform_aaa) ! Flag AAAAAAA 361 WRITE(numout,*) "~~~~~~~~ ****************" 362 WRITE(numout,*) 385 363 ! 386 364 ENDIF … … 415 393 416 394 ! Initialise time level indices 417 Nbb = 1 ; Nnn = 2; Naa = 3;Nrhs = Naa395 Nbb = 1 ; Nnn = 2 ; Naa = 3 ; Nrhs = Naa 418 396 #if defined key_agrif 419 Kbb_a = Nbb ; Kmm_a = Nnn;Krhs_a = Nrhs ! agrif_oce module copies of time level indices397 Kbb_a = Nbb ; Kmm_a = Nnn ; Krhs_a = Nrhs ! agrif_oce module copies of time level indices 420 398 #endif 421 399 ! !-------------------------------! … … 423 401 ! !-------------------------------! 424 402 425 CALL nemo_ctl ! Control prints 403 CALL nemo_ctl ! Control prints of namctl and namcfg 426 404 ! 427 405 ! ! General initialization … … 437 415 CALL Agrif_Declare_Var_ini ! " " " " " DOM 438 416 #endif 439 CALL dom_init( Nbb, Nnn, Naa ) ! Domain440 IF( ln_crs ) CALL crs_init( Nnn )! coarsened grid: domain initialization417 CALL dom_init( Nbb, Nnn, Naa ) ! Domain 418 IF( ln_crs ) CALL crs_init( Nnn ) ! coarsened grid: domain initialization 441 419 IF( sn_cfctl%l_prtctl ) & 442 420 & CALL prt_ctl_init ! Print control -
NEMO/trunk/src/OCE/oce.F90
r13237 r14053 16 16 PRIVATE 17 17 18 PUBLIC oce_alloc ! routine called by nemo_init in nemogcm.F90 18 PUBLIC oce_alloc ! routine called by nemo_init in nemogcm.F90 19 PUBLIC oce_SWE_alloc ! routine called by nemo_init in SWE/nemogcm.F90 (Shallow Water Eq. case) 19 20 20 21 !! dynamics and tracer fields … … 68 69 INTEGER, PUBLIC, DIMENSION(2) :: noce_array !: unused array but seems to be needed to prevent agrif from creating an empty module 69 70 71 !! Shallow Water Eq. case (SWE) 72 LOGICAL, PUBLIC :: lk_SWE = .FALSE. !: shallow water flag =T in SWE configurations only 73 74 !! Stand-Alone Surface module (SAS) 75 LOGICAL, PUBLIC :: l_SAS = .FALSE. !: SAS flag =T in SAS configurations only 76 77 70 78 !!---------------------------------------------------------------------- 71 79 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 115 123 END FUNCTION oce_alloc 116 124 125 126 INTEGER FUNCTION oce_SWE_alloc() 127 !!---------------------------------------------------------------------- 128 !! *** FUNCTION oce_SWE_alloc *** 129 !!---------------------------------------------------------------------- 130 INTEGER :: ierr(2) 131 !!---------------------------------------------------------------------- 132 ! 133 lk_SWE = .TRUE. ! =T SWE case 134 ! 135 ierr(:) = 0 136 ALLOCATE( uu(jpi,jpj,jpk,jpt) , vv (jpi,jpj,jpk,jpt) , & 137 & ww(jpi,jpj,jpk) , hdiv(jpi,jpj,jpk) , ssh(jpi,jpj,jpt) , STAT=ierr(1) ) 138 ! 139 ALLOCATE( ts(jpi,jpj,jpk,jpts,jpt) , fraqsr_1lev(jpi,jpj) , & 140 & uu_b(jpi,jpj,jpt) , vv_b(jpi,jpj,jpt) , rn2(jpi,jpj,jpk) , STAT=ierr(2) ) 141 ! 142 oce_SWE_alloc = MAXVAL( ierr ) 143 IF( oce_SWE_alloc /= 0 ) CALL ctl_stop( 'STOP', 'oce_SWE_alloc: failed to allocate arrays' ) 144 ! 145 END FUNCTION oce_SWE_alloc 146 117 147 !!====================================================================== 118 148 END MODULE oce -
NEMO/trunk/src/OCE/step.F90
r14010 r14053 42 42 !!---------------------------------------------------------------------- 43 43 USE step_oce ! time stepping definition modules 44 !45 USE iom ! xIOs server46 44 47 45 IMPLICIT NONE -
NEMO/trunk/src/OCE/step_oce.F90
r14010 r14053 3 3 !! *** MODULE step_oce *** 4 4 !! Ocean time-stepping : module used in both initialisation phase and time stepping 5 !! (i.e. nemo_init and stp or stp_MLF routines) 5 6 !!====================================================================== 6 7 !! History : 3.3 ! 2010-08 (C. Ethe) Original code - reorganisation of the initial phase … … 9 10 USE oce ! ocean dynamics and tracers variables 10 11 USE dom_oce ! ocean space and time domain variables 11 USE domain, ONLY : dom_tile 12 USE zdf_oce ! ocean vertical physics variables 13 USE zdfdrg , ONLY : ln_drgimp ! implicit top/bottom friction 12 USE domain , ONLY : dom_tile 14 13 15 14 USE daymod ! calendar (day routine) … … 20 19 USE sbccpl ! surface boundary condition: coupled formulation (call send at end of step) 21 20 USE sbcapr ! surface boundary condition: atmospheric pressure 22 USE tide_mod, ONLY : ln_tide, tide_update23 21 USE sbcwave ! Wave intialisation 22 USE tide_mod ! tides 23 24 USE bdy_oce , ONLY : ln_bdy 25 USE bdydta ! open boundary condition data (bdy_dta routine) 26 USE bdytra ! bdy cond. for tracers (bdy_tra routine) 27 USE bdydyn3d ! bdy cond. for baroclinic vel. (bdy_dyn3d routine) 24 28 25 29 USE isf_oce ! ice shelf boundary condition 26 30 USE isfstp ! ice shelf boundary condition (isf_stp routine) 31 32 USE sshwzv ! vertical velocity and ssh (ssh_nxt routine) 33 ! (ssh_swp routine) 34 ! (wzv routine) 35 USE domvvl ! variable vertical scale factors (dom_vvl_sf_nxt routine) 36 ! (dom_vvl_sf_swp routine) 37 38 USE divhor ! horizontal divergence (div_hor routine) 39 USE dynadv ! advection (dyn_adv routine) 40 USE dynvor ! vorticity term (dyn_vor routine) 41 USE dynhpg ! hydrostatic pressure grad. (dyn_hpg routine) 42 USE dynldf ! lateral momentum diffusion (dyn_ldf routine) 43 USE dynzdf ! vertical diffusion (dyn_zdf routine) 44 USE dynspg ! surface pressure gradient (dyn_spg routine) 45 USE dynatf ! time-filtering (dyn_atf routine) 27 46 28 47 USE traqsr ! solar radiation penetration (tra_qsr routine) … … 40 59 USE eosbn2 ! equation of state (eos_bn2 routine) 41 60 42 USE divhor ! horizontal divergence (div_hor routine)43 USE dynadv ! advection (dyn_adv routine)44 USE dynvor ! vorticity term (dyn_vor routine)45 USE dynhpg ! hydrostatic pressure grad. (dyn_hpg routine)46 USE dynldf ! lateral momentum diffusion (dyn_ldf routine)47 USE dynzdf ! vertical diffusion (dyn_zdf routine)48 USE dynspg ! surface pressure gradient (dyn_spg routine)49 50 USE dynatf ! time-filtering (dyn_atf routine)51 52 61 USE stopar ! Stochastic parametrization (sto_par routine) 53 62 USE stopts 54 55 USE bdy_oce , ONLY : ln_bdy56 USE bdydta ! open boundary condition data (bdy_dta routine)57 USE bdytra ! bdy cond. for tracers (bdy_tra routine)58 USE bdydyn3d ! bdy cond. for baroclinic vel. (bdy_dyn3d routine)59 60 USE sshwzv ! vertical velocity and ssh (ssh_nxt routine)61 ! (ssh_swp routine)62 ! (wzv routine)63 USE domvvl ! variable vertical scale factors (dom_vvl_sf_nxt routine)64 ! (dom_vvl_sf_swp routine)65 63 66 64 USE ldfslp ! iso-neutral slopes (ldf_slp routine) … … 68 66 USE ldftra ! lateral eddy diffusive coef. (ldf_tra routine) 69 67 68 USE zdf_oce ! ocean vertical physics variables 70 69 USE zdfphy ! vertical physics manager (zdf_phy_init routine) 71 USE zdfosm , ONLY : osm_rst, dyn_osm, tra_osm ! OSMOSIS routines used in step.F90 70 USE zdfdrg , ONLY : ln_drgimp ! implicit top/bottom friction 71 USE zdfosm , ONLY : osm_rst, dyn_osm, tra_osm ! OSMOSIS routines used in step.F90 72 72 USE zdfmfc ! Mass FLux Convection routine used in step.F90 73 73 … … 83 83 USE diahth ! thermocline depth (dia_hth routine) 84 84 USE diahsb ! heat, salt and volume budgets (dia_hsb routine) 85 USE diacfl 86 USE diaobs ! Observation operator 85 USE diacfl ! CFL diagnostics (dia_cfl routine) 86 USE diaobs ! Observation operator (dia_obs routine) 87 87 USE diadetide ! Weights computation for daily detiding of model diagnostics 88 88 USE diamlr ! IOM context management for multiple-linear-regression analysis … … 94 94 USE asminc ! assimilation increments (tra_asm_inc routine) 95 95 ! (dyn_asm_inc routine) 96 USE asmbkg 96 USE asmbkg ! writing out state trajectory 97 97 USE stpctl ! time stepping control (stp_ctl routine) 98 98 USE restart ! ocean restart (rst_wri routine) -
NEMO/trunk/src/OCE/stpctl.F90
r13616 r14053 26 26 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 27 27 USE lib_mpp ! distributed memory computing 28 !29 28 USE netcdf ! NetCDF library 29 30 30 IMPLICIT NONE 31 31 PRIVATE … … 71 71 CHARACTER(len=20) :: clname 72 72 !!---------------------------------------------------------------------- 73 ! 73 74 IF( nstop > 0 .AND. ngrdstop > -1 ) RETURN ! stpctl was already called by a child grid 74 75 ! … … 179 180 END DO 180 181 IF( kt == nitend ) istatus = NF90_CLOSE(nrunid) 181 END 182 ENDIF 182 183 ! !== error handling ==! 183 184 ! !== done by all processes at every time step ==!
Note: See TracChangeset
for help on using the changeset viewer.
