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 6808 for branches/NERC/dev_r5549_BDY_ZEROGRAD/NEMOGCM/NEMO/OFF_SRC/domrea.F90 – NEMO

Ignore:
Timestamp:
2016-07-19T10:38:35+02:00 (8 years ago)
Author:
jamesharle
Message:

merge with trunk@6232 for consistency with SSB code

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/NERC/dev_r5549_BDY_ZEROGRAD/NEMOGCM/NEMO/OFF_SRC/domrea.F90

    r5504 r6808  
    44   !! Ocean initialization : domain initialization 
    55   !!============================================================================== 
     6   !! History :  OPA  ! 1990-10  (C. Levy - G. Madec)  Original code 
     7   !!                 ! 1992-01  (M. Imbard) insert time step initialization 
     8   !!                 ! 1996-06  (G. Madec) generalized vertical coordinate  
     9   !!                 ! 1997-02  (G. Madec) creation of domwri.F 
     10   !!                 ! 2001-05  (E.Durand - G. Madec) insert closed sea 
     11   !!  NEMO      1.0  ! 2002-08  (G. Madec)  F90: Free form and module 
     12   !!---------------------------------------------------------------------- 
    613 
    714   !!---------------------------------------------------------------------- 
     
    1017   !!   dom_ctl        : control print for the ocean domain 
    1118   !!---------------------------------------------------------------------- 
    12    !! * Modules used 
    1319   USE oce             !  
     20   USE trc_oce         ! shared ocean/biogeochemical variables 
    1421   USE dom_oce         ! ocean space and time domain 
    1522   USE phycst          ! physical constants 
     23   USE domstp          ! domain: set the time-step 
     24   ! 
    1625   USE in_out_manager  ! I/O manager 
    1726   USE lib_mpp         ! distributed memory computing library 
    18  
    19    USE domstp          ! domain: set the time-step 
    20  
    2127   USE lbclnk          ! lateral boundary condition - MPP exchanges 
    22    USE trc_oce         ! shared ocean/biogeochemical variables 
    2328   USE wrk_nemo   
    2429    
     
    2631   PRIVATE 
    2732 
    28    !! * Routine accessibility 
    29    PUBLIC dom_rea       ! called by opa.F90 
     33   PUBLIC   dom_rea    ! called by nemogcm.F90 
    3034 
    3135   !! * Substitutions 
    32 #  include "domzgr_substitute.h90" 
    3336#  include "vectopt_loop_substitute.h90" 
    3437   !!---------------------------------------------------------------------- 
    35    !! NEMO/OFF 3.3 , NEMO Consortium (2010) 
     38   !! NEMO/OFF 3.7 , NEMO Consortium (2015) 
    3639   !! $Id$ 
    3740   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    3841   !!---------------------------------------------------------------------- 
    39  
    4042CONTAINS 
    4143 
     
    5153      !!      - dom_stp: defined the model time step 
    5254      !!      - dom_rea: read the meshmask file if nmsh=1 
    53       !! 
    54       !! History : 
    55       !!        !  90-10  (C. Levy - G. Madec)  Original code 
    56       !!        !  91-11  (G. Madec) 
    57       !!        !  92-01  (M. Imbard) insert time step initialization 
    58       !!        !  96-06  (G. Madec) generalized vertical coordinate  
    59       !!        !  97-02  (G. Madec) creation of domwri.F 
    60       !!        !  01-05  (E.Durand - G. Madec) insert closed sea 
    61       !!   8.5  !  02-08  (G. Madec)  F90: Free form and module 
    62       !!---------------------------------------------------------------------- 
    63       !! * Local declarations 
    64       INTEGER ::   jk                ! dummy loop argument 
    65       INTEGER ::   iconf = 0         ! temporary integers 
    66       !!---------------------------------------------------------------------- 
    67  
     55      !!---------------------------------------------------------------------- 
     56      INTEGER ::   jk          ! dummy loop index 
     57      INTEGER ::   iconf = 0   ! local integers 
     58      !!---------------------------------------------------------------------- 
     59      ! 
    6860      IF(lwp) THEN 
    6961         WRITE(numout,*) 
     
    7163         WRITE(numout,*) '~~~~~~~~' 
    7264      ENDIF 
    73  
    74       CALL dom_nam      ! read namelist ( namrun, namdom, namcla ) 
     65      ! 
     66      CALL dom_nam      ! read namelist ( namrun, namdom ) 
    7567      CALL dom_zgr      ! Vertical mesh and bathymetry option 
    7668      CALL dom_grd      ! Create a domain file 
    77  
    78      ! 
    79       ! - ML - Used in dom_vvl_sf_nxt and lateral diffusion routines 
    80       !        but could be usefull in many other routines 
    81       e12t    (:,:) = e1t(:,:) * e2t(:,:) 
    82       e1e2t   (:,:) = e1t(:,:) * e2t(:,:) 
    83       e12u    (:,:) = e1u(:,:) * e2u(:,:) 
    84       e12v    (:,:) = e1v(:,:) * e2v(:,:) 
    85       e12f    (:,:) = e1f(:,:) * e2f(:,:) 
    86       r1_e12t (:,:) = 1._wp    / e12t(:,:) 
    87       r1_e12u (:,:) = 1._wp    / e12u(:,:) 
    88       r1_e12v (:,:) = 1._wp    / e12v(:,:) 
    89       r1_e12f (:,:) = 1._wp    / e12f(:,:) 
    90       re2u_e1u(:,:) = e2u(:,:) / e1u(:,:) 
    91       re1v_e2v(:,:) = e1v(:,:) / e2v(:,:) 
    92       ! 
    93       hu(:,:) = 0._wp                          ! Ocean depth at U- and V-points 
    94       hv(:,:) = 0._wp 
    95       DO jk = 1, jpk 
    96          hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk) 
    97          hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk) 
     69      ! 
     70      !                                      ! associated horizontal metrics 
     71      ! 
     72      r1_e1t(:,:) = 1._wp / e1t(:,:)   ;   r1_e2t (:,:) = 1._wp / e2t(:,:) 
     73      r1_e1u(:,:) = 1._wp / e1u(:,:)   ;   r1_e2u (:,:) = 1._wp / e2u(:,:) 
     74      r1_e1v(:,:) = 1._wp / e1v(:,:)   ;   r1_e2v (:,:) = 1._wp / e2v(:,:) 
     75      r1_e1f(:,:) = 1._wp / e1f(:,:)   ;   r1_e2f (:,:) = 1._wp / e2f(:,:) 
     76      ! 
     77!!gm BUG if scale factor reduction !!!! 
     78      e1e2t (:,:) = e1t(:,:) * e2t(:,:)   ;   r1_e1e2t(:,:) = 1._wp / e1e2t(:,:) 
     79      e1e2u (:,:) = e1u(:,:) * e2u(:,:)   ;   r1_e1e2u(:,:) = 1._wp / e1e2u(:,:) 
     80      e1e2v (:,:) = e1v(:,:) * e2v(:,:)   ;   r1_e1e2v(:,:) = 1._wp / e1e2v(:,:) 
     81      e1e2f (:,:) = e1f(:,:) * e2f(:,:)   ;   r1_e1e2f(:,:) = 1._wp / e1e2f(:,:) 
     82      !    
     83      e2_e1u(:,:) = e2u(:,:) / e1u(:,:) 
     84      e1_e2v(:,:) = e1v(:,:) / e2v(:,:) 
     85      ! 
     86      hu_n(:,:) = e3u_n(:,:,1) * umask(:,:,1)     ! Ocean depth at U- and V-points 
     87      hv_n(:,:) = e3v_n(:,:,1) * vmask(:,:,1) 
     88      DO jk = 2, jpk 
     89         hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk) 
     90         hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk) 
    9891      END DO 
    9992      !                                        ! Inverse of the local depth 
    100       hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask(:,:,1) ) * umask(:,:,1) 
    101       hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask(:,:,1) ) * vmask(:,:,1) 
    102  
     93      r1_hu_n(:,:) = 1._wp / ( hu_n(:,:) + 1._wp - umask(:,:,1) ) * umask(:,:,1) 
     94      r1_hv_n(:,:) = 1._wp / ( hv_n(:,:) + 1._wp - vmask(:,:,1) ) * vmask(:,:,1) 
     95      ! 
    10396      CALL dom_stp      ! Time step 
    10497      CALL dom_msk      ! Masks 
    10598      CALL dom_ctl      ! Domain control 
    106  
     99      ! 
    107100   END SUBROUTINE dom_rea 
     101 
    108102 
    109103   SUBROUTINE dom_nam 
     
    115109      !! ** input   : - namrun namelist 
    116110      !!              - namdom namelist 
    117       !!              - namcla namelist 
    118111      !!---------------------------------------------------------------------- 
    119112      USE ioipsl 
    120       INTEGER  ::   ios                 ! Local integer output status for namelist read 
    121       NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,               & 
    122          &             nn_no   , cn_exp    , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl,   & 
    123          &             nn_it000, nn_itend  , nn_date0    , nn_leapy     , nn_istate , nn_stock ,   & 
    124          &             nn_write, ln_dimgnnn, ln_mskland  , ln_cfmeta    , ln_clobber, nn_chunksz, nn_euler 
    125       NAMELIST/namdom/ nn_bathy , rn_bathy, rn_e3zps_min, rn_e3zps_rat, nn_msh    , rn_hmin,   & 
    126          &             nn_acc   , rn_atfp     , rn_rdt      , rn_rdtmin ,            & 
    127          &             rn_rdtmax, rn_rdth     , nn_baro     , nn_closea , ln_crs, & 
    128          &             jphgr_msh, & 
     113      INTEGER  ::   ios   ! Local integer output status for namelist read 
     114      ! 
     115      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                         & 
     116         &             nn_no   , cn_exp    , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl,             & 
     117         &             nn_it000, nn_itend  , nn_date0    , nn_time0, nn_leapy     , nn_istate , nn_stock ,   & 
     118         &             nn_write, ln_iscpl  , ln_mskland  , ln_cfmeta    , ln_clobber, nn_chunksz, nn_euler 
     119      NAMELIST/namdom/ nn_bathy , rn_bathy , rn_e3zps_min, rn_e3zps_rat , nn_msh    , rn_hmin   , rn_isfhmin,& 
     120         &             rn_atfp  , rn_rdt   , nn_baro     , nn_closea    , ln_crs    , jphgr_msh,             & 
    129121         &             ppglam0, ppgphi0, ppe1_deg, ppe2_deg, ppe1_m, ppe2_m, & 
    130122         &             ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, & 
    131123         &             ppa2, ppkth2, ppacr2 
    132       NAMELIST/namcla/ nn_cla 
    133124#if defined key_netcdf4 
    134125      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip 
     
    161152         WRITE(numout,*) '      frequency of restart file       nn_stock   = ', nn_stock 
    162153         WRITE(numout,*) '      frequency of output file        nn_write   = ', nn_write 
    163          WRITE(numout,*) '      multi file dimgout              ln_dimgnnn = ', ln_dimgnnn 
    164154         WRITE(numout,*) '      mask land points                ln_mskland = ', ln_mskland 
    165155         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta  = ', ln_cfmeta 
     
    178168      nstocklist = nn_stocklist 
    179169      nwrite = nn_write 
    180  
    181  
     170      ! 
    182171      !                             ! control of output frequency 
    183172      IF ( nstock == 0 .OR. nstock > nitend ) THEN 
     
    194183      ! parameters correspondting to nit000 - 1 (as we start the step loop with a call to day) 
    195184      ndastp = ndate0 - 1        ! ndate0 read in the namelist in dom_nam, we assume that we start run at 00:00 
    196       adatrj = ( REAL( nit000-1, wp ) * rdttra(1) ) / rday 
     185      adatrj = ( REAL( nit000-1, wp ) * rdt ) / rday 
    197186 
    198187#if defined key_agrif 
     
    239228         WRITE(numout,*) '      asselin time filter parameter        rn_atfp   = ', rn_atfp 
    240229         WRITE(numout,*) '      time-splitting: nb of sub time-step  nn_baro   = ', nn_baro 
    241          WRITE(numout,*) '      acceleration of converge             nn_acc    = ', nn_acc 
    242          WRITE(numout,*) '        nn_acc=1: surface tracer rdt       rn_rdtmin = ', rn_rdtmin 
    243          WRITE(numout,*) '                  bottom  tracer rdt       rdtmax    = ', rn_rdtmax 
    244          WRITE(numout,*) '                  depth of transition      rn_rdth   = ', rn_rdth 
    245230         WRITE(numout,*) '      suppression of closed seas (=0)      nn_closea = ', nn_closea 
    246231         WRITE(numout,*) '      type of horizontal mesh jphgr_msh           = ', jphgr_msh 
     
    268253      e3zps_rat = rn_e3zps_rat 
    269254      nmsh      = nn_msh 
    270       nacc      = nn_acc 
    271255      atfp      = rn_atfp 
    272256      rdt       = rn_rdt 
    273       rdtmin    = rn_rdtmin 
    274       rdtmax    = rn_rdtmin 
    275       rdth      = rn_rdth 
    276  
    277       REWIND( numnam_ref )              ! Namelist namcla in reference namelist : Cross land advection 
    278       READ  ( numnam_ref, namcla, IOSTAT = ios, ERR = 905) 
    279 905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcla in reference namelist', lwp ) 
    280  
    281       REWIND( numnam_cfg )              ! Namelist namcla in configuration namelist : Cross land advection 
    282       READ  ( numnam_cfg, namcla, IOSTAT = ios, ERR = 906 ) 
    283 906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcla in configuration namelist', lwp ) 
    284       IF(lwm) WRITE( numond, namcla ) 
    285  
    286       IF(lwp) THEN 
    287          WRITE(numout,*) 
    288          WRITE(numout,*) '   Namelist namcla' 
    289          WRITE(numout,*) '      cross land advection                 nn_cla    = ', nn_cla 
    290       ENDIF 
    291  
    292257#if defined key_netcdf4 
    293258      !                             ! NetCDF 4 case   ("key_netcdf4" defined) 
     
    321286   END SUBROUTINE dom_nam 
    322287 
     288 
    323289   SUBROUTINE dom_zgr 
    324290      !!---------------------------------------------------------------------- 
     
    341307      INTEGER ::   ios 
    342308      !! 
    343       NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco, ln_isfcav 
     309      NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco, ln_isfcav, ln_linssh 
    344310      !!---------------------------------------------------------------------- 
    345311 
     
    362328         WRITE(numout,*) '             s- or hybrid z-s-coordinate    ln_sco    = ', ln_sco 
    363329         WRITE(numout,*) '             ice shelf cavity               ln_isfcav = ', ln_isfcav 
     330         WRITE(numout,*) '             Linear free surface            ln_linssh = ', ln_linssh 
    364331      ENDIF 
    365332 
     
    374341   END SUBROUTINE dom_zgr 
    375342 
     343 
    376344   SUBROUTINE dom_ctl 
    377345      !!---------------------------------------------------------------------- 
     
    382350      !! ** Method  :   compute and print extrema of masked scale factors 
    383351      !! 
    384       !! History : 
    385       !!   8.5  !  02-08  (G. Madec)    Original code 
    386       !!---------------------------------------------------------------------- 
    387       !! * Local declarations 
     352      !!---------------------------------------------------------------------- 
    388353      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2 
    389354      INTEGER, DIMENSION(2) ::   iloc      !  
     
    421386         ijma2 = iloc(2) + njmpp - 1 
    422387      ENDIF 
    423  
     388      ! 
    424389      IF(lwp) THEN 
    425390         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1 
     
    428393         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2 
    429394      ENDIF 
    430  
     395      ! 
    431396   END SUBROUTINE dom_ctl 
     397 
    432398 
    433399   SUBROUTINE dom_grd 
     
    538504         CALL iom_get( inum2, jpdom_data, 'facvolt', facvol ) 
    539505#endif 
    540  
    541506         !                                                         ! horizontal mesh (inum3) 
    542507         CALL iom_get( inum3, jpdom_data, 'glamt', glamt ) 
     
    578543            CALL iom_get( inum4, jpdom_unknown, 'esigw', esigw ) 
    579544 
    580             CALL iom_get( inum4, jpdom_data, 'e3t_0', fse3t_n(:,:,:) ) ! scale factors 
    581             CALL iom_get( inum4, jpdom_data, 'e3u_0', fse3u_n(:,:,:) ) 
    582             CALL iom_get( inum4, jpdom_data, 'e3v_0', fse3v_n(:,:,:) ) 
    583             CALL iom_get( inum4, jpdom_data, 'e3w_0', fse3w_n(:,:,:) ) 
     545            CALL iom_get( inum4, jpdom_data, 'e3t_0', e3t_0(:,:,:) ) ! scale factors 
     546            CALL iom_get( inum4, jpdom_data, 'e3u_0', e3u_0(:,:,:) ) 
     547            CALL iom_get( inum4, jpdom_data, 'e3v_0', e3v_0(:,:,:) ) 
     548            CALL iom_get( inum4, jpdom_data, 'e3w_0', e3w_0(:,:,:) ) 
    584549 
    585550            CALL iom_get( inum4, jpdom_unknown, 'gdept_1d', gdept_1d ) ! depth 
     
    595560            ! 
    596561            IF( nmsh <= 6 ) THEN                                        ! 3D vertical scale factors 
    597                CALL iom_get( inum4, jpdom_data, 'e3t_0', fse3t_n(:,:,:) ) 
    598                CALL iom_get( inum4, jpdom_data, 'e3u_0', fse3u_n(:,:,:) ) 
    599                CALL iom_get( inum4, jpdom_data, 'e3v_0', fse3v_n(:,:,:) ) 
    600                CALL iom_get( inum4, jpdom_data, 'e3w_0', fse3w_n(:,:,:) ) 
     562               CALL iom_get( inum4, jpdom_data, 'e3t_0', e3t_0(:,:,:) ) 
     563               CALL iom_get( inum4, jpdom_data, 'e3u_0', e3u_0(:,:,:) ) 
     564               CALL iom_get( inum4, jpdom_data, 'e3v_0', e3v_0(:,:,:) ) 
     565               CALL iom_get( inum4, jpdom_data, 'e3w_0', e3w_0(:,:,:) ) 
    601566            ELSE                                                        ! 2D bottom scale factors 
    602567               CALL iom_get( inum4, jpdom_data, 'e3t_ps', e3tp ) 
     
    604569               !                                                        ! deduces the 3D scale factors 
    605570               DO jk = 1, jpk 
    606                   fse3t_n(:,:,jk) = e3t_1d(jk)                                    ! set to the ref. factors 
    607                   fse3u_n(:,:,jk) = e3t_1d(jk) 
    608                   fse3v_n(:,:,jk) = e3t_1d(jk) 
    609                   fse3w_n(:,:,jk) = e3w_1d(jk) 
     571                  e3t_0(:,:,jk) = e3t_1d(jk)                                    ! set to the ref. factors 
     572                  e3u_0(:,:,jk) = e3t_1d(jk) 
     573                  e3v_0(:,:,jk) = e3t_1d(jk) 
     574                  e3w_0(:,:,jk) = e3w_1d(jk) 
    610575               END DO 
    611576               DO jj = 1,jpj                                                  ! adjust the deepest values 
    612577                  DO ji = 1,jpi 
    613578                     ik = mbkt(ji,jj) 
    614                      fse3t_n(ji,jj,ik) = e3tp(ji,jj) * tmask(ji,jj,1) + e3t_1d(1) * ( 1._wp - tmask(ji,jj,1) ) 
    615                      fse3w_n(ji,jj,ik) = e3wp(ji,jj) * tmask(ji,jj,1) + e3w_1d(1) * ( 1._wp - tmask(ji,jj,1) ) 
     579                     e3t_0(ji,jj,ik) = e3tp(ji,jj) * tmask(ji,jj,1) + e3t_1d(1) * ( 1._wp - tmask(ji,jj,1) ) 
     580                     e3w_0(ji,jj,ik) = e3wp(ji,jj) * tmask(ji,jj,1) + e3w_1d(1) * ( 1._wp - tmask(ji,jj,1) ) 
    616581                  END DO 
    617582               END DO 
     
    619584                  DO jj = 1, jpjm1 
    620585                     DO ji = 1, jpim1 
    621                         fse3u_n(ji,jj,jk) = MIN( fse3t_n(ji,jj,jk), fse3t_n(ji+1,jj,jk) ) 
    622                         fse3v_n(ji,jj,jk) = MIN( fse3t_n(ji,jj,jk), fse3t_n(ji,jj+1,jk) ) 
     586                        e3u_0(ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji+1,jj,jk) ) 
     587                        e3v_0(ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji,jj+1,jk) ) 
    623588                     END DO 
    624589                  END DO 
    625590               END DO 
    626                CALL lbc_lnk( fse3u_n(:,:,:) , 'U', 1._wp )   ;   CALL lbc_lnk( fse3uw_n(:,:,:), 'U', 1._wp )   ! lateral boundary conditions 
    627                CALL lbc_lnk( fse3v_n(:,:,:) , 'V', 1._wp )   ;   CALL lbc_lnk( fse3vw_n(:,:,:), 'V', 1._wp ) 
     591               CALL lbc_lnk( e3u_0(:,:,:) , 'U', 1._wp )   ;   CALL lbc_lnk( e3uw_0(:,:,:), 'U', 1._wp )   ! lateral boundary conditions 
     592               CALL lbc_lnk( e3v_0(:,:,:) , 'V', 1._wp )   ;   CALL lbc_lnk( e3vw_0(:,:,:), 'V', 1._wp ) 
    628593               ! 
    629594               DO jk = 1, jpk                        ! set to z-scale factor if zero (i.e. along closed boundaries) 
    630                   WHERE( fse3u_n(:,:,jk) == 0._wp )   fse3u_n(:,:,jk) = e3t_1d(jk) 
    631                   WHERE( fse3v_n(:,:,jk) == 0._wp )   fse3v_n(:,:,jk) = e3t_1d(jk) 
     595                  WHERE( e3u_0(:,:,jk) == 0._wp )   e3u_0(:,:,jk) = e3t_1d(jk) 
     596                  WHERE( e3v_0(:,:,jk) == 0._wp )   e3v_0(:,:,jk) = e3t_1d(jk) 
    632597               END DO 
    633598            END IF 
    634599 
    635600            IF( iom_varid( inum4, 'gdept_0', ldstop = .FALSE. ) > 0 ) THEN   ! 3D depth of t- and w-level 
    636                CALL iom_get( inum4, jpdom_data, 'gdept_0', fsdept_n(:,:,:) ) 
    637                CALL iom_get( inum4, jpdom_data, 'gdepw_0', fsdepw_n(:,:,:) ) 
     601               CALL iom_get( inum4, jpdom_data, 'gdept_0', gdept_0(:,:,:) ) 
     602               CALL iom_get( inum4, jpdom_data, 'gdepw_0', gdepw_0(:,:,:) ) 
    638603            ELSE                                                           ! 2D bottom depth 
    639604               CALL iom_get( inum4, jpdom_data, 'hdept', zprt ) 
     
    641606               ! 
    642607               DO jk = 1, jpk                                              ! deduces the 3D depth 
    643                   fsdept_n(:,:,jk) = gdept_1d(jk) 
    644                   fsdepw_n(:,:,jk) = gdepw_1d(jk) 
     608                  gdept_0(:,:,jk) = gdept_1d(jk) 
     609                  gdepw_0(:,:,jk) = gdepw_1d(jk) 
    645610               END DO 
    646611               DO jj = 1, jpj 
     
    648613                     ik = mbkt(ji,jj) 
    649614                     IF( ik > 0 ) THEN 
    650                         fsdepw_n(ji,jj,ik+1) = zprw(ji,jj) 
    651                         fsdept_n(ji,jj,ik  ) = zprt(ji,jj) 
    652                         fsdept_n(ji,jj,ik+1) = fsdept_n(ji,jj,ik) + fse3t_n(ji,jj,ik) 
     615                        gdepw_0(ji,jj,ik+1) = zprw(ji,jj) 
     616                        gdept_0(ji,jj,ik  ) = zprt(ji,jj) 
     617                        gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + e3t_0(ji,jj,ik) 
    653618                     ENDIF 
    654619                  END DO 
     
    664629            CALL iom_get( inum4, jpdom_unknown, 'e3w_1d'  , e3w_1d   ) 
    665630            DO jk = 1, jpk 
    666                fse3t_n(:,:,jk) = e3t_1d(jk)                              ! set to the ref. factors 
    667                fse3u_n(:,:,jk) = e3t_1d(jk) 
    668                fse3v_n(:,:,jk) = e3t_1d(jk) 
    669                fse3w_n(:,:,jk) = e3w_1d(jk) 
    670                fsdept_n(:,:,jk) = gdept_1d(jk) 
    671                fsdepw_n(:,:,jk) = gdepw_1d(jk) 
     631               e3t_0(:,:,jk) = e3t_1d(jk)                              ! set to the ref. factors 
     632               e3u_0(:,:,jk) = e3t_1d(jk) 
     633               e3v_0(:,:,jk) = e3t_1d(jk) 
     634               e3w_0(:,:,jk) = e3w_1d(jk) 
     635               gdept_0(:,:,jk) = gdept_1d(jk) 
     636               gdepw_0(:,:,jk) = gdepw_1d(jk) 
    672637            END DO 
    673638         ENDIF 
     639 
     640      ! 
     641      !              !==  time varying part of coordinate system  ==! 
     642      ! 
     643      !       before        !          now          !       after         ! 
     644      ;  gdept_b = gdept_0  ;   gdept_n = gdept_0   !        ---          ! depth of grid-points 
     645      ;  gdepw_b = gdepw_0  ;   gdepw_n = gdepw_0   !        ---          ! 
     646      ;                     ;   gde3w_n = gde3w_0   !        ---          ! 
     647      ! 
     648      ;    e3t_b =   e3t_0  ;     e3t_n =   e3t_0   ;   e3t_a =  e3t_0    ! scale factors 
     649      ;    e3u_b =   e3u_0  ;     e3u_n =   e3u_0   ;   e3u_a =  e3u_0    ! 
     650      ;    e3v_b =   e3v_0  ;     e3v_n =   e3v_0   ;   e3v_a =  e3v_0    ! 
     651      ;                     ;     e3f_n =   e3f_0   !        ---          ! 
     652      ;    e3w_b =   e3w_0  ;     e3w_n =   e3w_0   !        ---          ! 
     653      ;   e3uw_b =  e3uw_0  ;    e3uw_n =  e3uw_0   !        ---          ! 
     654      ;   e3vw_b =  e3vw_0  ;    e3vw_n =  e3vw_0   !        ---          ! 
     655      ! 
    674656 
    675657!!gm BUG in s-coordinate this does not work! 
     
    701683            &                     e2t  (1,jj), e2u  (1,jj),   & 
    702684            &                     e2v  (1,jj), jj = 1, jpj, 10 ) 
    703       ENDIF 
    704  
    705  
    706       IF( nprint == 1 .AND. lwp ) THEN 
    707          WRITE(numout,*) '          e1u e2u ' 
    708          CALL prihre( e1u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    709          CALL prihre( e2u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    710          WRITE(numout,*) '          e1v e2v  ' 
    711          CALL prihre( e1v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    712          CALL prihre( e2v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    713685      ENDIF 
    714686 
     
    756728      !!                                     (min value = 1 over land) 
    757729      !!---------------------------------------------------------------------- 
    758       ! 
    759730      INTEGER ::   ji, jj   ! dummy loop indices 
    760731      REAL(wp), POINTER, DIMENSION(:,:) :: zmbk 
     
    785756   END SUBROUTINE zgr_bot_level 
    786757 
     758 
    787759   SUBROUTINE dom_msk 
    788760      !!--------------------------------------------------------------------- 
     
    799771      !!               tpol     : ??? 
    800772      !!---------------------------------------------------------------------- 
    801       ! 
    802       INTEGER  ::   ji, jj, jk                   ! dummy loop indices 
    803       INTEGER  ::   iif, iil, ijf, ijl       ! local integers 
     773      INTEGER  ::   ji, jj, jk           ! dummy loop indices 
     774      INTEGER  ::   iif, iil, ijf, ijl   ! local integers 
    804775      INTEGER, POINTER, DIMENSION(:,:) ::  imsk  
    805       ! 
    806776      !!--------------------------------------------------------------------- 
    807777       
     
    839809      DO jj = 1, jpjm1 
    840810         DO ji = 1, fs_jpim1   ! vector loop 
    841             umask_i(ji,jj)  = ssmask(ji,jj) * ssmask(ji+1,jj  )  * MIN(1._wp,SUM(umask(ji,jj,:))) 
    842             vmask_i(ji,jj)  = ssmask(ji,jj) * ssmask(ji  ,jj+1)  * MIN(1._wp,SUM(vmask(ji,jj,:))) 
     811            ssumask(ji,jj)  = ssmask(ji,jj) * ssmask(ji+1,jj  )  * MIN(1._wp,SUM(umask(ji,jj,:))) 
     812            ssvmask(ji,jj)  = ssmask(ji,jj) * ssmask(ji  ,jj+1)  * MIN(1._wp,SUM(vmask(ji,jj,:))) 
    843813         END DO 
    844814         DO ji = 1, jpim1      ! NO vector opt. 
    845             fmask_i(ji,jj) =  ssmask(ji,jj  ) * ssmask(ji+1,jj  )   & 
     815            ssfmask(ji,jj) =  ssmask(ji,jj  ) * ssmask(ji+1,jj  )   & 
    846816               &            * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:))) 
    847817         END DO 
    848818      END DO 
    849       CALL lbc_lnk( umask_i, 'U', 1._wp )      ! Lateral boundary conditions 
    850       CALL lbc_lnk( vmask_i, 'V', 1._wp ) 
    851       CALL lbc_lnk( fmask_i, 'F', 1._wp ) 
     819      CALL lbc_lnk( ssumask, 'U', 1._wp )      ! Lateral boundary conditions 
     820      CALL lbc_lnk( ssvmask, 'V', 1._wp ) 
     821      CALL lbc_lnk( ssfmask, 'F', 1._wp ) 
    852822 
    853823      ! 3. Ocean/land mask at wu-, wv- and w points  
    854824      !---------------------------------------------- 
    855       wmask (:,:,1) = tmask(:,:,1) ! ???????? 
    856       wumask(:,:,1) = umask(:,:,1) ! ???????? 
    857       wvmask(:,:,1) = vmask(:,:,1) ! ???????? 
    858       DO jk=2,jpk 
    859          wmask (:,:,jk)=tmask(:,:,jk) * tmask(:,:,jk-1) 
    860          wumask(:,:,jk)=umask(:,:,jk) * umask(:,:,jk-1)    
    861          wvmask(:,:,jk)=vmask(:,:,jk) * vmask(:,:,jk-1) 
     825      wmask (:,:,1) = tmask(:,:,1)        ! surface value 
     826      wumask(:,:,1) = umask(:,:,1)  
     827      wvmask(:,:,1) = vmask(:,:,1) 
     828      DO jk = 2, jpk                      ! deeper value 
     829         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1) 
     830         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)    
     831         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1) 
    862832      END DO 
    863       ! 
    864       IF( nprint == 1 .AND. lwp ) THEN    ! Control print 
    865          imsk(:,:) = INT( tmask_i(:,:) ) 
    866          WRITE(numout,*) ' tmask_i : ' 
    867          CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, 1, jpj, 1, 1, numout) 
    868          WRITE (numout,*) 
    869          WRITE (numout,*) ' dommsk: tmask for each level' 
    870          WRITE (numout,*) ' ----------------------------' 
    871          DO jk = 1, jpk 
    872             imsk(:,:) = INT( tmask(:,:,jk) ) 
    873             WRITE(numout,*) 
    874             WRITE(numout,*) ' level = ',jk 
    875             CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, 1, jpj, 1, 1, numout) 
    876          END DO 
    877       ENDIF 
    878833      ! 
    879834      CALL wrk_dealloc( jpi, jpj, imsk ) 
Note: See TracChangeset for help on using the changeset viewer.