Changeset 14053 for NEMO/trunk/src/OCE/DOM
- Timestamp:
- 2020-12-03T14:48:38+01:00 (4 years ago)
- Location:
- NEMO/trunk/src/OCE/DOM
- Files:
-
- 8 edited
-
dom_oce.F90 (modified) (7 diffs)
-
domain.F90 (modified) (18 diffs)
-
dommsk.F90 (modified) (1 diff)
-
domqco.F90 (modified) (11 diffs)
-
domvvl.F90 (modified) (3 diffs)
-
domzgr_substitute.h90 (modified) (2 diffs)
-
istate.F90 (modified) (5 diffs)
-
phycst.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
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)
Note: See TracChangeset
for help on using the changeset viewer.
