New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 14053 for NEMO/trunk/src/OCE/DOM/domvvl.F90 – NEMO

Ignore:
Timestamp:
2020-12-03T14:48:38+01:00 (4 years ago)
Author:
techene
Message:

#2385 added to the trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/DOM/domvvl.F90

    r13982 r14053  
    99   !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability 
    1010   !!            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 ratio 
     11   !!             -   ! 2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio 
    1212   !!---------------------------------------------------------------------- 
    1313 
     
    766766      !! ** Purpose :   Read or write VVL file in restart file 
    767767      !! 
    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 
    773782      !!---------------------------------------------------------------------- 
    774783      INTEGER         , INTENT(in) ::   kt        ! ocean time-step 
     
    776785      CHARACTER(len=*), INTENT(in) ::   cdrw      ! "READ"/"WRITE" flag 
    777786      ! 
    778       INTEGER ::   ji, jj, jk 
    779       INTEGER ::   id1, id2, id3, id4, id5     ! local integers 
    780       !!---------------------------------------------------------------------- 
    781       ! 
    782       IF( TRIM(cdrw) == 'READ' ) THEN        ! Read/initialise  
    783          !                                   ! =============== 
    784          IF( ln_rstart ) THEN                   !* Read the restart file 
    785             CALL rst_read_open                  !  open the restart file if necessary 
    786             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  ==! 
    787796            ! 
    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 
    791803            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. ) 
    793805            ! 
    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' 
    799817               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 computed  
    802                IF(lwp) write(numout,*) 'dom_vvl_rst : e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) found in restart files' 
    803818               WHERE ( tmask(:,:,:) == 0.0_wp )  
    804                   e3t(:,:,:,Kmm) = e3t_0(:,:,:) 
    805819                  e3t(:,:,:,Kbb) = e3t_0(:,:,:) 
    806820               END WHERE 
    807                IF( l_1st_euler ) THEN 
    808                   e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    809                ENDIF 
    810             ELSE IF( id1 > 0 ) THEN 
    811                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 ) THEN 
    818                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             ELSE 
    825                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, jpk 
    829                   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 DO 
    833                e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    834                l_1st_euler = .true. 
    835821            ENDIF 
    836             !                             ! ----------- ! 
    837             IF( ln_vvl_zstar ) THEN       ! z_star case ! 
    838                !                          ! ----------- ! 
     822            !                                         ! ------------ ! 
     823            IF( ln_vvl_zstar ) THEN                   ! z_star case ! 
     824               !                                      ! ------------ ! 
    839825               IF( MIN( id3, id4 ) > 0 ) THEN 
    840826                  CALL ctl_stop( 'dom_vvl_rst: z_star cannot restart from a z_tilde or layer run' ) 
    841827               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' 
    847834                  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  
    849843                  tilde_e3t_b(:,:,:) = 0.0_wp 
    850844                  tilde_e3t_n(:,:,:) = 0.0_wp 
    851845               ENDIF 
    852                !                          ! ------------ ! 
    853                IF( ln_vvl_ztilde ) THEN   ! z_tilde case ! 
    854                   !                       ! ------------ ! 
     846               !                                      ! ------------ ! 
     847               IF( ln_vvl_ztilde ) THEN               ! z_tilde case ! 
     848                  !                                   ! ------------ ! 
    855849                  IF( id5 > 0 ) THEN  ! required array exists 
    856850                     CALL iom_get( numror, jpdom_auto, 'hdiv_lf', hdiv_lf(:,:,:) ) 
    857851                  ELSE                ! array is missing 
    858                      hdiv_lf(:,:,:) = 0.0_wp 
     852                     hdiv_lf(:,:,:) = 0.0_wp  
    859853                  ENDIF 
    860854               ENDIF 
    861855            ENDIF 
    862856            ! 
    863          ELSE                                   !* Initialize at "rest" 
     857         ELSE                                   !==  Initialize at "rest" with ssh  ==! 
    864858            ! 
    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            ! 
    919864            IF( ln_vvl_ztilde .OR. ln_vvl_layer) THEN 
    920865               tilde_e3t_b(:,:,:) = 0._wp 
    921866               tilde_e3t_n(:,:,:) = 0._wp 
    922867               IF( ln_vvl_ztilde ) hdiv_lf(:,:,:) = 0._wp 
    923             END IF 
     868            ENDIF 
    924869         ENDIF 
    925          ! 
    926       ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN   ! Create restart file 
    927          !                                   ! =================== 
     870         !                                       !=======================! 
     871      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN       !  Create restart file  ! 
     872         !                                       !=======================! 
     873         ! 
    928874         IF(lwp) WRITE(numout,*) '---- dom_vvl_rst ----' 
    929875         !                                           ! --------- ! 
Note: See TracChangeset for help on using the changeset viewer.