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 5836 for trunk/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

Ignore:
Timestamp:
2015-10-26T15:49:40+01:00 (9 years ago)
Author:
cetlod
Message:

merge the simplification branch onto the trunk, see ticket #1612

Location:
trunk/NEMOGCM/NEMO/OPA_SRC/DOM
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90

    r5506 r5836  
    7171      !!                                   =2 put at location runoff 
    7272      !!---------------------------------------------------------------------- 
    73       INTEGER ::   jc            ! dummy loop indices 
    74       INTEGER :: isrow           ! local index 
    75       !!---------------------------------------------------------------------- 
    76        
     73      INTEGER ::   jc      ! dummy loop indices 
     74      INTEGER ::   isrow   ! local index 
     75      !!---------------------------------------------------------------------- 
     76      ! 
    7777      IF(lwp) WRITE(numout,*) 
    7878      IF(lwp) WRITE(numout,*)'dom_clo : closed seas ' 
    7979      IF(lwp) WRITE(numout,*)'~~~~~~~' 
    80  
     80      ! 
    8181      ! initial values 
    8282      ncsnr(:) = 1  ;  ncsi1(:) = 1  ;  ncsi2(:) = 1  ;  ncsir(:,:) = 1 
    8383      ncstt(:) = 0  ;  ncsj1(:) = 1  ;  ncsj2(:) = 1  ;  ncsjr(:,:) = 1 
    84  
     84      ! 
    8585      ! set the closed seas (in data domain indices) 
    8686      ! ------------------- 
    87  
     87      ! 
    8888      IF( cp_cfg == "orca" ) THEN 
    8989         ! 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90

    r5123 r5836  
    77   !! History :  1.0  ! 2005-10  (A. Beckmann, G. Madec)  reactivate s-coordinate  
    88   !!            3.3  ! 2010-11  (G. Madec) add mbk. arrays associated to the deepest ocean level 
    9    !!            4.0  ! 2011-01  (A. R. Porter, STFC Daresbury) dynamical allocation 
     9   !!            3.4  ! 2011-01  (A. R. Porter, STFC Daresbury) dynamical allocation 
    1010   !!            3.5  ! 2012     (S. Mocavero, I. Epicoco) Add arrays associated 
    1111   !!                             to the optimization of BDY communications 
     12   !!            3.7  ! 2015-11  (G. Madec) introduce surface and scale factor ratio 
    1213   !!---------------------------------------------------------------------- 
    1314 
     
    2021 
    2122   IMPLICIT NONE 
    22    PUBLIC             ! allows the acces to par_oce when dom_oce is used 
    23    !                  ! exception to coding rules... to be suppressed ??? 
     23   PUBLIC             ! allows the acces to par_oce when dom_oce is used (exception to coding rules) 
    2424 
    2525   PUBLIC dom_oce_alloc  ! Called from nemogcm.F90 
     
    107107   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   r2dtra  !: = 2*rdttra except at nit000 (=rdttra) if neuler=0 
    108108 
    109    !                                         !!* Namelist namcla : cross land advection 
    110    INTEGER, PUBLIC ::   nn_cla               !: =1 cross land advection for exchanges through some straits (ORCA2) 
    111  
    112109   !!---------------------------------------------------------------------- 
    113110   !! space domain parameters 
     
    158155   !! horizontal curvilinear coordinate and scale factors 
    159156   !! --------------------------------------------------------------------- 
    160    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  glamt, glamu   !: longitude of t-, u-, v- and f-points (degre) 
    161    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  glamv, glamf   !: 
    162    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  gphit, gphiu   !: latitude  of t-, u-, v- and f-points (degre) 
    163    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  gphiv, gphif   !: 
    164    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::  e1t, e2t, r1_e1t, r1_e2t   !: horizontal scale factors and inverse at t-point (m) 
    165    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::  e1u, e2u, r1_e1u, r1_e2u   !: horizontal scale factors and inverse at u-point (m) 
    166    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::  e1v, e2v, r1_e1v, r1_e2v   !: horizontal scale factors and inverse at v-point (m) 
    167    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::  e1f, e2f, r1_e1f, r1_e2f   !: horizontal scale factors and inverse at f-point (m) 
    168    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  e1e2t          !: surface at t-point (m2) 
    169    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  ff             !: coriolis factor (2.*omega*sin(yphi) ) (s-1) 
     157   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   glamt , glamu, glamv , glamf    !: longitude at t, u, v, f-points [degree] 
     158   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   gphit , gphiu, gphiv , gphif    !: latitude  at t, u, v, f-points [degree] 
     159   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::   e1t   , e2t  , r1_e1t, r1_e2t   !: t-point horizontal scale factors    [m] 
     160   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::   e1u   , e2u  , r1_e1u, r1_e2u   !: horizontal scale factors at u-point [m] 
     161   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::   e1v   , e2v  , r1_e1v, r1_e2v   !: horizontal scale factors at v-point [m] 
     162   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) ::   e1f   , e2f  , r1_e1f, r1_e2f   !: horizontal scale factors at f-point [m] 
     163   ! 
     164   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e1e2t , r1_e1e2t                !: associated metrics at t-point 
     165   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e1e2u , r1_e1e2u , e2_e1u       !: associated metrics at u-point 
     166   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e1e2v , r1_e1e2v , e1_e2v       !: associated metrics at v-point 
     167   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e1e2f , r1_e1e2f                !: associated metrics at f-point 
     168   ! 
     169   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ff                              !: coriolis factor                   [1/s] 
    170170 
    171171   !!---------------------------------------------------------------------- 
     
    216216   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ht_0           !: reference depth at t-       points (meters) 
    217217   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hu_0 , hv_0    !: reference depth at u- and v-points (meters) 
    218    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   re2u_e1u       !: scale factor coeffs at u points (e2u/e1u) 
    219    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   re1v_e2v       !: scale factor coeffs at v points (e1v/e2v) 
    220    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e12t , r1_e12t !: horizontal cell surface and inverse at t points 
    221    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e12u , r1_e12u !: horizontal cell surface and inverse at u points 
    222    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e12v , r1_e12v !: horizontal cell surface and inverse at v points 
    223    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e12f , r1_e12f !: horizontal cell surface and inverse at f points 
    224218 
    225219   INTEGER, PUBLIC ::   nla10              !: deepest    W level Above  ~10m (nlb10 - 1) 
     
    265259 
    266260   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   tpol, fpol          !: north fold mask (jperio= 3 or 4) 
    267  
    268 #if defined key_noslip_accurate 
    269    INTEGER, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:  )  :: npcoa              !: ??? 
    270    INTEGER, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  :: nicoa, njcoa       !: ??? 
    271 #endif 
    272261 
    273262   !!---------------------------------------------------------------------- 
     
    333322   INTEGER FUNCTION dom_oce_alloc() 
    334323      !!---------------------------------------------------------------------- 
    335       INTEGER, DIMENSION(12) :: ierr 
     324      INTEGER, DIMENSION(13) :: ierr 
    336325      !!---------------------------------------------------------------------- 
    337326      ierr(:) = 0 
     
    346335         &      tpol(jpiglo)  , fpol(jpiglo)                               , STAT=ierr(2) ) 
    347336         ! 
    348       ALLOCATE( glamt(jpi,jpj) , gphit(jpi,jpj) , e1t(jpi,jpj) , e2t(jpi,jpj) , r1_e1t(jpi,jpj) , r1_e2t(jpi,jpj) ,   &  
    349          &      glamu(jpi,jpj) , gphiu(jpi,jpj) , e1u(jpi,jpj) , e2u(jpi,jpj) , r1_e1u(jpi,jpj) , r1_e2u(jpi,jpj) ,   &   
    350          &      glamv(jpi,jpj) , gphiv(jpi,jpj) , e1v(jpi,jpj) , e2v(jpi,jpj) , r1_e1v(jpi,jpj) , r1_e2v(jpi,jpj) ,   &   
    351          &      glamf(jpi,jpj) , gphif(jpi,jpj) , e1f(jpi,jpj) , e2f(jpi,jpj) , r1_e1f(jpi,jpj) , r1_e2f(jpi,jpj) ,   & 
    352          &      e1e2t(jpi,jpj) , ff   (jpi,jpj) , STAT=ierr(3) )      
     337      ALLOCATE( glamt(jpi,jpj) ,    glamu(jpi,jpj) ,  glamv(jpi,jpj) ,  glamf(jpi,jpj) ,     & 
     338         &      gphit(jpi,jpj) ,    gphiu(jpi,jpj) ,  gphiv(jpi,jpj) ,  gphif(jpi,jpj) ,     & 
     339         &       e1t (jpi,jpj) ,     e2t (jpi,jpj) , r1_e1t(jpi,jpj) , r1_e2t(jpi,jpj) ,     & 
     340         &       e1u (jpi,jpj) ,     e2u (jpi,jpj) , r1_e1u(jpi,jpj) , r1_e2u(jpi,jpj) ,     & 
     341         &       e1v (jpi,jpj) ,     e2v (jpi,jpj) , r1_e1v(jpi,jpj) , r1_e2v(jpi,jpj) ,     & 
     342         &       e1f (jpi,jpj) ,     e2f (jpi,jpj) , r1_e1f(jpi,jpj) , r1_e2f(jpi,jpj) ,     & 
     343         &      e1e2t(jpi,jpj) , r1_e1e2t(jpi,jpj)                                     ,     & 
     344         &      e1e2u(jpi,jpj) , r1_e1e2u(jpi,jpj) , e2_e1u(jpi,jpj)                   ,     & 
     345         &      e1e2v(jpi,jpj) , r1_e1e2v(jpi,jpj) , e1_e2v(jpi,jpj)                   ,     & 
     346         &      e1e2f(jpi,jpj) , r1_e1e2f(jpi,jpj)                                     ,     & 
     347         &        ff (jpi,jpj)                                                         , STAT=ierr(3) ) 
    353348         ! 
    354349      ALLOCATE( gdep3w_0(jpi,jpj,jpk) , e3v_0(jpi,jpj,jpk) , e3f_0 (jpi,jpj,jpk) ,                         & 
     
    364359         &      gdept_b (jpi,jpj,jpk) ,gdepw_b(jpi,jpj,jpk) , e3w_b (jpi,jpj,jpk) ,                           & 
    365360         &      e3t_a   (jpi,jpj,jpk) , e3u_a (jpi,jpj,jpk) , e3v_a (jpi,jpj,jpk) ,                           & 
    366          &      ehu_a    (jpi,jpj)    , ehv_a (jpi,jpj),                                                     & 
    367          &      ehur_a   (jpi,jpj)    , ehvr_a (jpi,jpj),                                                     & 
    368          &      ehu_b    (jpi,jpj)    , ehv_b (jpi,jpj),                                                     & 
    369          &      ehur_b   (jpi,jpj)    , ehvr_b (jpi,jpj),                                  STAT=ierr(5) )                           
     361         &      ehu_a   (jpi,jpj)     , ehv_a (jpi,jpj),                                                     & 
     362         &      ehur_a  (jpi,jpj)     , ehvr_a(jpi,jpj),                                                     & 
     363         &      ehu_b   (jpi,jpj)     , ehv_b (jpi,jpj),                                                     & 
     364         &      ehur_b  (jpi,jpj)     , ehvr_b(jpi,jpj),                                  STAT=ierr(5) )                           
    370365#endif 
    371366         ! 
    372       ALLOCATE( hu      (jpi,jpj) , hur     (jpi,jpj) , hu_0(jpi,jpj) , ht_0  (jpi,jpj) ,     & 
    373          &      hv      (jpi,jpj) , hvr     (jpi,jpj) , hv_0(jpi,jpj) , ht    (jpi,jpj) ,     & 
    374          &      re2u_e1u(jpi,jpj) , re1v_e2v(jpi,jpj) ,                                       & 
    375          &      e12t    (jpi,jpj) , r1_e12t (jpi,jpj) ,                                       & 
    376          &      e12u    (jpi,jpj) , r1_e12u (jpi,jpj) ,                                       & 
    377          &      e12v    (jpi,jpj) , r1_e12v (jpi,jpj) ,                                       & 
    378          &      e12f    (jpi,jpj) , r1_e12f (jpi,jpj) ,                                   STAT=ierr(6)  ) 
     367      ALLOCATE( hu(jpi,jpj) , hur(jpi,jpj) , hu_0(jpi,jpj) , ht_0(jpi,jpj) ,     & 
     368         &      hv(jpi,jpj) , hvr(jpi,jpj) , hv_0(jpi,jpj) , ht  (jpi,jpj) , STAT=ierr(6)  ) 
    379369         ! 
    380370      ALLOCATE( gdept_1d(jpk) , gdepw_1d(jpk) ,                                     & 
     
    387377         &      scosrf(jpi,jpj) , scobot(jpi,jpj) ,     & 
    388378         &      hifv  (jpi,jpj) , hiff  (jpi,jpj) ,     & 
    389          &      hift  (jpi,jpj) , hifu  (jpi,jpj) , rx1 (jpi,jpj) , STAT=ierr(8) ) 
     379         &      hift  (jpi,jpj) , hifu  (jpi,jpj) , rx1(jpi,jpj) , STAT=ierr(8) ) 
    390380 
    391381      ALLOCATE( mbathy(jpi,jpj) , bathy(jpi,jpj) ,                                      & 
    392382         &     tmask_i(jpi,jpj) , umask_i(jpi,jpj), vmask_i(jpi,jpj), fmask_i(jpi,jpj), & 
    393          &     bmask(jpi,jpj)  ,                                                       & 
     383         &     bmask  (jpi,jpj) ,                                                       & 
    394384         &     mbkt   (jpi,jpj) , mbku (jpi,jpj) , mbkv(jpi,jpj) , STAT=ierr(9) ) 
    395385 
    396386! (ISF) Allocation of basic array    
    397       ALLOCATE( misfdep(jpi,jpj) , risfdep(jpi,jpj),     & 
    398          &     mikt(jpi,jpj), miku(jpi,jpj), mikv(jpi,jpj) ,           & 
    399          &     mikf(jpi,jpj), ssmask(jpi,jpj), STAT=ierr(10) ) 
     387      ALLOCATE( misfdep(jpi,jpj) , risfdep(jpi,jpj),                   & 
     388         &      mikt(jpi,jpj), miku(jpi,jpj), mikv(jpi,jpj) ,           & 
     389         &      mikf(jpi,jpj), ssmask(jpi,jpj), STAT=ierr(10) ) 
    400390 
    401391      ALLOCATE( tmask(jpi,jpj,jpk) , umask(jpi,jpj,jpk),     &  
     
    403393 
    404394      ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(12) ) 
    405  
    406 #if defined key_noslip_accurate 
    407       ALLOCATE( npcoa(4,jpk), nicoa(2*(jpi+jpj),4,jpk), njcoa(2*(jpi+jpj),4,jpk), STAT=ierr(12) ) 
    408 #endif 
    409395      ! 
    410396      dom_oce_alloc = MAXVAL(ierr) 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r5363 r5836  
    1919   !!   dom_nam        : read and contral domain namelists 
    2020   !!   dom_ctl        : control print for the ocean domain 
     21   !!   dom_stiff      : diagnose maximum grid stiffness/hydrostatic consistency (s-coordinate) 
    2122   !!---------------------------------------------------------------------- 
    2223   USE oce             ! ocean variables 
     
    2526   USE phycst          ! physical constants 
    2627   USE closea          ! closed seas 
    27    USE in_out_manager  ! I/O manager 
    28    USE lib_mpp         ! distributed memory computing library 
    29  
    3028   USE domhgr          ! domain: set the horizontal mesh 
    3129   USE domzgr          ! domain: set the vertical mesh 
     
    3634   USE c1d             ! 1D vertical configuration 
    3735   USE dyncor_c1d      ! Coriolis term (c1d case)         (cor_c1d routine) 
     36   ! 
     37   USE in_out_manager  ! I/O manager 
     38   USE lib_mpp         ! distributed memory computing library 
     39   USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    3840   USE timing          ! Timing 
    39    USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    4041 
    4142   IMPLICIT NONE 
     
    8182      ENDIF 
    8283      ! 
    83                              CALL dom_nam      ! read namelist ( namrun, namdom, namcla ) 
     84                             CALL dom_nam      ! read namelist ( namrun, namdom ) 
    8485                             CALL dom_clo      ! Closed seas and lake 
    8586                             CALL dom_hgr      ! Horizontal mesh 
     
    8889      IF( ln_sco )           CALL dom_stiff    ! Maximum stiffness ratio/hydrostatic consistency 
    8990      ! 
    90       ht_0(:,:) = 0.0_wp                       ! Reference ocean depth at T-points 
    91       hu_0(:,:) = 0.0_wp                       ! Reference ocean depth at U-points 
    92       hv_0(:,:) = 0.0_wp                       ! Reference ocean depth at V-points 
     91      ht_0(:,:) = 0._wp                        ! Reference ocean depth at T-points 
     92      hu_0(:,:) = 0._wp                        ! Reference ocean depth at U-points 
     93      hv_0(:,:) = 0._wp                        ! Reference ocean depth at V-points 
    9394      DO jk = 1, jpk 
    9495         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk) 
     
    9798      END DO 
    9899      ! 
    99       IF( lk_vvl )           CALL dom_vvl_init ! Vertical variable mesh 
     100      IF( lk_vvl         )   CALL dom_vvl_init ! Vertical variable mesh 
    100101      ! 
    101102      IF( lk_c1d         )   CALL cor_c1d      ! 1D configuration: Coriolis set at T-point 
     
    131132      !! ** input   : - namrun namelist 
    132133      !!              - namdom namelist 
    133       !!              - namcla namelist 
    134134      !!              - namnc4 namelist   ! "key_netcdf4" only 
    135135      !!---------------------------------------------------------------------- 
     
    146146         &             ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, & 
    147147         &             ppa2, ppkth2, ppacr2 
    148       NAMELIST/namcla/ nn_cla 
    149148#if defined key_netcdf4 
    150149      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip 
     
    155154      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run 
    156155      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901) 
    157 901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in reference namelist', lwp ) 
     156901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in reference namelist', lwp ) 
    158157 
    159158      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run 
    160159      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 ) 
    161 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp ) 
     160902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp ) 
    162161      IF(lwm) WRITE ( numond, namrun ) 
    163162      ! 
     
    251250904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp ) 
    252251      IF(lwm) WRITE ( numond, namdom ) 
    253  
     252      ! 
    254253      IF(lwp) THEN 
    255254         WRITE(numout,*) 
     
    293292         WRITE(numout,*) '                                      ppacr2            = ', ppacr2 
    294293      ENDIF 
    295  
     294      ! 
    296295      ntopo     = nn_bathy          ! conversion DOCTOR names into model names (this should disappear soon) 
    297296      e3zps_min = rn_e3zps_min 
     
    304303      rdtmax    = rn_rdtmin 
    305304      rdth      = rn_rdth 
    306  
    307       REWIND( numnam_ref )              ! Namelist namcla in reference namelist : Cross land advection 
    308       READ  ( numnam_ref, namcla, IOSTAT = ios, ERR = 905) 
    309 905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcla in reference namelist', lwp ) 
    310  
    311       REWIND( numnam_cfg )              ! Namelist namcla in configuration namelist : Cross land advection 
    312       READ  ( numnam_cfg, namcla, IOSTAT = ios, ERR = 906 ) 
    313 906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcla in configuration namelist', lwp ) 
    314       IF(lwm) WRITE( numond, namcla ) 
    315  
    316       IF(lwp) THEN 
    317          WRITE(numout,*) 
    318          WRITE(numout,*) '   Namelist namcla' 
    319          WRITE(numout,*) '      cross land advection                 nn_cla    = ', nn_cla 
    320       ENDIF 
    321       IF ( nn_cla .EQ. 1 ) THEN 
    322          IF  ( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA R2  
    323             CONTINUE 
    324          ELSE 
    325             CALL ctl_stop( 'STOP', 'Cross land advation iplemented only for ORCA2 configuration: cp_cfg = "orca" and jp_cfg = 2 ' ) 
    326          ENDIF 
    327       ENDIF 
    328305 
    329306#if defined key_netcdf4 
     
    409386   END SUBROUTINE dom_ctl 
    410387 
     388 
    411389   SUBROUTINE dom_stiff 
    412390      !!---------------------------------------------------------------------- 
     
    427405      REAL(wp), DIMENSION(4) :: zr1 
    428406      !!---------------------------------------------------------------------- 
    429       rx1(:,:) = 0.e0 
    430       zrxmax   = 0.e0 
    431       zr1(:)   = 0.e0 
    432        
     407      rx1(:,:) = 0._wp 
     408      zrxmax   = 0._wp 
     409      zr1(:)   = 0._wp 
     410      ! 
    433411      DO ji = 2, jpim1 
    434412         DO jj = 2, jpjm1 
     
    455433         END DO 
    456434      END DO 
    457  
    458435      CALL lbc_lnk( rx1, 'T', 1. ) 
    459  
    460       zrxmax = MAXVAL(rx1) 
    461  
     436      ! 
     437      zrxmax = MAXVAL( rx1 ) 
     438      ! 
    462439      IF( lk_mpp )   CALL mpp_max( zrxmax ) ! max over the global domain 
    463  
     440      ! 
    464441      IF(lwp) THEN 
    465442         WRITE(numout,*) 
     
    467444         WRITE(numout,*) '~~~~~~~~~' 
    468445      ENDIF 
    469  
     446      ! 
    470447   END SUBROUTINE dom_stiff 
    471  
    472  
    473448 
    474449   !!====================================================================== 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domhgr.F90

    r5656 r5836  
    1414   !!                            use of parameters in par_CONFIG-Rxx.h90, not in namelist 
    1515   !!             -   ! 2004-05  (A. Koch-Larrouy) Add Gyre configuration  
    16    !!            4.0  ! 2011-02  (G. Madec) add cell surface (e1e2t) 
     16   !!            3.7  ! 2015-09  (G. Madec, S. Flavoni) add cell surface and their inverse 
     17   !!                                       add optional read of e1e2u & e1e2v 
    1718   !!---------------------------------------------------------------------- 
    1819 
     
    2324   USE dom_oce        ! ocean space and time domain 
    2425   USE phycst         ! physical constants 
     26   USE domwri         ! write 'meshmask.nc' & 'coordinate_e1e2u_v.nc' files 
     27   ! 
    2528   USE in_out_manager ! I/O manager 
    2629   USE lib_mpp        ! MPP library 
     
    3538 
    3639   !!---------------------------------------------------------------------- 
    37    !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
     40   !! NEMO/OPA 3.7 , NEMO Consortium (2014) 
    3841   !! $Id$  
    3942   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    106109      REAL(wp) ::   zphi1, zsin_alpha, zim05, zjm05 
    107110      INTEGER  ::   isrow                ! index for ORCA1 starting row 
    108  
     111      INTEGER  ::   ie1e2u_v             ! fag for u- & v-surface read in coordinate file or not 
    109112      !!---------------------------------------------------------------------- 
    110113      ! 
     
    122125         WRITE(numout,*) '             meridional grid-spacing (meters)  ppe2_m   = ', ppe2_m   
    123126      ENDIF 
    124  
    125  
    126       SELECT CASE( jphgr_msh )   ! type of horizontal mesh 
    127  
    128       CASE ( 0 )                     !  curvilinear coordinate on the sphere read in coordinate.nc file 
    129  
     127      ! 
     128      ! 
     129      SELECT CASE( jphgr_msh )   !  type of horizontal mesh   
     130      ! 
     131      CASE ( 0 )                     !==  read in coordinate.nc file  ==! 
     132         ! 
    130133         IF(lwp) WRITE(numout,*) 
    131134         IF(lwp) WRITE(numout,*) '          curvilinear coordinate on the sphere read in "coordinate" file' 
    132  
    133          CALL hgr_read           ! Defaultl option  :   NetCDF file 
    134  
    135          !                                                ! ===================== 
    136          IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN    ! ORCA R2 configuration 
    137             !                                             ! ===================== 
    138             IF( nn_cla == 0 ) THEN 
    139                ! 
    140                ii0 = 139   ;   ii1 = 140        ! Gibraltar Strait (e2u = 20 km) 
    141                ij0 = 102   ;   ij1 = 102   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  20.e3 
    142                IF(lwp) WRITE(numout,*) 
    143                IF(lwp) WRITE(numout,*) '             orca_r2: Gibraltar    : e2u reduced to 20 km' 
    144                ! 
    145                ii0 = 160   ;   ii1 = 160        ! Bab el Mandeb (e2u = 18 km) 
    146                ij0 =  88   ;   ij1 =  88   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  18.e3 
    147                                                e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  30.e3 
    148                IF(lwp) WRITE(numout,*) 
    149                IF(lwp) WRITE(numout,*) '             orca_r2: Bab el Mandeb: e2u reduced to 30 km' 
    150                IF(lwp) WRITE(numout,*) '                                     e1v reduced to 18 km' 
    151             ENDIF 
    152  
    153             ii0 = 145   ;   ii1 = 146        ! Danish Straits (e2u = 10 km) 
    154             ij0 = 116   ;   ij1 = 116   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
    155             IF(lwp) WRITE(numout,*) 
    156             IF(lwp) WRITE(numout,*) '             orca_r2: Danish Straits : e2u reduced to 10 km' 
    157             ! 
     135         ! 
     136         ie1e2u_v = 0                  ! set to unread e1e2u and e1e2v 
     137         ! 
     138         CALL hgr_read( ie1e2u_v )     ! read the coordinate.nc file 
     139         ! 
     140         IF( ie1e2u_v == 0 ) THEN      ! e1e2u and e1e2v have not been read: compute them 
     141            !                          ! e2u and e1v does not include a reduction in some strait: apply reduction 
     142            e1e2u (:,:) = e1u(:,:) * e2u(:,:)    
     143            e1e2v (:,:) = e1v(:,:) * e2v(:,:)  
    158144         ENDIF 
    159  
    160             !                                             ! ===================== 
    161          IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN    ! ORCA R1 configuration 
    162             !                                             ! ===================== 
    163             ! This dirty section will be suppressed by simplification process: all this will come back in input files 
    164             ! Currently these hard-wired indices relate to configuration with 
    165             ! extend grid (jpjglo=332) 
    166             ! which had a grid-size of 362x292. 
    167             !  
    168             isrow = 332 - jpjglo 
    169             ! 
    170             ii0 = 282           ;   ii1 = 283        ! Gibraltar Strait (e2u = 20 km) 
    171             ij0 = 241 - isrow   ;   ij1 = 241 - isrow   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  20.e3 
    172             IF(lwp) WRITE(numout,*) 
    173             IF(lwp) WRITE(numout,*) '             orca_r1: Gibraltar : e2u reduced to 20 km' 
    174  
    175             ii0 = 314           ;   ii1 = 315        ! Bhosporus Strait (e2u = 10 km) 
    176             ij0 = 248 - isrow   ;   ij1 = 248 - isrow   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
    177             IF(lwp) WRITE(numout,*) 
    178             IF(lwp) WRITE(numout,*) '             orca_r1: Bhosporus : e2u reduced to 10 km' 
    179  
    180             ii0 =  44           ;   ii1 =  44        ! Lombok Strait (e1v = 13 km) 
    181             ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  13.e3 
    182             IF(lwp) WRITE(numout,*) 
    183             IF(lwp) WRITE(numout,*) '             orca_r1: Lombok : e1v reduced to 10 km' 
    184  
    185             ii0 =  48           ;   ii1 =  48        ! Sumba Strait (e1v = 8 km) [closed from bathy_11 on] 
    186             ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  8.e3 
    187             IF(lwp) WRITE(numout,*) 
    188             IF(lwp) WRITE(numout,*) '             orca_r1: Sumba : e1v reduced to 8 km' 
    189  
    190             ii0 =  53           ;   ii1 =  53        ! Ombai Strait (e1v = 13 km) 
    191             ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3 
    192             IF(lwp) WRITE(numout,*) 
    193             IF(lwp) WRITE(numout,*) '             orca_r1: Ombai : e1v reduced to 13 km' 
    194  
    195             ii0 =  56           ;   ii1 =  56        ! Timor Passage (e1v = 20 km) 
    196             ij0 = 164 - isrow   ;   ij1 = 145 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3 
    197             IF(lwp) WRITE(numout,*) 
    198             IF(lwp) WRITE(numout,*) '             orca_r1: Timor Passage : e1v reduced to 20 km' 
    199  
    200             ii0 =  55           ;   ii1 =  55        ! West Halmahera Strait (e1v = 30 km) 
    201             ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 30.e3 
    202             IF(lwp) WRITE(numout,*) 
    203             IF(lwp) WRITE(numout,*) '             orca_r1: W Halmahera : e1v reduced to 30 km' 
    204  
    205             ii0 =  58           ;   ii1 =  58        ! East Halmahera Strait (e1v = 50 km) 
    206             ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 50.e3 
    207             IF(lwp) WRITE(numout,*) 
    208             IF(lwp) WRITE(numout,*) '             orca_r1: E Halmahera : e1v reduced to 50 km' 
    209             ! 
    210             ! 
    211          ENDIF 
    212  
    213          !                                                ! ====================== 
    214          IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN   ! ORCA R05 configuration 
    215             !                                             ! ====================== 
    216             ii0 = 563   ;   ii1 = 564        ! Gibraltar Strait (e2u = 20 km) 
    217             ij0 = 327   ;   ij1 = 327   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  20.e3 
    218             IF(lwp) WRITE(numout,*) 
    219             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Gibraltar Strait' 
    220             ! 
    221             ii0 = 627   ;   ii1 = 628        ! Bosphore Strait (e2u = 10 km) 
    222             ij0 = 343   ;   ij1 = 343   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
    223             IF(lwp) WRITE(numout,*) 
    224             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Bosphore Strait' 
    225             ! 
    226             ii0 =  93   ;   ii1 =  94        ! Sumba Strait (e2u = 40 km) 
    227             ij0 = 232   ;   ij1 = 232   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  40.e3 
    228             IF(lwp) WRITE(numout,*) 
    229             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Sumba Strait' 
    230             ! 
    231             ii0 = 103   ;   ii1 = 103        ! Ombai Strait (e2u = 15 km) 
    232             ij0 = 232   ;   ij1 = 232   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  15.e3 
    233             IF(lwp) WRITE(numout,*) 
    234             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Ombai Strait' 
    235             ! 
    236             ii0 =  15   ;   ii1 =  15        ! Palk Strait (e2u = 10 km) 
    237             ij0 = 270   ;   ij1 = 270   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
    238             IF(lwp) WRITE(numout,*) 
    239             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Palk Strait' 
    240             ! 
    241             ii0 =  87   ;   ii1 =  87        ! Lombok Strait (e1v = 10 km) 
    242             ij0 = 232   ;   ij1 = 233   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
    243             IF(lwp) WRITE(numout,*) 
    244             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e1v at the Lombok Strait' 
    245             ! 
    246             ! 
    247             ii0 = 662   ;   ii1 = 662        ! Bab el Mandeb (e1v = 25 km) 
    248             ij0 = 276   ;   ij1 = 276   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  25.e3 
    249             IF(lwp) WRITE(numout,*) 
    250             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e1v at the Bab el Mandeb' 
    251             ! 
    252          ENDIF 
    253  
    254  
    255          ! N.B. :  General case, lat and long function of both i and j indices: 
    256          !     e1t(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphit(ji,jj) ) * fsdila( zti, ztj ) )**2   & 
    257          !                                  + (                           fsdiph( zti, ztj ) )**2  ) 
    258          !     e1u(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphiu(ji,jj) ) * fsdila( zui, zuj ) )**2   & 
    259          !                                  + (                           fsdiph( zui, zuj ) )**2  ) 
    260          !     e1v(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphiv(ji,jj) ) * fsdila( zvi, zvj ) )**2   & 
    261          !                                  + (                           fsdiph( zvi, zvj ) )**2  ) 
    262          !     e1f(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphif(ji,jj) ) * fsdila( zfi, zfj ) )**2   & 
    263          !                                  + (                           fsdiph( zfi, zfj ) )**2  ) 
    264          ! 
    265          !     e2t(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphit(ji,jj) ) * fsdjla( zti, ztj ) )**2   & 
    266          !                                  + (                           fsdjph( zti, ztj ) )**2  ) 
    267          !     e2u(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphiu(ji,jj) ) * fsdjla( zui, zuj ) )**2   & 
    268          !                                  + (                           fsdjph( zui, zuj ) )**2  ) 
    269          !     e2v(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphiv(ji,jj) ) * fsdjla( zvi, zvj ) )**2   & 
    270          !                                  + (                           fsdjph( zvi, zvj ) )**2  ) 
    271          !     e2f(ji,jj) = ra * rad * SQRT(  ( cos( rad*gphif(ji,jj) ) * fsdjla( zfi, zfj ) )**2   & 
    272          !                                  + (                           fsdjph( zfi, zfj ) )**2  ) 
    273  
    274  
    275       CASE ( 1 )                     ! geographical mesh on the sphere with regular grid-spacing 
    276  
     145         ! 
     146      CASE ( 1 )                     !==  geographical mesh on the sphere with regular (in degree) grid-spacing  ==! 
     147         ! 
    277148         IF(lwp) WRITE(numout,*) 
    278149         IF(lwp) WRITE(numout,*) '          geographical mesh on the sphere with regular grid-spacing' 
    279150         IF(lwp) WRITE(numout,*) '          given by ppe1_deg and ppe2_deg'  
    280  
     151         ! 
    281152         DO jj = 1, jpj 
    282153            DO ji = 1, jpi 
    283                zti = FLOAT( ji - 1 + nimpp - 1 )         ;   ztj = FLOAT( jj - 1 + njmpp - 1 ) 
    284                zui = FLOAT( ji - 1 + nimpp - 1 ) + 0.5   ;   zuj = FLOAT( jj - 1 + njmpp - 1 ) 
    285                zvi = FLOAT( ji - 1 + nimpp - 1 )         ;   zvj = FLOAT( jj - 1 + njmpp - 1 ) + 0.5 
    286                zfi = FLOAT( ji - 1 + nimpp - 1 ) + 0.5   ;   zfj = FLOAT( jj - 1 + njmpp - 1 ) + 0.5 
     154               zti = REAL( ji - 1 + nimpp - 1 )         ;   ztj = REAL( jj - 1 + njmpp - 1 ) 
     155               zui = REAL( ji - 1 + nimpp - 1 ) + 0.5   ;   zuj = REAL( jj - 1 + njmpp - 1 ) 
     156               zvi = REAL( ji - 1 + nimpp - 1 )         ;   zvj = REAL( jj - 1 + njmpp - 1 ) + 0.5 
     157               zfi = REAL( ji - 1 + nimpp - 1 ) + 0.5   ;   zfj = REAL( jj - 1 + njmpp - 1 ) + 0.5 
    287158         ! Longitude 
    288159               glamt(ji,jj) = ppglam0 + ppe1_deg * zti 
     
    307178            END DO 
    308179         END DO 
    309  
    310  
    311       CASE ( 2:3 )                   ! f- or beta-plane with regular grid-spacing 
    312  
     180         ! 
     181      CASE ( 2:3 )                   !==  f- or beta-plane with regular grid-spacing  ==! 
     182         ! 
    313183         IF(lwp) WRITE(numout,*) 
    314184         IF(lwp) WRITE(numout,*) '          f- or beta-plane with regular grid-spacing' 
    315185         IF(lwp) WRITE(numout,*) '          given by ppe1_m and ppe2_m'  
    316  
     186         ! 
    317187         ! Position coordinates (in kilometers) 
    318188         !                          ========== 
    319          glam0 = 0.e0 
     189         glam0 = 0._wp 
    320190         gphi0 = - ppe2_m * 1.e-3 
    321           
     191         ! 
    322192#if defined key_agrif  
    323193         IF ( cp_cfg == 'eel' .AND. jp_cfg == 6 ) THEN    ! for EEL6 configuration only 
     
    332202         DO jj = 1, jpj 
    333203            DO ji = 1, jpi 
    334                glamt(ji,jj) = glam0 + ppe1_m * 1.e-3 * ( FLOAT( ji - 1 + nimpp - 1 )       ) 
    335                glamu(ji,jj) = glam0 + ppe1_m * 1.e-3 * ( FLOAT( ji - 1 + nimpp - 1 ) + 0.5 ) 
     204               glamt(ji,jj) = glam0 + ppe1_m * 1.e-3 * ( REAL( ji - 1 + nimpp - 1 )       ) 
     205               glamu(ji,jj) = glam0 + ppe1_m * 1.e-3 * ( REAL( ji - 1 + nimpp - 1 ) + 0.5 ) 
    336206               glamv(ji,jj) = glamt(ji,jj) 
    337207               glamf(ji,jj) = glamu(ji,jj) 
    338     
    339                gphit(ji,jj) = gphi0 + ppe2_m * 1.e-3 * ( FLOAT( jj - 1 + njmpp - 1 )       ) 
     208               ! 
     209               gphit(ji,jj) = gphi0 + ppe2_m * 1.e-3 * ( REAL( jj - 1 + njmpp - 1 )       ) 
    340210               gphiu(ji,jj) = gphit(ji,jj) 
    341                gphiv(ji,jj) = gphi0 + ppe2_m * 1.e-3 * ( FLOAT( jj - 1 + njmpp - 1 ) + 0.5 ) 
     211               gphiv(ji,jj) = gphi0 + ppe2_m * 1.e-3 * ( REAL( jj - 1 + njmpp - 1 ) + 0.5 ) 
    342212               gphif(ji,jj) = gphiv(ji,jj) 
    343213            END DO 
    344214         END DO 
    345  
     215         ! 
    346216         ! Horizontal scale factors (in meters) 
    347217         !                              ====== 
     
    350220         e1v(:,:) = ppe1_m      ;      e2v(:,:) = ppe2_m 
    351221         e1f(:,:) = ppe1_m      ;      e2f(:,:) = ppe2_m 
    352  
    353       CASE ( 4 )                     ! geographical mesh on the sphere, isotropic MERCATOR type 
    354  
     222         ! 
     223      CASE ( 4 )                     !==  geographical mesh on the sphere, isotropic MERCATOR type  ==! 
     224         ! 
    355225         IF(lwp) WRITE(numout,*) 
    356226         IF(lwp) WRITE(numout,*) '          geographical mesh on the sphere, MERCATOR type' 
    357227         IF(lwp) WRITE(numout,*) '          longitudinal/latitudinal spacing given by ppe1_deg' 
    358228         IF ( ppgphi0 == -90 ) CALL ctl_stop( ' Mercator grid cannot start at south pole !!!! ' ) 
    359  
     229         ! 
    360230         !  Find index corresponding to the equator, given the grid spacing e1_deg 
    361231         !  and the (approximate) southern latitude ppgphi0. 
     
    365235         ijeq = ABS( 180./rpi * LOG( COS( zarg ) / SIN( zarg ) ) / ppe1_deg ) 
    366236         IF(  ppgphi0 > 0 )  ijeq = -ijeq 
    367  
     237         ! 
    368238         IF(lwp) WRITE(numout,*) '          Index of the equator on the MERCATOR grid:', ijeq 
    369  
     239         ! 
    370240         DO jj = 1, jpj 
    371241            DO ji = 1, jpi 
    372                zti = FLOAT( ji - 1 + nimpp - 1 )         ;   ztj = FLOAT( jj - ijeq + njmpp - 1 ) 
    373                zui = FLOAT( ji - 1 + nimpp - 1 ) + 0.5   ;   zuj = FLOAT( jj - ijeq + njmpp - 1 ) 
    374                zvi = FLOAT( ji - 1 + nimpp - 1 )         ;   zvj = FLOAT( jj - ijeq + njmpp - 1 ) + 0.5 
    375                zfi = FLOAT( ji - 1 + nimpp - 1 ) + 0.5   ;   zfj = FLOAT( jj - ijeq + njmpp - 1 ) + 0.5 
     242               zti = REAL( ji - 1 + nimpp - 1 )         ;   ztj = REAL( jj - ijeq + njmpp - 1 ) 
     243               zui = REAL( ji - 1 + nimpp - 1 ) + 0.5   ;   zuj = REAL( jj - ijeq + njmpp - 1 ) 
     244               zvi = REAL( ji - 1 + nimpp - 1 )         ;   zvj = REAL( jj - ijeq + njmpp - 1 ) + 0.5 
     245               zfi = REAL( ji - 1 + nimpp - 1 ) + 0.5   ;   zfj = REAL( jj - ijeq + njmpp - 1 ) + 0.5 
    376246         ! Longitude 
    377247               glamt(ji,jj) = ppglam0 + ppe1_deg * zti 
     
    396266            END DO 
    397267         END DO 
    398  
    399       CASE ( 5 )                   ! beta-plane with regular grid-spacing and rotated domain (GYRE configuration) 
    400  
     268         ! 
     269      CASE ( 5 )                   !==  beta-plane with regular grid-spacing and rotated domain ==! (GYRE configuration) 
     270         ! 
    401271         IF(lwp) WRITE(numout,*) 
    402272         IF(lwp) WRITE(numout,*) '          beta-plane with regular grid-spacing and rotated domain (GYRE configuration)' 
    403273         IF(lwp) WRITE(numout,*) '          given by ppe1_m and ppe2_m' 
    404  
     274         ! 
    405275         ! Position coordinates (in kilometers) 
    406276         !                          ========== 
    407  
     277         ! 
    408278         ! angle 45deg and ze1=106.e+3 / jp_cfg forced -> zlam1 = -85deg, zphi1 = 29degN 
    409          zlam1 = -85 
    410          zphi1 = 29 
     279         zlam1 = -85._wp 
     280         zphi1 =  29._wp 
    411281         ! resolution in meters 
    412          ze1 = 106000. / FLOAT(jp_cfg)             
     282         ze1 = 106000. / REAL( jp_cfg , wp )             
    413283         ! benchmark: forced the resolution to be about 100 km 
    414          IF( nbench /= 0 )   ze1 = 106000.e0      
    415          zsin_alpha = - SQRT( 2. ) / 2. 
    416          zcos_alpha =   SQRT( 2. ) / 2. 
     284         IF( nbench /= 0 )   ze1 = 106000._wp      
     285         zsin_alpha = - SQRT( 2._wp ) * 0.5_wp 
     286         zcos_alpha =   SQRT( 2._wp ) * 0.5_wp 
    417287         ze1deg = ze1 / (ra * rad) 
    418          IF( nbench /= 0 )   ze1deg = ze1deg / FLOAT(jp_cfg)        ! benchmark: keep the lat/+lon 
    419          !                                                          ! at the right jp_cfg resolution 
    420          glam0 = zlam1 + zcos_alpha * ze1deg * FLOAT( jpjglo-2 ) 
    421          gphi0 = zphi1 + zsin_alpha * ze1deg * FLOAT( jpjglo-2 ) 
    422  
     288         IF( nbench /= 0 )   ze1deg = ze1deg / REAL( jp_cfg , wp )   ! benchmark: keep the lat/+lon 
     289         !                                                           ! at the right jp_cfg resolution 
     290         glam0 = zlam1 + zcos_alpha * ze1deg * REAL( jpjglo-2 , wp ) 
     291         gphi0 = zphi1 + zsin_alpha * ze1deg * REAL( jpjglo-2 , wp ) 
     292         ! 
    423293         IF( nprint==1 .AND. lwp )   THEN 
    424294            WRITE(numout,*) '          ze1', ze1, 'cosalpha', zcos_alpha, 'sinalpha', zsin_alpha 
    425295            WRITE(numout,*) '          ze1deg', ze1deg, 'glam0', glam0, 'gphi0', gphi0 
    426296         ENDIF 
    427  
     297         ! 
    428298         DO jj = 1, jpj 
    429            DO ji = 1, jpi 
    430              zim1 = FLOAT( ji + nimpp - 1 ) - 1.   ;   zim05 = FLOAT( ji + nimpp - 1 ) - 1.5 
    431              zjm1 = FLOAT( jj + njmpp - 1 ) - 1.   ;   zjm05 = FLOAT( jj + njmpp - 1 ) - 1.5 
    432  
    433              glamf(ji,jj) = glam0 + zim1  * ze1deg * zcos_alpha + zjm1  * ze1deg * zsin_alpha 
    434              gphif(ji,jj) = gphi0 - zim1  * ze1deg * zsin_alpha + zjm1  * ze1deg * zcos_alpha 
    435  
    436              glamt(ji,jj) = glam0 + zim05 * ze1deg * zcos_alpha + zjm05 * ze1deg * zsin_alpha 
    437              gphit(ji,jj) = gphi0 - zim05 * ze1deg * zsin_alpha + zjm05 * ze1deg * zcos_alpha 
    438  
    439              glamu(ji,jj) = glam0 + zim1  * ze1deg * zcos_alpha + zjm05 * ze1deg * zsin_alpha 
    440              gphiu(ji,jj) = gphi0 - zim1  * ze1deg * zsin_alpha + zjm05 * ze1deg * zcos_alpha 
    441  
    442              glamv(ji,jj) = glam0 + zim05 * ze1deg * zcos_alpha + zjm1  * ze1deg * zsin_alpha 
    443              gphiv(ji,jj) = gphi0 - zim05 * ze1deg * zsin_alpha + zjm1  * ze1deg * zcos_alpha 
    444            END DO 
    445           END DO 
    446  
     299            DO ji = 1, jpi 
     300               zim1 = REAL( ji + nimpp - 1 ) - 1.   ;   zim05 = REAL( ji + nimpp - 1 ) - 1.5 
     301               zjm1 = REAL( jj + njmpp - 1 ) - 1.   ;   zjm05 = REAL( jj + njmpp - 1 ) - 1.5 
     302               ! 
     303               glamf(ji,jj) = glam0 + zim1  * ze1deg * zcos_alpha + zjm1  * ze1deg * zsin_alpha 
     304               gphif(ji,jj) = gphi0 - zim1  * ze1deg * zsin_alpha + zjm1  * ze1deg * zcos_alpha 
     305               ! 
     306               glamt(ji,jj) = glam0 + zim05 * ze1deg * zcos_alpha + zjm05 * ze1deg * zsin_alpha 
     307               gphit(ji,jj) = gphi0 - zim05 * ze1deg * zsin_alpha + zjm05 * ze1deg * zcos_alpha 
     308               ! 
     309               glamu(ji,jj) = glam0 + zim1  * ze1deg * zcos_alpha + zjm05 * ze1deg * zsin_alpha 
     310               gphiu(ji,jj) = gphi0 - zim1  * ze1deg * zsin_alpha + zjm05 * ze1deg * zcos_alpha 
     311               ! 
     312               glamv(ji,jj) = glam0 + zim05 * ze1deg * zcos_alpha + zjm1  * ze1deg * zsin_alpha 
     313               gphiv(ji,jj) = gphi0 - zim05 * ze1deg * zsin_alpha + zjm1  * ze1deg * zcos_alpha 
     314            END DO 
     315         END DO 
     316         ! 
    447317         ! Horizontal scale factors (in meters) 
    448318         !                              ====== 
     
    451321         e1v(:,:) =  ze1     ;      e2v(:,:) = ze1 
    452322         e1f(:,:) =  ze1     ;      e2f(:,:) = ze1 
    453  
     323         ! 
    454324      CASE DEFAULT 
    455325         WRITE(ctmp1,*) '          bad flag value for jphgr_msh = ', jphgr_msh 
    456326         CALL ctl_stop( ctmp1 ) 
    457  
     327         ! 
    458328      END SELECT 
    459329       
    460       ! T-cell surface 
    461       ! -------------- 
    462       e1e2t(:,:) = e1t(:,:) * e2t(:,:) 
    463      
    464       ! Useful shortcuts (JC: note the duplicated e2e2t array ! Need some cleaning) 
    465       ! --------------------------------------------------------------------------- 
    466       e12t    (:,:) = e1t(:,:) * e2t(:,:) 
    467       e12u    (:,:) = e1u(:,:) * e2u(:,:) 
    468       e12v    (:,:) = e1v(:,:) * e2v(:,:) 
    469       e12f    (:,:) = e1f(:,:) * e2f(:,:) 
    470       r1_e12t (:,:) = 1._wp    / e12t(:,:) 
    471       r1_e12u (:,:) = 1._wp    / e12u(:,:) 
    472       r1_e12v (:,:) = 1._wp    / e12v(:,:) 
    473       r1_e12f (:,:) = 1._wp    / e12f(:,:) 
    474       re2u_e1u(:,:) = e2u(:,:) / e1u(:,:) 
    475       re1v_e2v(:,:) = e1v(:,:) / e2v(:,:) 
    476       r1_e1t  (:,:) = 1._wp    / e1t(:,:) 
    477       r1_e1u  (:,:) = 1._wp    / e1u(:,:) 
    478       r1_e1v  (:,:) = 1._wp    / e1v(:,:) 
    479       r1_e1f  (:,:) = 1._wp    / e1f(:,:) 
    480       r1_e2t  (:,:) = 1._wp    / e2t(:,:) 
    481       r1_e2u  (:,:) = 1._wp    / e2u(:,:) 
    482       r1_e2v  (:,:) = 1._wp    / e2v(:,:) 
    483       r1_e2f  (:,:) = 1._wp    / e2f(:,:) 
    484  
    485       ! Control printing : Grid informations (if not restart) 
    486       ! ---------------- 
    487  
    488       IF( lwp .AND. .NOT.ln_rstart ) THEN 
     330      ! associated horizontal metrics 
     331      ! ----------------------------- 
     332      ! 
     333      r1_e1t(:,:) = 1._wp / e1t(:,:)   ;   r1_e2t (:,:) = 1._wp / e2t(:,:) 
     334      r1_e1u(:,:) = 1._wp / e1u(:,:)   ;   r1_e2u (:,:) = 1._wp / e2u(:,:) 
     335      r1_e1v(:,:) = 1._wp / e1v(:,:)   ;   r1_e2v (:,:) = 1._wp / e2v(:,:) 
     336      r1_e1f(:,:) = 1._wp / e1f(:,:)   ;   r1_e2f (:,:) = 1._wp / e2f(:,:) 
     337      ! 
     338      e1e2t (:,:) = e1t(:,:) * e2t(:,:)   ;   r1_e1e2t(:,:) = 1._wp / e1e2t(:,:) 
     339      e1e2f (:,:) = e1f(:,:) * e2f(:,:)   ;   r1_e1e2f(:,:) = 1._wp / e1e2f(:,:) 
     340      IF( jphgr_msh /= 0 ) THEN               ! e1e2u and e1e2v have not been set: compute them 
     341         e1e2u (:,:) = e1u(:,:) * e2u(:,:)    
     342         e1e2v (:,:) = e1v(:,:) * e2v(:,:)  
     343      ENDIF 
     344      r1_e1e2u(:,:) = 1._wp / e1e2u(:,:)     ! compute their invert in both cases 
     345      r1_e1e2v(:,:) = 1._wp / e1e2v(:,:) 
     346      !    
     347      e2_e1u(:,:) = e2u(:,:) / e1u(:,:) 
     348      e1_e2v(:,:) = e1v(:,:) / e2v(:,:) 
     349 
     350      IF( lwp .AND. .NOT.ln_rstart ) THEN      ! Control print : Grid informations (if not restart) 
    489351         WRITE(numout,*) 
    490352         WRITE(numout,*) '          longitude and e1 scale factors' 
     
    4963589300     FORMAT( 1x, i4, f8.2,1x, f8.2,1x, f8.2,1x, f8.2, 1x,    & 
    497359            f19.10, 1x, f19.10, 1x, f19.10, 1x, f19.10 ) 
    498           
     360            ! 
    499361         WRITE(numout,*) 
    500362         WRITE(numout,*) '          latitude and e2 scale factors' 
     
    506368      ENDIF 
    507369 
    508        
    509       IF( nprint == 1 .AND. lwp ) THEN 
    510          WRITE(numout,*) '          e1u e2u ' 
    511          CALL prihre( e1u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    512          CALL prihre( e2u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    513          WRITE(numout,*) '          e1v e2v  ' 
    514          CALL prihre( e1v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    515          CALL prihre( e2v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    516          WRITE(numout,*) '          e1f e2f  ' 
    517          CALL prihre( e1f,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    518          CALL prihre( e2f,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout ) 
    519       ENDIF 
    520  
    521370 
    522371      ! ================= ! 
     
    525374 
    526375      SELECT CASE( jphgr_msh )   ! type of horizontal mesh 
    527  
     376      ! 
    528377      CASE ( 0, 1, 4 )               ! mesh on the sphere 
    529  
     378         ! 
    530379         ff(:,:) = 2. * omega * SIN( rad * gphif(:,:) )  
    531  
     380         ! 
    532381      CASE ( 2 )                     ! f-plane at ppgphi0  
    533  
     382         ! 
    534383         ff(:,:) = 2. * omega * SIN( rad * ppgphi0 ) 
    535  
     384         ! 
    536385         IF(lwp) WRITE(numout,*) '          f-plane: Coriolis parameter = constant = ', ff(1,1) 
    537  
     386         ! 
    538387      CASE ( 3 )                     ! beta-plane 
    539  
     388         ! 
    540389         zbeta   = 2. * omega * COS( rad * ppgphi0 ) / ra                       ! beta at latitude ppgphi0 
    541          zphi0   = ppgphi0 - FLOAT( jpjglo/2) * ppe2_m / ( ra * rad )           ! latitude of the first row F-points 
    542           
     390         zphi0   = ppgphi0 - REAL( jpjglo/2) * ppe2_m / ( ra * rad )           ! latitude of the first row F-points 
     391         ! 
    543392#if defined key_agrif 
    544393         IF ( cp_cfg == 'eel' .AND. jp_cfg == 6 ) THEN    ! for EEL6 configuration only 
    545394            IF( .NOT. Agrif_Root() ) THEN 
    546               zphi0 = ppgphi0 - FLOAT( Agrif_Parent(jpjglo)/2)*Agrif_Parent(ppe2_m)   &  
    547                     &           / (ra * rad) 
     395              zphi0 = ppgphi0 - REAL( Agrif_Parent(jpjglo)/2)*Agrif_Parent(ppe2_m) / (ra * rad) 
    548396            ENDIF 
    549397         ENDIF 
    550398#endif          
    551399         zf0     = 2. * omega * SIN( rad * zphi0 )                              ! compute f0 1st point south 
    552  
     400         ! 
    553401         ff(:,:) = ( zf0  + zbeta * gphif(:,:) * 1.e+3 )                        ! f = f0 +beta* y ( y=0 at south) 
    554           
     402         ! 
    555403         IF(lwp) THEN 
    556404            WRITE(numout,*)  
     
    565413            IF(lwp) WRITE(numout,*) '          Coriolis parameter varies globally from ', zminff,' to ', zmaxff 
    566414         END IF 
    567  
     415         ! 
    568416      CASE ( 5 )                     ! beta-plane and rotated domain (gyre configuration) 
    569  
     417         ! 
    570418         zbeta = 2. * omega * COS( rad * ppgphi0 ) / ra                     ! beta at latitude ppgphi0 
    571          zphi0 = 15.e0                                                      ! latitude of the first row F-points 
     419         zphi0 = 15._wp                                                     ! latitude of the first row F-points 
    572420         zf0   = 2. * omega * SIN( rad * zphi0 )                            ! compute f0 1st point south 
    573  
     421         ! 
    574422         ff(:,:) = ( zf0 + zbeta * ABS( gphif(:,:) - zphi0 ) * rad * ra )   ! f = f0 +beta* y ( y=0 at south) 
    575  
     423         ! 
    576424         IF(lwp) THEN 
    577425            WRITE(numout,*)  
     
    579427            WRITE(numout,*) '          Coriolis parameter varies in this processor from ', ff(nldi,nldj),' to ', ff(nldi,nlej) 
    580428         ENDIF 
    581  
     429         ! 
    582430         IF( lk_mpp ) THEN  
    583431            zminff=ff(nldi,nldj) 
     
    587435            IF(lwp) WRITE(numout,*) '          Coriolis parameter varies globally from ', zminff,' to ', zmaxff 
    588436         END IF 
    589  
     437         ! 
    590438      END SELECT 
    591439 
     
    596444 
    597445      IF( nperio == 2 ) THEN 
    598          znorme = SQRT( SUM( gphiu(:,2) * gphiu(:,2) ) ) / FLOAT( jpi ) 
     446         znorme = SQRT( SUM( gphiu(:,2) * gphiu(:,2) ) ) / REAL( jpi ) 
    599447         IF( znorme > 1.e-13 ) CALL ctl_stop( ' ===>>>> : symmetrical condition: rerun with good equator line' ) 
    600448      ENDIF 
     
    605453 
    606454 
    607    SUBROUTINE hgr_read 
     455   SUBROUTINE hgr_read( ke1e2u_v ) 
    608456      !!--------------------------------------------------------------------- 
    609457      !!              ***  ROUTINE hgr_read  *** 
    610458      !! 
    611       !! ** Purpose :   Read a coordinate file in NetCDF format  
    612       !! 
    613       !! ** Method  :   The mesh file has been defined trough a analytical  
    614       !!      or semi-analytical method. It is read in a NetCDF file.  
    615       !!      
     459      !! ** Purpose :   Read a coordinate file in NetCDF format using IOM 
     460      !! 
    616461      !!---------------------------------------------------------------------- 
    617462      USE iom 
    618  
     463      !! 
     464      INTEGER, INTENT( inout ) ::   ke1e2u_v   ! fag: e1e2u & e1e2v read in coordinate file (=1) or not (=0) 
     465      ! 
    619466      INTEGER ::   inum   ! temporary logical unit 
    620467      !!---------------------------------------------------------------------- 
    621  
     468      ! 
    622469      IF(lwp) THEN 
    623470         WRITE(numout,*) 
     
    625472         WRITE(numout,*) '~~~~~~~~      jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpk = ', jpk 
    626473      ENDIF 
    627        
     474      ! 
    628475      CALL iom_open( 'coordinates', inum ) 
    629        
     476      ! 
    630477      CALL iom_get( inum, jpdom_data, 'glamt', glamt, lrowattr=ln_use_jattr ) 
    631478      CALL iom_get( inum, jpdom_data, 'glamu', glamu, lrowattr=ln_use_jattr ) 
    632479      CALL iom_get( inum, jpdom_data, 'glamv', glamv, lrowattr=ln_use_jattr ) 
    633480      CALL iom_get( inum, jpdom_data, 'glamf', glamf, lrowattr=ln_use_jattr ) 
    634        
     481      ! 
    635482      CALL iom_get( inum, jpdom_data, 'gphit', gphit, lrowattr=ln_use_jattr ) 
    636483      CALL iom_get( inum, jpdom_data, 'gphiu', gphiu, lrowattr=ln_use_jattr ) 
    637484      CALL iom_get( inum, jpdom_data, 'gphiv', gphiv, lrowattr=ln_use_jattr ) 
    638485      CALL iom_get( inum, jpdom_data, 'gphif', gphif, lrowattr=ln_use_jattr ) 
    639        
    640       CALL iom_get( inum, jpdom_data, 'e1t', e1t, lrowattr=ln_use_jattr ) 
    641       CALL iom_get( inum, jpdom_data, 'e1u', e1u, lrowattr=ln_use_jattr ) 
    642       CALL iom_get( inum, jpdom_data, 'e1v', e1v, lrowattr=ln_use_jattr ) 
    643       CALL iom_get( inum, jpdom_data, 'e1f', e1f, lrowattr=ln_use_jattr ) 
    644        
    645       CALL iom_get( inum, jpdom_data, 'e2t', e2t, lrowattr=ln_use_jattr ) 
    646       CALL iom_get( inum, jpdom_data, 'e2u', e2u, lrowattr=ln_use_jattr ) 
    647       CALL iom_get( inum, jpdom_data, 'e2v', e2v, lrowattr=ln_use_jattr ) 
    648       CALL iom_get( inum, jpdom_data, 'e2f', e2f, lrowattr=ln_use_jattr ) 
    649        
     486      ! 
     487      CALL iom_get( inum, jpdom_data, 'e1t'  , e1t  , lrowattr=ln_use_jattr ) 
     488      CALL iom_get( inum, jpdom_data, 'e1u'  , e1u  , lrowattr=ln_use_jattr ) 
     489      CALL iom_get( inum, jpdom_data, 'e1v'  , e1v  , lrowattr=ln_use_jattr ) 
     490      CALL iom_get( inum, jpdom_data, 'e1f'  , e1f  , lrowattr=ln_use_jattr ) 
     491      ! 
     492      CALL iom_get( inum, jpdom_data, 'e2t'  , e2t  , lrowattr=ln_use_jattr ) 
     493      CALL iom_get( inum, jpdom_data, 'e2u'  , e2u  , lrowattr=ln_use_jattr ) 
     494      CALL iom_get( inum, jpdom_data, 'e2v'  , e2v  , lrowattr=ln_use_jattr ) 
     495      CALL iom_get( inum, jpdom_data, 'e2f'  , e2f  , lrowattr=ln_use_jattr ) 
     496      ! 
     497      IF( iom_varid( inum, 'e1e2u', ldstop = .FALSE. ) > 0 ) THEN 
     498         IF(lwp) WRITE(numout,*) 'hgr_read : e1e2u & e1e2v read in coordinates file' 
     499         CALL iom_get( inum, jpdom_data, 'e1e2u'  , e1e2u  , lrowattr=ln_use_jattr ) 
     500         CALL iom_get( inum, jpdom_data, 'e1e2v'  , e1e2v  , lrowattr=ln_use_jattr ) 
     501         ke1e2u_v = 1 
     502      ELSE 
     503         ke1e2u_v = 0 
     504      ENDIF 
     505      ! 
    650506      CALL iom_close( inum ) 
    651507       
     508!!gm   THIS is TO BE REMOVED !!!!!!! 
     509 
    652510! need to be define for the extended grid south of -80S 
    653511! some point are undefined but you need to have e1 and e2 .NE. 0 
     
    676534         e2f=1.0e2 
    677535      END WHERE 
     536!!gm end 
    678537        
    679538    END SUBROUTINE hgr_read 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90

    r5552 r5836  
    1717   !!             -   ! 2005-11  (V. Garnier) Surface pressure gradient organization 
    1818   !!            3.2  ! 2009-07  (R. Benshila) Suppression of rigid-lid option 
     19   !!            3.6  ! 2015-05  (P. Mathiot) ISF: add wmask,wumask and wvmask 
    1920   !!---------------------------------------------------------------------- 
    2021 
    2122   !!---------------------------------------------------------------------- 
    2223   !!   dom_msk        : compute land/ocean mask 
    23    !!   dom_msk_nsa    : update land/ocean mask when no-slip accurate option is used. 
    2424   !!---------------------------------------------------------------------- 
    2525   USE oce             ! ocean dynamics and tracers 
     
    3636 
    3737   PUBLIC   dom_msk         ! routine called by inidom.F90 
    38    PUBLIC   dom_msk_alloc   ! routine called by nemogcm.F90 
    3938 
    4039   !                            !!* Namelist namlbc : lateral boundary condition * 
     
    4241   LOGICAL, PUBLIC :: ln_vorlat  !  consistency of vorticity boundary condition  
    4342   !                                            with analytical eqs. 
    44  
    45  
    46    INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  icoord ! Workspace for dom_msk_nsa() 
    4743 
    4844   !! * Substitutions 
     
    5450   !!---------------------------------------------------------------------- 
    5551CONTAINS 
    56     
    57    INTEGER FUNCTION dom_msk_alloc() 
    58       !!--------------------------------------------------------------------- 
    59       !!                 ***  FUNCTION dom_msk_alloc  *** 
    60       !!--------------------------------------------------------------------- 
    61       dom_msk_alloc = 0 
    62 #if defined key_noslip_accurate 
    63       ALLOCATE(icoord(jpi*jpj*jpk,3), STAT=dom_msk_alloc) 
    64 #endif 
    65       IF( dom_msk_alloc /= 0 )   CALL ctl_warn('dom_msk_alloc: failed to allocate icoord array') 
    66       ! 
    67    END FUNCTION dom_msk_alloc 
    68  
    6952 
    7053   SUBROUTINE dom_msk 
     
    129112      !!               tmask_i  : interior ocean mask 
    130113      !!---------------------------------------------------------------------- 
    131       ! 
    132       INTEGER  ::   ji, jj, jk      ! dummy loop indices 
     114      INTEGER  ::   ji, jj, jk               ! dummy loop indices 
    133115      INTEGER  ::   iif, iil, ii0, ii1, ii   ! local integers 
    134116      INTEGER  ::   ijf, ijl, ij0, ij1       !   -       - 
     
    199181      END DO   
    200182 
    201 !!gm  ???? 
    202 #if defined key_zdfkpp 
    203       IF( cp_cfg == 'orca' ) THEN 
    204          IF( jp_cfg == 2 )   THEN       ! land point on Bab el Mandeb zonal section 
    205             ij0 =  87   ;   ij1 =  88 
    206             ii0 = 160   ;   ii1 = 161 
    207             tmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0._wp 
    208          ELSE 
    209             IF(lwp) WRITE(numout,*) 
    210             IF(lwp) WRITE(numout,cform_war) 
    211             IF(lwp) WRITE(numout,*) 
    212             IF(lwp) WRITE(numout,*)'          A mask must be applied on Bab el Mandeb strait' 
    213             IF(lwp) WRITE(numout,*)'          in case of ORCAs configurations' 
    214             IF(lwp) WRITE(numout,*)'          This is a problem which is not yet solved' 
    215             IF(lwp) WRITE(numout,*) 
    216          ENDIF 
    217       ENDIF 
    218 #endif 
    219 !!gm end 
    220  
    221183      ! Interior domain mask (used for global sum) 
    222184      ! -------------------- 
     
    284246      ! 3. Ocean/land mask at wu-, wv- and w points  
    285247      !---------------------------------------------- 
    286       wmask (:,:,1) = tmask(:,:,1) ! ???????? 
    287       wumask(:,:,1) = umask(:,:,1) ! ???????? 
    288       wvmask(:,:,1) = vmask(:,:,1) ! ???????? 
    289       DO jk=2,jpk 
    290          wmask (:,:,jk)=tmask(:,:,jk) * tmask(:,:,jk-1) 
    291          wumask(:,:,jk)=umask(:,:,jk) * umask(:,:,jk-1)    
    292          wvmask(:,:,jk)=vmask(:,:,jk) * vmask(:,:,jk-1) 
     248      wmask (:,:,1) = tmask(:,:,1)     ! surface 
     249      wumask(:,:,1) = umask(:,:,1) 
     250      wvmask(:,:,1) = vmask(:,:,1) 
     251      DO jk = 2, jpk                   ! interior values 
     252         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1) 
     253         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)    
     254         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1) 
    293255      END DO 
    294256 
     
    339301      ENDIF 
    340302 
    341  
    342       ! mask for second order calculation of vorticity 
    343       ! ---------------------------------------------- 
    344       CALL dom_msk_nsa 
    345  
    346        
    347303      ! Lateral boundary conditions on velocity (modify fmask) 
    348304      ! ---------------------------------------      
     
    377333      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA_R2 configuration 
    378334         !                                                 ! Increased lateral friction near of some straits 
    379          IF( nn_cla == 0 ) THEN 
    380             !                                ! Gibraltar strait  : partial slip (fmask=0.5) 
    381             ij0 = 101   ;   ij1 = 101 
    382             ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
    383             ij0 = 102   ;   ij1 = 102 
    384             ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
    385             ! 
    386             !                                ! Bab el Mandeb : partial slip (fmask=1) 
    387             ij0 =  87   ;   ij1 =  88 
    388             ii0 = 160   ;   ii1 = 160   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
    389             ij0 =  88   ;   ij1 =  88 
    390             ii0 = 159   ;   ii1 = 159   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
    391             ! 
    392          ENDIF 
     335         !                                ! Gibraltar strait  : partial slip (fmask=0.5) 
     336         ij0 = 101   ;   ij1 = 101 
     337         ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
     338         ij0 = 102   ;   ij1 = 102 
     339         ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
     340         ! 
     341         !                                ! Bab el Mandeb : partial slip (fmask=1) 
     342         ij0 =  87   ;   ij1 =  88 
     343         ii0 = 160   ;   ii1 = 160   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
     344         ij0 =  88   ;   ij1 =  88 
     345         ii0 = 159   ;   ii1 = 159   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
     346         ! 
    393347         !                                ! Danish straits  : strong slip (fmask > 2) 
    394348! We keep this as an example but it is instable in this case  
     
    500454      ! 
    501455   END SUBROUTINE dom_msk 
    502  
    503 #if defined key_noslip_accurate 
    504    !!---------------------------------------------------------------------- 
    505    !!   'key_noslip_accurate' :         accurate no-slip boundary condition 
    506    !!---------------------------------------------------------------------- 
    507     
    508    SUBROUTINE dom_msk_nsa 
    509       !!--------------------------------------------------------------------- 
    510       !!                 ***  ROUTINE dom_msk_nsa  *** 
    511       !!  
    512       !! ** Purpose : 
    513       !! 
    514       !! ** Method  : 
    515       !! 
    516       !! ** Action : 
    517       !!---------------------------------------------------------------------- 
    518       INTEGER  ::   ji, jj, jk, jl      ! dummy loop indices 
    519       INTEGER  ::   ine, inw, ins, inn, itest, ierror, iind, ijnd 
    520       REAL(wp) ::   zaa 
    521       !!--------------------------------------------------------------------- 
    522       ! 
    523       IF( nn_timing == 1 )  CALL timing_start('dom_msk_nsa') 
    524       ! 
    525       IF(lwp) WRITE(numout,*) 
    526       IF(lwp) WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition' 
    527       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~   using Schchepetkin and O Brian scheme' 
    528       IF( lk_mpp )   CALL ctl_stop( ' mpp version is not yet implemented' ) 
    529  
    530       ! mask for second order calculation of vorticity 
    531       ! ---------------------------------------------- 
    532       ! noslip boundary condition: fmask=1  at convex corner, store 
    533       ! index of straight coast meshes ( 'west', refering to a coast, 
    534       ! means west of the ocean, aso) 
    535        
    536       DO jk = 1, jpk 
    537          DO jl = 1, 4 
    538             npcoa(jl,jk) = 0 
    539             DO ji = 1, 2*(jpi+jpj) 
    540                nicoa(ji,jl,jk) = 0 
    541                njcoa(ji,jl,jk) = 0 
    542             END DO 
    543          END DO 
    544       END DO 
    545        
    546       IF( jperio == 2 ) THEN 
    547          WRITE(numout,*) ' ' 
    548          WRITE(numout,*) ' symetric boundary conditions need special' 
    549          WRITE(numout,*) ' treatment not implemented. we stop.' 
    550          STOP 
    551       ENDIF 
    552        
    553       ! convex corners 
    554        
    555       DO jk = 1, jpkm1 
    556          DO jj = 1, jpjm1 
    557             DO ji = 1, jpim1 
    558                zaa = tmask(ji  ,jj,jk) + tmask(ji  ,jj+1,jk)   & 
    559                   &+ tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk) 
    560                IF( ABS(zaa-3._wp) <= 0.1_wp )   fmask(ji,jj,jk) = 1._wp 
    561             END DO 
    562          END DO 
    563       END DO 
    564  
    565       ! north-south straight coast 
    566  
    567       DO jk = 1, jpkm1 
    568          inw = 0 
    569          ine = 0 
    570          DO jj = 2, jpjm1 
    571             DO ji = 2, jpim1 
    572                zaa = tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk) 
    573                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    574                   inw = inw + 1 
    575                   nicoa(inw,1,jk) = ji 
    576                   njcoa(inw,1,jk) = jj 
    577                   IF( nprint == 1 ) WRITE(numout,*) ' west  : ', jk, inw, ji, jj 
    578                ENDIF 
    579                zaa = tmask(ji,jj,jk) + tmask(ji,jj+1,jk) 
    580                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    581                   ine = ine + 1 
    582                   nicoa(ine,2,jk) = ji 
    583                   njcoa(ine,2,jk) = jj 
    584                   IF( nprint == 1 ) WRITE(numout,*) ' east  : ', jk, ine, ji, jj 
    585                ENDIF 
    586             END DO 
    587          END DO 
    588          npcoa(1,jk) = inw 
    589          npcoa(2,jk) = ine 
    590       END DO 
    591  
    592       ! west-east straight coast 
    593  
    594       DO jk = 1, jpkm1 
    595          ins = 0 
    596          inn = 0 
    597          DO jj = 2, jpjm1 
    598             DO ji =2, jpim1 
    599                zaa = tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) 
    600                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    601                   ins = ins + 1 
    602                   nicoa(ins,3,jk) = ji 
    603                   njcoa(ins,3,jk) = jj 
    604                   IF( nprint == 1 ) WRITE(numout,*) ' south : ', jk, ins, ji, jj 
    605                ENDIF 
    606                zaa = tmask(ji+1,jj,jk) + tmask(ji,jj,jk) 
    607                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    608                   inn = inn + 1 
    609                   nicoa(inn,4,jk) = ji 
    610                   njcoa(inn,4,jk) = jj 
    611                   IF( nprint == 1 ) WRITE(numout,*) ' north : ', jk, inn, ji, jj 
    612                ENDIF 
    613             END DO 
    614          END DO 
    615          npcoa(3,jk) = ins 
    616          npcoa(4,jk) = inn 
    617       END DO 
    618  
    619       itest = 2 * ( jpi + jpj ) 
    620       DO jk = 1, jpk 
    621          IF( npcoa(1,jk) > itest .OR. npcoa(2,jk) > itest .OR.   & 
    622              npcoa(3,jk) > itest .OR. npcoa(4,jk) > itest ) THEN 
    623              
    624             WRITE(ctmp1,*) ' level jk = ',jk 
    625             WRITE(ctmp2,*) ' straight coast index arraies are too small.:' 
    626             WRITE(ctmp3,*) ' npe, npw, nps, npn = ', npcoa(1,jk), npcoa(2,jk),   & 
    627                 &                                     npcoa(3,jk), npcoa(4,jk) 
    628             WRITE(ctmp4,*) ' 2*(jpi+jpj) = ',itest,'. we stop.' 
    629             CALL ctl_stop( ctmp1, ctmp2, ctmp3, ctmp4 ) 
    630         ENDIF 
    631       END DO 
    632  
    633       ierror = 0 
    634       iind = 0 
    635       ijnd = 0 
    636       IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 )   iind = 2 
    637       IF( nperio == 3 .OR. nperio == 4 .OR. nperio == 5 .OR. nperio == 6 )   ijnd = 2 
    638       DO jk = 1, jpk 
    639          DO jl = 1, npcoa(1,jk) 
    640             IF( nicoa(jl,1,jk)+3 > jpi+iind ) THEN 
    641                ierror = ierror+1 
    642                icoord(ierror,1) = nicoa(jl,1,jk) 
    643                icoord(ierror,2) = njcoa(jl,1,jk) 
    644                icoord(ierror,3) = jk 
    645             ENDIF 
    646          END DO 
    647          DO jl = 1, npcoa(2,jk) 
    648             IF(nicoa(jl,2,jk)-2 < 1-iind ) THEN 
    649                ierror = ierror + 1 
    650                icoord(ierror,1) = nicoa(jl,2,jk) 
    651                icoord(ierror,2) = njcoa(jl,2,jk) 
    652                icoord(ierror,3) = jk 
    653             ENDIF 
    654          END DO 
    655          DO jl = 1, npcoa(3,jk) 
    656             IF( njcoa(jl,3,jk)+3 > jpj+ijnd ) THEN 
    657                ierror = ierror + 1 
    658                icoord(ierror,1) = nicoa(jl,3,jk) 
    659                icoord(ierror,2) = njcoa(jl,3,jk) 
    660                icoord(ierror,3) = jk 
    661             ENDIF 
    662          END DO 
    663          DO jl = 1, npcoa(4,jk) 
    664             IF( njcoa(jl,4,jk)-2 < 1) THEN 
    665                ierror=ierror + 1 
    666                icoord(ierror,1) = nicoa(jl,4,jk) 
    667                icoord(ierror,2) = njcoa(jl,4,jk) 
    668                icoord(ierror,3) = jk 
    669             ENDIF 
    670          END DO 
    671       END DO 
    672        
    673       IF( ierror > 0 ) THEN 
    674          IF(lwp) WRITE(numout,*) 
    675          IF(lwp) WRITE(numout,*) '              Problem on lateral conditions' 
    676          IF(lwp) WRITE(numout,*) '                 Bad marking off at points:' 
    677          DO jl = 1, ierror 
    678             IF(lwp) WRITE(numout,*) 'Level:',icoord(jl,3),   & 
    679                &                  '  Point(',icoord(jl,1),',',icoord(jl,2),')' 
    680          END DO 
    681          CALL ctl_stop( 'We stop...' ) 
    682       ENDIF 
    683       ! 
    684       IF( nn_timing == 1 )  CALL timing_stop('dom_msk_nsa') 
    685       ! 
    686    END SUBROUTINE dom_msk_nsa 
    687  
    688 #else 
    689    !!---------------------------------------------------------------------- 
    690    !!   Default option :                                      Empty routine 
    691    !!---------------------------------------------------------------------- 
    692    SUBROUTINE dom_msk_nsa        
    693    END SUBROUTINE dom_msk_nsa 
    694 #endif 
    695456    
    696457   !!====================================================================== 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domngb.F90

    r3294 r5836  
    3737      !!                -> not good if located at too high latitude... 
    3838      !!---------------------------------------------------------------------- 
    39       ! 
    4039      REAL(wp)        , INTENT(in   ) ::   plon, plat   ! longitude,latitude of the point 
    4140      INTEGER         , INTENT(  out) ::   kii, kjj     ! i-,j-index of the closes grid point 
     
    4948      IF( nn_timing == 1 )  CALL timing_start('dom_ngb') 
    5049      ! 
    51       CALL wrk_alloc( jpi, jpj, zglam, zgphi, zmask, zdist ) 
     50      CALL wrk_alloc( jpi,jpj,  zglam, zgphi, zmask, zdist ) 
    5251      ! 
    5352      zmask(:,:) = 0._wp 
     
    7675      ENDIF 
    7776      ! 
    78       CALL wrk_dealloc( jpi, jpj, zglam, zgphi, zmask, zdist ) 
     77      CALL wrk_dealloc( jpi,jpj,  zglam, zgphi, zmask, zdist ) 
    7978      ! 
    8079      IF( nn_timing == 1 )  CALL timing_stop('dom_ngb') 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90

    r5506 r5836  
    1010   !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability 
    1111   !!---------------------------------------------------------------------- 
    12    !!   'key_vvl'                              variable volume 
    13    !!---------------------------------------------------------------------- 
     12 
    1413   !!---------------------------------------------------------------------- 
    1514   !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
     
    1918   !!   dom_vvl_rst      : read/write restart file 
    2019   !!   dom_vvl_ctl      : Check the vvl options 
    21    !!   dom_vvl_orca_fix : Recompute some area-weighted interpolations of vertical scale factors  
    22    !!                    : to account for manual changes to e[1,2][u,v] in some Straits  
    2320   !!---------------------------------------------------------------------- 
    24    !! * Modules used 
    2521   USE oce             ! ocean dynamics and tracers 
    2622   USE dom_oce         ! ocean space and time domain 
     
    3733   PRIVATE 
    3834 
    39    !! * Routine accessibility 
    4035   PUBLIC  dom_vvl_init       ! called by domain.F90 
    4136   PUBLIC  dom_vvl_sf_nxt     ! called by step.F90 
    4237   PUBLIC  dom_vvl_sf_swp     ! called by step.F90 
    4338   PUBLIC  dom_vvl_interpol   ! called by dynnxt.F90 
    44    PRIVATE dom_vvl_orca_fix   ! called by dom_vvl_interpol 
    45  
    46    !!* Namelist nam_vvl 
    47    LOGICAL , PUBLIC                                      :: ln_vvl_zstar = .FALSE.              ! zstar  vertical coordinate 
    48    LOGICAL , PUBLIC                                      :: ln_vvl_ztilde = .FALSE.             ! ztilde vertical coordinate 
    49    LOGICAL , PUBLIC                                      :: ln_vvl_layer = .FALSE.              ! level  vertical coordinate 
    50    LOGICAL , PUBLIC                                      :: ln_vvl_ztilde_as_zstar = .FALSE.    ! ztilde vertical coordinate 
    51    LOGICAL , PUBLIC                                      :: ln_vvl_zstar_at_eqtor = .FALSE.     ! ztilde vertical coordinate 
    52    LOGICAL , PUBLIC                                      :: ln_vvl_kepe = .FALSE.               ! kinetic/potential energy transfer 
    53    !                                                                                            ! conservation: not used yet 
    54    REAL(wp)                                              :: rn_ahe3                   ! thickness diffusion coefficient 
    55    REAL(wp)                                              :: rn_rst_e3t                ! ztilde to zstar restoration timescale [days] 
    56    REAL(wp)                                              :: rn_lf_cutoff              ! cutoff frequency for low-pass filter  [days] 
    57    REAL(wp)                                              :: rn_zdef_max               ! maximum fractional e3t deformation 
    58    LOGICAL , PUBLIC                                      :: ln_vvl_dbg = .FALSE.      ! debug control prints 
    59  
    60    !! * Module variables 
    61    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_td, vn_td                       ! thickness diffusion transport 
    62    REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdiv_lf                            ! low frequency part of hz divergence 
    63    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_b, tilde_e3t_n           ! baroclinic scale factors 
    64    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a, dtilde_e3t_a          ! baroclinic scale factors 
    65    REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_e3t                        ! retoring period for scale factors 
    66    REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_hdv                        ! retoring period for low freq. divergence 
     39 
     40   !                                                      !!* Namelist nam_vvl 
     41   LOGICAL , PUBLIC :: ln_vvl_zstar           = .FALSE.    ! zstar  vertical coordinate 
     42   LOGICAL , PUBLIC :: ln_vvl_ztilde          = .FALSE.    ! ztilde vertical coordinate 
     43   LOGICAL , PUBLIC :: ln_vvl_layer           = .FALSE.    ! level  vertical coordinate 
     44   LOGICAL , PUBLIC :: ln_vvl_ztilde_as_zstar = .FALSE.    ! ztilde vertical coordinate 
     45   LOGICAL , PUBLIC :: ln_vvl_zstar_at_eqtor  = .FALSE.    ! ztilde vertical coordinate 
     46   LOGICAL , PUBLIC :: ln_vvl_kepe            = .FALSE.    ! kinetic/potential energy transfer 
     47   !                                                       ! conservation: not used yet 
     48   REAL(wp)         :: rn_ahe3                             ! thickness diffusion coefficient 
     49   REAL(wp)         :: rn_rst_e3t                          ! ztilde to zstar restoration timescale [days] 
     50   REAL(wp)         :: rn_lf_cutoff                        ! cutoff frequency for low-pass filter  [days] 
     51   REAL(wp)         :: rn_zdef_max                         ! maximum fractional e3t deformation 
     52   LOGICAL , PUBLIC :: ln_vvl_dbg = .FALSE.                ! debug control prints 
     53 
     54   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_td, vn_td                ! thickness diffusion transport 
     55   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdiv_lf                     ! low frequency part of hz divergence 
     56   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_b, tilde_e3t_n    ! baroclinic scale factors 
     57   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a, dtilde_e3t_a   ! baroclinic scale factors 
     58   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_e3t                 ! retoring period for scale factors 
     59   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_hdv                 ! retoring period for low freq. divergence 
    6760 
    6861   !! * Substitutions 
     
    372365            DO jj = 1, jpjm1 
    373366               DO ji = 1, fs_jpim1   ! vector opt. 
    374                   un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * re2u_e1u(ji,jj) & 
    375                                   & * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
    376                   vn_td(ji,jj,jk) = rn_ahe3 * vmask(ji,jj,jk) * re1v_e2v(ji,jj) &  
    377                                   & * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji  ,jj+1,jk) ) 
     367                  un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * e2_e1u(ji,jj)          & 
     368                     &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
     369                  vn_td(ji,jj,jk) = rn_ahe3 * vmask(ji,jj,jk) * e1_e2v(ji,jj)          &  
     370                     &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji  ,jj+1,jk) ) 
    378371                  zwu(ji,jj) = zwu(ji,jj) + un_td(ji,jj,jk) 
    379372                  zwv(ji,jj) = zwv(ji,jj) + vn_td(ji,jj,jk) 
     
    394387                  tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + (   un_td(ji-1,jj  ,jk) - un_td(ji,jj,jk)    & 
    395388                     &                                          +     vn_td(ji  ,jj-1,jk) - vn_td(ji,jj,jk)    & 
    396                      &                                            ) * r1_e12t(ji,jj) 
     389                     &                                            ) * r1_e1e2t(ji,jj) 
    397390               END DO 
    398391            END DO 
     
    695688      !!                - vertical interpolation: simple averaging 
    696689      !!---------------------------------------------------------------------- 
    697       !! * Arguments 
    698       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in    ) ::  pe3_in     ! input e3 to be interpolated 
    699       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) ::  pe3_out    ! output interpolated e3 
    700       CHARACTER(LEN=*), INTENT( in )                    ::  pout       ! grid point of out scale factors 
    701       !                                                                !   =  'U', 'V', 'W, 'F', 'UW' or 'VW' 
    702       !! * Local declarations 
     690      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::  pe3_in    ! input e3 to be interpolated 
     691      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::  pe3_out   ! output interpolated e3 
     692      CHARACTER(LEN=*)                , INTENT(in   ) ::  pout      ! grid point of out scale factors 
     693      !                                                             !   =  'U', 'V', 'W, 'F', 'UW' or 'VW' 
     694      ! 
    703695      INTEGER ::   ji, jj, jk                                          ! dummy loop indices 
    704       LOGICAL ::   l_is_orca                                           ! local logical 
    705       !!---------------------------------------------------------------------- 
     696      !!---------------------------------------------------------------------- 
     697      ! 
    706698      IF( nn_timing == 1 )  CALL timing_start('dom_vvl_interpol') 
    707          ! 
    708       l_is_orca = .FALSE. 
    709       IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) l_is_orca = .TRUE.      ! ORCA R2 configuration - will need to correct some locations 
    710  
    711       SELECT CASE ( pout ) 
    712          !               ! ------------------------------------- ! 
    713       CASE( 'U' )        ! interpolation from T-point to U-point ! 
    714          !               ! ------------------------------------- ! 
    715          ! horizontal surface weighted interpolation 
     699      ! 
     700      SELECT CASE ( pout )    !==  type of interpolation  ==! 
     701         ! 
     702      CASE( 'U' )                   !* from T- to U-point : hor. surface weighted mean 
    716703         DO jk = 1, jpk 
    717704            DO jj = 1, jpjm1 
    718705               DO ji = 1, fs_jpim1   ! vector opt. 
    719                   pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * r1_e12u(ji,jj)                                   & 
    720                      &                       * (   e12t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) )     & 
    721                      &                           + e12t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) ) 
     706                  pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * r1_e1e2u(ji,jj)                                   & 
     707                     &                       * (   e1e2t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) )     & 
     708                     &                           + e1e2t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) ) 
    722709               END DO 
    723710            END DO 
    724711         END DO 
    725          ! 
    726          IF( l_is_orca ) CALL dom_vvl_orca_fix( pe3_in, pe3_out, pout ) 
    727          ! boundary conditions 
    728712         CALL lbc_lnk( pe3_out(:,:,:), 'U', 1._wp ) 
    729713         pe3_out(:,:,:) = pe3_out(:,:,:) + e3u_0(:,:,:) 
    730          !               ! ------------------------------------- ! 
    731       CASE( 'V' )        ! interpolation from T-point to V-point ! 
    732          !               ! ------------------------------------- ! 
    733          ! horizontal surface weighted interpolation 
     714         ! 
     715      CASE( 'V' )                   !* from T- to V-point : hor. surface weighted mean 
    734716         DO jk = 1, jpk 
    735717            DO jj = 1, jpjm1 
    736718               DO ji = 1, fs_jpim1   ! vector opt. 
    737                   pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) * r1_e12v(ji,jj)                                   & 
    738                      &                       * (   e12t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) )     & 
    739                      &                           + e12t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) ) 
     719                  pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) * r1_e1e2v(ji,jj)                                   & 
     720                     &                       * (   e1e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) )     & 
     721                     &                           + e1e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) ) 
    740722               END DO 
    741723            END DO 
    742724         END DO 
    743          ! 
    744          IF( l_is_orca ) CALL dom_vvl_orca_fix( pe3_in, pe3_out, pout ) 
    745          ! boundary conditions 
    746725         CALL lbc_lnk( pe3_out(:,:,:), 'V', 1._wp ) 
    747726         pe3_out(:,:,:) = pe3_out(:,:,:) + e3v_0(:,:,:) 
    748          !               ! ------------------------------------- ! 
    749       CASE( 'F' )        ! interpolation from U-point to F-point ! 
    750          !               ! ------------------------------------- ! 
    751          ! horizontal surface weighted interpolation 
     727         ! 
     728      CASE( 'F' )                   !* from U-point to F-point : hor. surface weighted mean 
    752729         DO jk = 1, jpk 
    753730            DO jj = 1, jpjm1 
    754731               DO ji = 1, fs_jpim1   ! vector opt. 
    755                   pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) * r1_e12f(ji,jj)               & 
    756                      &                       * (   e12u(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3u_0(ji,jj  ,jk) )     & 
    757                      &                           + e12u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) ) 
     732                  pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) * r1_e1e2f(ji,jj)               & 
     733                     &                       * (   e1e2u(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3u_0(ji,jj  ,jk) )     & 
     734                     &                           + e1e2u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) ) 
    758735               END DO 
    759736            END DO 
    760737         END DO 
    761          ! 
    762          IF( l_is_orca ) CALL dom_vvl_orca_fix( pe3_in, pe3_out, pout ) 
    763          ! boundary conditions 
    764738         CALL lbc_lnk( pe3_out(:,:,:), 'F', 1._wp ) 
    765739         pe3_out(:,:,:) = pe3_out(:,:,:) + e3f_0(:,:,:) 
    766          !               ! ------------------------------------- ! 
    767       CASE( 'W' )        ! interpolation from T-point to W-point ! 
    768          !               ! ------------------------------------- ! 
    769          ! vertical simple interpolation 
     740         ! 
     741      CASE( 'W' )                   !* from T- to W-point : vertical simple mean 
     742         ! 
    770743         pe3_out(:,:,1) = e3w_0(:,:,1) + pe3_in(:,:,1) - e3t_0(:,:,1) 
    771          ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 
     744         ! - ML - The use of mask in this formulea enables the special treatment of the last w-point without indirect adressing 
     745!!gm BUG? use here wmask in case of ISF ?  to be checked 
    772746         DO jk = 2, jpk 
    773747            pe3_out(:,:,jk) = e3w_0(:,:,jk) + ( 1.0_wp - 0.5_wp * tmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3t_0(:,:,jk-1) )   & 
    774748               &                            +            0.5_wp * tmask(:,:,jk)   * ( pe3_in(:,:,jk  ) - e3t_0(:,:,jk  ) ) 
    775749         END DO 
    776          !               ! -------------------------------------- ! 
    777       CASE( 'UW' )       ! interpolation from U-point to UW-point ! 
    778          !               ! -------------------------------------- ! 
    779          ! vertical simple interpolation 
     750         ! 
     751      CASE( 'UW' )                  !* from U- to UW-point : vertical simple mean 
     752         ! 
    780753         pe3_out(:,:,1) = e3uw_0(:,:,1) + pe3_in(:,:,1) - e3u_0(:,:,1) 
    781754         ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 
     755!!gm BUG? use here wumask in case of ISF ?  to be checked 
    782756         DO jk = 2, jpk 
    783757            pe3_out(:,:,jk) = e3uw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * umask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3u_0(:,:,jk-1) )   & 
    784758               &                             +            0.5_wp * umask(:,:,jk)   * ( pe3_in(:,:,jk  ) - e3u_0(:,:,jk  ) ) 
    785759         END DO 
    786          !               ! -------------------------------------- ! 
    787       CASE( 'VW' )       ! interpolation from V-point to VW-point ! 
    788          !               ! -------------------------------------- ! 
    789          ! vertical simple interpolation 
     760         ! 
     761      CASE( 'VW' )                  !* from V- to VW-point : vertical simple mean 
     762         ! 
    790763         pe3_out(:,:,1) = e3vw_0(:,:,1) + pe3_in(:,:,1) - e3v_0(:,:,1) 
    791764         ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 
     765!!gm BUG? use here wvmask in case of ISF ?  to be checked 
    792766         DO jk = 2, jpk 
    793767            pe3_out(:,:,jk) = e3vw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * vmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3v_0(:,:,jk-1) )   & 
     
    796770      END SELECT 
    797771      ! 
    798  
    799772      IF( nn_timing == 1 )  CALL timing_stop('dom_vvl_interpol') 
    800  
     773      ! 
    801774   END SUBROUTINE dom_vvl_interpol 
     775 
    802776 
    803777   SUBROUTINE dom_vvl_rst( kt, cdrw ) 
     
    813787      !!                they are set to 0. 
    814788      !!---------------------------------------------------------------------- 
    815       !! * Arguments 
    816789      INTEGER         , INTENT(in) ::   kt     ! ocean time-step 
    817790      CHARACTER(len=*), INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag 
    818       !! * Local declarations 
     791      ! 
    819792      INTEGER ::   jk 
    820793      INTEGER ::   id1, id2, id3, id4, id5     ! local integers 
     
    911884            END IF 
    912885         ENDIF 
    913  
     886         ! 
    914887      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN   ! Create restart file 
    915888         !                                   ! =================== 
     
    931904            CALL iom_rstput( kt, nitrst, numrow, 'hdiv_lf', hdiv_lf(:,:,:) ) 
    932905         ENDIF 
    933  
    934       ENDIF 
     906         ! 
     907      ENDIF 
     908      ! 
    935909      IF( nn_timing == 1 )  CALL timing_stop('dom_vvl_rst') 
    936  
     910      ! 
    937911   END SUBROUTINE dom_vvl_rst 
    938912 
     
    945919      !!                for vertical coordinate 
    946920      !!---------------------------------------------------------------------- 
    947       INTEGER ::   ioptio 
    948       INTEGER ::   ios 
    949  
     921      INTEGER ::   ioptio, ios 
     922      !! 
    950923      NAMELIST/nam_vvl/ ln_vvl_zstar, ln_vvl_ztilde, ln_vvl_layer, ln_vvl_ztilde_as_zstar, & 
    951                       & ln_vvl_zstar_at_eqtor      , rn_ahe3     , rn_rst_e3t            , & 
    952                       & rn_lf_cutoff               , rn_zdef_max , ln_vvl_dbg                ! not yet implemented: ln_vvl_kepe 
     924         &              ln_vvl_zstar_at_eqtor      , rn_ahe3     , rn_rst_e3t            , & 
     925         &              rn_lf_cutoff               , rn_zdef_max , ln_vvl_dbg                ! not yet implemented: ln_vvl_kepe 
    953926      !!----------------------------------------------------------------------  
    954  
     927      ! 
    955928      REWIND( numnam_ref )              ! Namelist nam_vvl in reference namelist :  
    956929      READ  ( numnam_ref, nam_vvl, IOSTAT = ios, ERR = 901) 
    957930901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_vvl in reference namelist', lwp ) 
    958  
     931      ! 
    959932      REWIND( numnam_cfg )              ! Namelist nam_vvl in configuration namelist : Parameters of the run 
    960933      READ  ( numnam_cfg, nam_vvl, IOSTAT = ios, ERR = 902 ) 
    961934902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_vvl in configuration namelist', lwp ) 
    962935      IF(lwm) WRITE ( numond, nam_vvl ) 
    963  
     936      ! 
    964937      IF(lwp) THEN                    ! Namelist print 
    965938         WRITE(numout,*) 
     
    994967         WRITE(numout,*) '                                         ln_vvl_dbg     = ', ln_vvl_dbg 
    995968      ENDIF 
    996  
     969      ! 
    997970      ioptio = 0                      ! Parameter control 
    998       IF( ln_vvl_ztilde_as_zstar ) ln_vvl_ztilde = .true. 
    999       IF( ln_vvl_zstar           )        ioptio = ioptio + 1 
    1000       IF( ln_vvl_ztilde          )        ioptio = ioptio + 1 
    1001       IF( ln_vvl_layer           )        ioptio = ioptio + 1 
    1002  
     971      IF( ln_vvl_ztilde_as_zstar )   ln_vvl_ztilde = .true. 
     972      IF( ln_vvl_zstar           )   ioptio = ioptio + 1 
     973      IF( ln_vvl_ztilde          )   ioptio = ioptio + 1 
     974      IF( ln_vvl_layer           )   ioptio = ioptio + 1 
     975      ! 
    1003976      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE vertical coordinate in namelist nam_vvl' ) 
    1004977      IF( .NOT. ln_vvl_zstar .AND. nn_isf .NE. 0) CALL ctl_stop( 'Only vvl_zstar has been tested with ice shelf cavity' ) 
    1005  
     978      ! 
    1006979      IF(lwp) THEN                   ! Print the choice 
    1007980         WRITE(numout,*) 
     
    1014987         ! IF( .NOT. ln_vvl_kepe ) WRITE(numout,*) '              kinetic to potential energy transfer : option not used' 
    1015988      ENDIF 
    1016  
     989      ! 
    1017990#if defined key_agrif 
    1018991      IF (.NOT.Agrif_Root()) CALL ctl_stop( 'AGRIF not implemented with non-linear free surface (key_vvl)' ) 
    1019992#endif 
    1020  
     993      ! 
    1021994   END SUBROUTINE dom_vvl_ctl 
    1022  
    1023    SUBROUTINE dom_vvl_orca_fix( pe3_in, pe3_out, pout ) 
    1024       !!--------------------------------------------------------------------- 
    1025       !!                   ***  ROUTINE dom_vvl_orca_fix  *** 
    1026       !!                      
    1027       !! ** Purpose :   Correct surface weighted, horizontally interpolated,  
    1028       !!                scale factors at locations that have been individually 
    1029       !!                modified in domhgr. Such modifications break the 
    1030       !!                relationship between e12t and e1u*e2u etc. 
    1031       !!                Recompute some scale factors ignoring the modified metric. 
    1032       !!---------------------------------------------------------------------- 
    1033       !! * Arguments 
    1034       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in    ) ::  pe3_in     ! input e3 to be interpolated 
    1035       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) ::  pe3_out    ! output interpolated e3 
    1036       CHARACTER(LEN=*), INTENT( in )                    ::  pout       ! grid point of out scale factors 
    1037       !                                                                !   =  'U', 'V', 'W, 'F', 'UW' or 'VW' 
    1038       !! * Local declarations 
    1039       INTEGER ::   ji, jj, jk                                          ! dummy loop indices 
    1040       INTEGER ::   ij0, ij1, ii0, ii1                                  ! dummy loop indices 
    1041       INTEGER ::   isrow                                               ! index for ORCA1 starting row 
    1042       !! acc 
    1043       !! Hmm with the time splitting these "fixes" seem to do more harm than good. Temporarily disabled for 
    1044       !! the ORCA2 tests (by changing jp_cfg test from 2 to 3) pending further investigations 
    1045       !!  
    1046       !                                                ! ===================== 
    1047       IF( cp_cfg == "orca" .AND. jp_cfg == 3 ) THEN    ! ORCA R2 configuration 
    1048          !                                             ! ===================== 
    1049       !! acc 
    1050          IF( nn_cla == 0 ) THEN 
    1051             ! 
    1052             ii0 = 139   ;   ii1 = 140        ! Gibraltar Strait (e2u was modified) 
    1053             ij0 = 102   ;   ij1 = 102 
    1054             DO jk = 1, jpkm1 
    1055                DO jj = mj0(ij0), mj1(ij1) 
    1056                   DO ji = mi0(ii0), mi1(ii1) 
    1057                      SELECT CASE ( pout ) 
    1058                      CASE( 'U' ) 
    1059                         pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1060                        &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1061                        &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1062                        &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1063                      CASE( 'F' ) 
    1064                         pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1065                        &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1066                        &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1067                        &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1068                      END SELECT 
    1069                   END DO 
    1070                END DO 
    1071             END DO 
    1072             ! 
    1073             ii0 = 160   ;   ii1 = 160        ! Bab el Mandeb (e2u and e1v were modified) 
    1074             ij0 =  88   ;   ij1 =  88 
    1075             DO jk = 1, jpkm1 
    1076                DO jj = mj0(ij0), mj1(ij1) 
    1077                   DO ji = mi0(ii0), mi1(ii1) 
    1078                      SELECT CASE ( pout ) 
    1079                      CASE( 'U' ) 
    1080                         pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1081                        &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1082                        &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1083                        &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1084                      CASE( 'V' ) 
    1085                         pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1086                        &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1087                        &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1088                        &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1089                      CASE( 'F' ) 
    1090                         pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1091                        &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1092                        &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1093                        &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1094                      END SELECT 
    1095                   END DO 
    1096                END DO 
    1097             END DO 
    1098          ENDIF 
    1099  
    1100          ii0 = 145   ;   ii1 = 146        ! Danish Straits (e2u was modified) 
    1101          ij0 = 116   ;   ij1 = 116 
    1102          DO jk = 1, jpkm1 
    1103             DO jj = mj0(ij0), mj1(ij1) 
    1104                DO ji = mi0(ii0), mi1(ii1) 
    1105                   SELECT CASE ( pout ) 
    1106                   CASE( 'U' ) 
    1107                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1108                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1109                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1110                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1111                   CASE( 'F' ) 
    1112                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1113                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1114                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1115                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1116                   END SELECT 
    1117                END DO 
    1118             END DO 
    1119          END DO 
    1120       ENDIF 
    1121       ! 
    1122          !                                             ! ===================== 
    1123       IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN    ! ORCA R1 configuration 
    1124          !                                             ! ===================== 
    1125          ! This dirty section will be suppressed by simplification process: 
    1126          ! all this will come back in input files 
    1127          ! Currently these hard-wired indices relate to configuration with 
    1128          ! extend grid (jpjglo=332) 
    1129          ! which had a grid-size of 362x292. 
    1130          isrow = 332 - jpjglo 
    1131          ! 
    1132          ii0 = 282           ;   ii1 = 283        ! Gibraltar Strait (e2u was modified) 
    1133          ij0 = 241 - isrow   ;   ij1 = 241 - isrow 
    1134          DO jk = 1, jpkm1 
    1135             DO jj = mj0(ij0), mj1(ij1) 
    1136                DO ji = mi0(ii0), mi1(ii1) 
    1137                   SELECT CASE ( pout ) 
    1138                   CASE( 'U' ) 
    1139                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1140                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1141                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1142                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1143                   CASE( 'F' ) 
    1144                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1145                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1146                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1147                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1148                   END SELECT 
    1149                END DO 
    1150             END DO 
    1151          END DO 
    1152          ! 
    1153          ii0 = 314           ;   ii1 = 315        ! Bhosporus Strait (e2u was modified) 
    1154          ij0 = 248 - isrow   ;   ij1 = 248 - isrow 
    1155          DO jk = 1, jpkm1 
    1156             DO jj = mj0(ij0), mj1(ij1) 
    1157                DO ji = mi0(ii0), mi1(ii1) 
    1158                   SELECT CASE ( pout ) 
    1159                   CASE( 'U' ) 
    1160                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        &   
    1161                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1162                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1163                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1164                   CASE( 'F' ) 
    1165                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    &   
    1166                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1167                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1168                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1169                   END SELECT 
    1170                END DO 
    1171             END DO 
    1172          END DO 
    1173          ! 
    1174          ii0 =  44           ;   ii1 =  44        ! Lombok Strait (e1v was modified) 
    1175          ij0 = 164 - isrow   ;   ij1 = 165 - isrow 
    1176          DO jk = 1, jpkm1 
    1177             DO jj = mj0(ij0), mj1(ij1) 
    1178                DO ji = mi0(ii0), mi1(ii1) 
    1179                   SELECT CASE ( pout ) 
    1180                   CASE( 'V' ) 
    1181                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1182                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1183                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1184                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1185                   END SELECT 
    1186                END DO 
    1187             END DO 
    1188          END DO 
    1189          ! 
    1190          ii0 =  48           ;   ii1 =  48        ! Sumba Strait (e1v was modified) [closed from bathy_11 on] 
    1191          ij0 = 164 - isrow   ;   ij1 = 165 - isrow 
    1192          DO jk = 1, jpkm1 
    1193             DO jj = mj0(ij0), mj1(ij1) 
    1194                DO ji = mi0(ii0), mi1(ii1) 
    1195                   SELECT CASE ( pout ) 
    1196                   CASE( 'V' ) 
    1197                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1198                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1199                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1200                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1201                   END SELECT 
    1202                END DO 
    1203             END DO 
    1204          END DO 
    1205          ! 
    1206          ii0 =  53          ;   ii1 =  53        ! Ombai Strait (e1v was modified) 
    1207          ij0 = 164 - isrow  ;   ij1 = 165  - isrow   
    1208          DO jk = 1, jpkm1 
    1209             DO jj = mj0(ij0), mj1(ij1) 
    1210                DO ji = mi0(ii0), mi1(ii1) 
    1211                   SELECT CASE ( pout ) 
    1212                   CASE( 'V' ) 
    1213                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1214                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1215                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1216                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1217                   END SELECT 
    1218                END DO 
    1219             END DO 
    1220          END DO 
    1221          ! 
    1222          ii0 =  56            ;   ii1 =  56        ! Timor Passage (e1v was modified) 
    1223          ij0 = 164 - isrow    ;   ij1 = 165  - isrow   
    1224          DO jk = 1, jpkm1 
    1225             DO jj = mj0(ij0), mj1(ij1) 
    1226                DO ji = mi0(ii0), mi1(ii1) 
    1227                   SELECT CASE ( pout ) 
    1228                   CASE( 'V' ) 
    1229                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1230                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1231                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1232                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1233                   END SELECT 
    1234                END DO 
    1235             END DO 
    1236          END DO 
    1237          ! 
    1238          ii0 =  55            ;   ii1 =  55        ! West Halmahera Strait (e1v was modified) 
    1239          ij0 = 181 - isrow    ;   ij1 = 182 - isrow   
    1240          DO jk = 1, jpkm1 
    1241             DO jj = mj0(ij0), mj1(ij1) 
    1242                DO ji = mi0(ii0), mi1(ii1) 
    1243                   SELECT CASE ( pout ) 
    1244                   CASE( 'V' ) 
    1245                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1246                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1247                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1248                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1249                   END SELECT 
    1250                END DO 
    1251             END DO 
    1252          END DO 
    1253          ! 
    1254          ii0 =  58            ;   ii1 =  58        ! East Halmahera Strait (e1v was modified) 
    1255          ij0 = 181 - isrow    ;   ij1 = 182 - isrow   
    1256          DO jk = 1, jpkm1 
    1257             DO jj = mj0(ij0), mj1(ij1) 
    1258                DO ji = mi0(ii0), mi1(ii1) 
    1259                   SELECT CASE ( pout ) 
    1260                   CASE( 'V' ) 
    1261                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1262                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1263                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1264                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1265                   END SELECT 
    1266                END DO 
    1267             END DO 
    1268          END DO 
    1269       ENDIF 
    1270          !                                             ! ===================== 
    1271       IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN   ! ORCA R05 configuration 
    1272          !                                             ! ===================== 
    1273          ! 
    1274          ii0 = 563   ;   ii1 = 564        ! Gibraltar Strait (e2u was modified) 
    1275          ij0 = 327   ;   ij1 = 327 
    1276          DO jk = 1, jpkm1 
    1277             DO jj = mj0(ij0), mj1(ij1) 
    1278                DO ji = mi0(ii0), mi1(ii1) 
    1279                   SELECT CASE ( pout ) 
    1280                   CASE( 'U' ) 
    1281                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1282                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1283                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1284                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1285                   CASE( 'F' ) 
    1286                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1287                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1288                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1289                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1290                   END SELECT 
    1291                END DO 
    1292             END DO 
    1293          END DO 
    1294          ! 
    1295          ii0 = 627   ;   ii1 = 628        ! Bosphorus Strait (e2u was modified) 
    1296          ij0 = 343   ;   ij1 = 343 
    1297          DO jk = 1, jpkm1 
    1298             DO jj = mj0(ij0), mj1(ij1) 
    1299                DO ji = mi0(ii0), mi1(ii1) 
    1300                   SELECT CASE ( pout ) 
    1301                   CASE( 'U' ) 
    1302                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        &   
    1303                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1304                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1305                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1306                   CASE( 'F' ) 
    1307                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    &   
    1308                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1309                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1310                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1311                   END SELECT 
    1312                END DO 
    1313             END DO 
    1314          END DO 
    1315          ! 
    1316          ii0 =  93   ;   ii1 =  94        ! Sumba Strait (e2u was modified) 
    1317          ij0 = 232   ;   ij1 = 232 
    1318          DO jk = 1, jpkm1 
    1319             DO jj = mj0(ij0), mj1(ij1) 
    1320                DO ji = mi0(ii0), mi1(ii1) 
    1321                   SELECT CASE ( pout ) 
    1322                   CASE( 'U' ) 
    1323                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1324                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1325                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1326                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1327                   CASE( 'F' ) 
    1328                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1329                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1330                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1331                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1332                   END SELECT 
    1333                END DO 
    1334             END DO 
    1335          END DO 
    1336          ! 
    1337          ii0 = 103   ;   ii1 = 103        ! Ombai Strait (e2u was modified) 
    1338          ij0 = 232   ;   ij1 = 232 
    1339          DO jk = 1, jpkm1 
    1340             DO jj = mj0(ij0), mj1(ij1) 
    1341                DO ji = mi0(ii0), mi1(ii1) 
    1342                   SELECT CASE ( pout ) 
    1343                   CASE( 'U' ) 
    1344                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1345                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1346                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1347                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1348                   CASE( 'F' ) 
    1349                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1350                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1351                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1352                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1353                   END SELECT 
    1354                END DO 
    1355             END DO 
    1356          END DO 
    1357          ! 
    1358          ii0 =  15   ;   ii1 =  15        ! Palk Strait (e2u was modified) 
    1359          ij0 = 270   ;   ij1 = 270 
    1360          DO jk = 1, jpkm1 
    1361             DO jj = mj0(ij0), mj1(ij1) 
    1362                DO ji = mi0(ii0), mi1(ii1) 
    1363                   SELECT CASE ( pout ) 
    1364                   CASE( 'U' ) 
    1365                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk)                                        & 
    1366                     &                    * (   e1t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) ) & 
    1367                     &                    +     e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) & 
    1368                     &                      ) / e1u(ji,jj)   +   e3u_0(ji,jj,jk) 
    1369                   CASE( 'F' ) 
    1370                      pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk)                    & 
    1371                     &                    * (   e1u(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3u_0(ji  ,jj,jk) ) & 
    1372                     &                    +     e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) & 
    1373                     &                      ) / e1f(ji,jj)   +   e3f_0(ji,jj,jk) 
    1374                   END SELECT 
    1375                END DO 
    1376             END DO 
    1377          END DO 
    1378          ! 
    1379          ii0 =  87   ;   ii1 =  87        ! Lombok Strait (e1v was modified) 
    1380          ij0 = 232   ;   ij1 = 233 
    1381          DO jk = 1, jpkm1 
    1382             DO jj = mj0(ij0), mj1(ij1) 
    1383                DO ji = mi0(ii0), mi1(ii1) 
    1384                   SELECT CASE ( pout ) 
    1385                   CASE( 'V' ) 
    1386                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1387                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1388                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1389                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1390                   END SELECT 
    1391                END DO 
    1392             END DO 
    1393          END DO 
    1394          ! 
    1395          ii0 = 662   ;   ii1 = 662        ! Bab el Mandeb (e1v was modified) 
    1396          ij0 = 276   ;   ij1 = 276 
    1397          DO jk = 1, jpkm1 
    1398             DO jj = mj0(ij0), mj1(ij1) 
    1399                DO ji = mi0(ii0), mi1(ii1) 
    1400                   SELECT CASE ( pout ) 
    1401                   CASE( 'V' ) 
    1402                      pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk)                                        & 
    1403                     &                    * (   e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) ) & 
    1404                     &                    +     e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) & 
    1405                     &                      ) / e2v(ji,jj)   +   e3v_0(ji,jj,jk) 
    1406                   END SELECT 
    1407                END DO 
    1408             END DO 
    1409          END DO 
    1410       ENDIF 
    1411    END SUBROUTINE dom_vvl_orca_fix 
    1412995 
    1413996   !!====================================================================== 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90

    r5603 r5836  
    77   !!            8.1  ! 1999-11  (M. Imbard)  NetCDF FORMAT with IOIPSL 
    88   !!   NEMO     1.0  ! 2002-08  (G. Madec)  F90 and several file 
    9    !!            3.0  ! 2008-01  (S. Masson) add dom_uniq  
    10    !!            4.0  ! 2011-01  (A. R. Porter, STFC Daresbury) dynamical allocation 
     9   !!            3.0  ! 2008-01  (S. Masson)  add dom_uniq  
    1110   !!---------------------------------------------------------------------- 
    1211 
    1312   !!---------------------------------------------------------------------- 
    1413   !!   dom_wri        : create and write mesh and mask file(s) 
    15    !!   dom_uniq       : 
     14   !!   dom_uniq       : identify unique point of a grid (TUVF) 
    1615   !!---------------------------------------------------------------------- 
    1716   USE dom_oce         ! ocean space and time domain 
     
    2625   PRIVATE 
    2726 
    28    PUBLIC dom_wri        ! routine called by inidom.F90 
    29  
     27   PUBLIC   dom_wri              ! routine called by inidom.F90 
     28   PUBLIC   dom_wri_coordinate   ! routine called by domhgr.F90 
    3029   !! * Substitutions 
    3130#  include "vectopt_loop_substitute.h90" 
     
    3635   !!---------------------------------------------------------------------- 
    3736CONTAINS 
     37 
     38   SUBROUTINE dom_wri_coordinate 
     39      !!---------------------------------------------------------------------- 
     40      !!                  ***  ROUTINE dom_wri_coordinate  *** 
     41      !!                    
     42      !! ** Purpose :   Create the NetCDF file which contains all the 
     43      !!              standard coordinate information plus the surface, 
     44      !!              e1e2u and e1e2v. By doing so, those surface will 
     45      !!              not be changed by the reduction of e1u or e2v scale  
     46      !!              factors in some straits.  
     47      !!                 NB: call just after the read of standard coordinate 
     48      !!              and the reduction of scale factors in some straits 
     49      !! 
     50      !! ** output file :   coordinate_e1e2u_v.nc 
     51      !!---------------------------------------------------------------------- 
     52      INTEGER           ::   inum0    ! temprary units for 'coordinate_e1e2u_v.nc' file 
     53      CHARACTER(len=21) ::   clnam0   ! filename (mesh and mask informations) 
     54      !                                   !  workspaces 
     55      REAL(wp), POINTER, DIMENSION(:,:  ) :: zprt, zprw  
     56      REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepu, zdepv 
     57      !!---------------------------------------------------------------------- 
     58      ! 
     59      IF( nn_timing == 1 )  CALL timing_start('dom_wri_coordinate') 
     60      ! 
     61      IF(lwp) WRITE(numout,*) 
     62      IF(lwp) WRITE(numout,*) 'dom_wri_coordinate : create NetCDF coordinate file' 
     63      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~~' 
     64       
     65      clnam0 = 'coordinate_e1e2u_v'  ! filename (mesh and mask informations) 
     66       
     67      !  create 'coordinate_e1e2u_v.nc' file 
     68      ! ============================ 
     69      ! 
     70      CALL iom_open( TRIM(clnam0), inum0, ldwrt = .TRUE., kiolib = jprstlib ) 
     71      ! 
     72      !                                                         ! horizontal mesh (inum3) 
     73      CALL iom_rstput( 0, 0, inum0, 'glamt', glamt, ktype = jp_r4 )     !    ! latitude 
     74      CALL iom_rstput( 0, 0, inum0, 'glamu', glamu, ktype = jp_r4 ) 
     75      CALL iom_rstput( 0, 0, inum0, 'glamv', glamv, ktype = jp_r4 ) 
     76      CALL iom_rstput( 0, 0, inum0, 'glamf', glamf, ktype = jp_r4 ) 
     77       
     78      CALL iom_rstput( 0, 0, inum0, 'gphit', gphit, ktype = jp_r4 )     !    ! longitude 
     79      CALL iom_rstput( 0, 0, inum0, 'gphiu', gphiu, ktype = jp_r4 ) 
     80      CALL iom_rstput( 0, 0, inum0, 'gphiv', gphiv, ktype = jp_r4 ) 
     81      CALL iom_rstput( 0, 0, inum0, 'gphif', gphif, ktype = jp_r4 ) 
     82       
     83      CALL iom_rstput( 0, 0, inum0, 'e1t', e1t, ktype = jp_r8 )         !    ! e1 scale factors 
     84      CALL iom_rstput( 0, 0, inum0, 'e1u', e1u, ktype = jp_r8 ) 
     85      CALL iom_rstput( 0, 0, inum0, 'e1v', e1v, ktype = jp_r8 ) 
     86      CALL iom_rstput( 0, 0, inum0, 'e1f', e1f, ktype = jp_r8 ) 
     87       
     88      CALL iom_rstput( 0, 0, inum0, 'e2t', e2t, ktype = jp_r8 )         !    ! e2 scale factors 
     89      CALL iom_rstput( 0, 0, inum0, 'e2u', e2u, ktype = jp_r8 ) 
     90      CALL iom_rstput( 0, 0, inum0, 'e2v', e2v, ktype = jp_r8 ) 
     91      CALL iom_rstput( 0, 0, inum0, 'e2f', e2f, ktype = jp_r8 ) 
     92       
     93      CALL iom_rstput( 0, 0, inum0, 'e1e2u', e1e2u, ktype = jp_r8 ) 
     94      CALL iom_rstput( 0, 0, inum0, 'e1e2v', e1e2v, ktype = jp_r8 ) 
     95 
     96      CALL iom_close( inum0 ) 
     97      ! 
     98      IF( nn_timing == 1 )  CALL timing_stop('dom_wri_coordinate') 
     99      ! 
     100   END SUBROUTINE dom_wri_coordinate 
     101 
    38102 
    39103   SUBROUTINE dom_wri 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r5656 r5836  
    501501            IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN    ! ORCA R2 configuration 
    502502               !                                             ! ===================== 
    503                IF( nn_cla == 0 ) THEN 
    504                   ii0 = 140   ;   ii1 = 140                  ! Gibraltar Strait open  
    505                   ij0 = 102   ;   ij1 = 102                  ! (Thomson, Ocean Modelling, 1995) 
    506                   DO ji = mi0(ii0), mi1(ii1) 
    507                      DO jj = mj0(ij0), mj1(ij1) 
    508                         mbathy(ji,jj) = 15 
    509                      END DO 
     503               ! 
     504               ii0 = 140   ;   ii1 = 140                  ! Gibraltar Strait open  
     505               ij0 = 102   ;   ij1 = 102                  ! (Thomson, Ocean Modelling, 1995) 
     506               DO ji = mi0(ii0), mi1(ii1) 
     507                  DO jj = mj0(ij0), mj1(ij1) 
     508                     mbathy(ji,jj) = 15 
    510509                  END DO 
    511                   IF(lwp) WRITE(numout,*) 
    512                   IF(lwp) WRITE(numout,*) '      orca_r2: Gibraltar strait open at i=',ii0,' j=',ij0 
    513                   ! 
    514                   ii0 = 160   ;   ii1 = 160                  ! Bab el mandeb Strait open 
    515                   ij0 = 88    ;   ij1 = 88                   ! (Thomson, Ocean Modelling, 1995) 
    516                   DO ji = mi0(ii0), mi1(ii1) 
    517                      DO jj = mj0(ij0), mj1(ij1) 
    518                         mbathy(ji,jj) = 12 
    519                      END DO 
     510               END DO 
     511               IF(lwp) WRITE(numout,*) 
     512               IF(lwp) WRITE(numout,*) '      orca_r2: Gibraltar strait open at i=',ii0,' j=',ij0 
     513               ! 
     514               ii0 = 160   ;   ii1 = 160                  ! Bab el mandeb Strait open 
     515               ij0 = 88    ;   ij1 = 88                   ! (Thomson, Ocean Modelling, 1995) 
     516               DO ji = mi0(ii0), mi1(ii1) 
     517                  DO jj = mj0(ij0), mj1(ij1) 
     518                     mbathy(ji,jj) = 12 
    520519                  END DO 
    521                   IF(lwp) WRITE(numout,*) 
    522                   IF(lwp) WRITE(numout,*) '      orca_r2: Bab el Mandeb strait open at i=',ii0,' j=',ij0 
    523                ENDIF 
     520               END DO 
     521               IF(lwp) WRITE(numout,*) 
     522               IF(lwp) WRITE(numout,*) '      orca_r2: Bab el Mandeb strait open at i=',ii0,' j=',ij0 
    524523               ! 
    525524            ENDIF 
     
    545544            !        
    546545            IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN    ! ORCA R2 configuration 
    547                ! 
    548               IF( nn_cla == 0 ) THEN 
    549                  ii0 = 140   ;   ii1 = 140                   ! Gibraltar Strait open  
    550                  ij0 = 102   ;   ij1 = 102                   ! (Thomson, Ocean Modelling, 1995) 
    551                  DO ji = mi0(ii0), mi1(ii1) 
    552                     DO jj = mj0(ij0), mj1(ij1) 
    553                        bathy(ji,jj) = 284._wp 
    554                     END DO 
     546            ! 
     547              ii0 = 140   ;   ii1 = 140                   ! Gibraltar Strait open  
     548              ij0 = 102   ;   ij1 = 102                   ! (Thomson, Ocean Modelling, 1995) 
     549              DO ji = mi0(ii0), mi1(ii1) 
     550                 DO jj = mj0(ij0), mj1(ij1) 
     551                    bathy(ji,jj) = 284._wp 
    555552                 END DO 
    556                  IF(lwp) WRITE(numout,*)      
    557                  IF(lwp) WRITE(numout,*) '      orca_r2: Gibraltar strait open at i=',ii0,' j=',ij0 
    558                  ! 
    559                  ii0 = 160   ;   ii1 = 160                   ! Bab el mandeb Strait open 
    560                  ij0 = 88    ;   ij1 = 88                    ! (Thomson, Ocean Modelling, 1995) 
    561                  DO ji = mi0(ii0), mi1(ii1) 
    562                     DO jj = mj0(ij0), mj1(ij1) 
    563                        bathy(ji,jj) = 137._wp 
    564                     END DO 
     553               END DO 
     554              IF(lwp) WRITE(numout,*)      
     555              IF(lwp) WRITE(numout,*) '      orca_r2: Gibraltar strait open at i=',ii0,' j=',ij0 
     556              ! 
     557              ii0 = 160   ;   ii1 = 160                   ! Bab el mandeb Strait open 
     558              ij0 = 88    ;   ij1 = 88                    ! (Thomson, Ocean Modelling, 1995) 
     559               DO ji = mi0(ii0), mi1(ii1) 
     560                 DO jj = mj0(ij0), mj1(ij1) 
     561                    bathy(ji,jj) = 137._wp 
    565562                 END DO 
    566                  IF(lwp) WRITE(numout,*) 
    567                  IF(lwp) WRITE(numout,*) '             orca_r2: Bab el Mandeb strait open at i=',ii0,' j=',ij0 
    568               ENDIF 
     563              END DO 
     564              IF(lwp) WRITE(numout,*) 
     565               IF(lwp) WRITE(numout,*) '             orca_r2: Bab el Mandeb strait open at i=',ii0,' j=',ij0 
    569566              ! 
    570567           ENDIF 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dtatsd.F90

    r5215 r5836  
    174174            END DO 
    175175         END DO 
    176          IF( nn_cla == 1 ) THEN                          ! Cross Land advection 
    177             il0 = 138   ;   il1 = 138                          ! set T & S profile at Gibraltar Strait 
    178             ij0 = 101   ;   ij1 = 102 
    179             ii0 = 139   ;   ii1 = 139 
    180             DO jl = mi0(il0), mi1(il1) 
    181                DO jj = mj0(ij0), mj1(ij1) 
    182                   DO ji = mi0(ii0), mi1(ii1) 
    183                      sf_tsd(jp_tem)%fnow(ji,jj,:) = sf_tsd(jp_tem)%fnow(jl,jj,:) 
    184                      sf_tsd(jp_sal)%fnow(ji,jj,:) = sf_tsd(jp_sal)%fnow(jl,jj,:) 
    185                   END DO 
    186                END DO 
    187             END DO 
    188             il0 = 164   ;   il1 = 164                          ! set T & S profile at Bab el Mandeb Strait 
    189             ij0 =  87   ;   ij1 =  88 
    190             ii0 = 161   ;   ii1 = 163 
    191             DO jl = mi0(il0), mi1(il1) 
    192                DO jj = mj0(ij0), mj1(ij1) 
    193                   DO ji = mi0(ii0), mi1(ii1) 
    194                      sf_tsd(jp_tem)%fnow(ji,jj,:) = sf_tsd(jp_tem)%fnow(jl,jj,:) 
    195                      sf_tsd(jp_sal)%fnow(ji,jj,:) = sf_tsd(jp_sal)%fnow(jl,jj,:) 
    196                   END DO 
    197                END DO 
    198             END DO 
    199          ELSE                                            ! No Cross Land advection 
    200             ij0 =  87   ;   ij1 =  96                          ! Reduced temperature in Red Sea 
    201             ii0 = 148   ;   ii1 = 160 
    202             sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ,  4:10 ) = 7.0_wp 
    203             sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 11:13 ) = 6.5_wp 
    204             sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 14:20 ) = 6.0_wp 
    205          ENDIF 
     176         ij0 =  87   ;   ij1 =  96                          ! Reduced temperature in Red Sea 
     177         ii0 = 148   ;   ii1 = 160 
     178         sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ,  4:10 ) = 7.0_wp 
     179         sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 11:13 ) = 6.5_wp 
     180         sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 14:20 ) = 6.0_wp 
    206181      ENDIF 
    207182      ! 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r5332 r5836  
    2929   USE daymod          ! calendar 
    3030   USE eosbn2          ! eq. of state, Brunt Vaisala frequency (eos     routine) 
    31    USE ldftra_oce      ! ocean active tracers: lateral physics 
     31   USE ldftra          ! lateral physics: ocean active tracers 
    3232   USE zdf_oce         ! ocean vertical physics 
    3333   USE phycst          ! physical constants 
    3434   USE dtatsd          ! data temperature and salinity   (dta_tsd routine) 
    3535   USE dtauvd          ! data: U & V current             (dta_uvd routine) 
    36    USE zpshde          ! partial step: hor. derivative (zps_hde routine) 
    37    USE eosbn2          ! equation of state            (eos bn2 routine) 
    3836   USE domvvl          ! varying vertical mesh 
    3937   USE dynspg_oce      ! pressure gradient schemes 
     
    7674      ! 
    7775 
    78       IF(lwp) WRITE(numout,*) ' ' 
     76      IF(lwp) WRITE(numout,*) 
    7977      IF(lwp) WRITE(numout,*) 'istate_ini : Initialization of the dynamics and tracers' 
    8078      IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    8179 
    82       CALL dta_tsd_init                       ! Initialisation of T & S input data 
    83       IF( lk_c1d ) CALL dta_uvd_init          ! Initialization of U & V input data 
     80                     CALL dta_tsd_init        ! Initialisation of T & S input data 
     81      IF( lk_c1d )   CALL dta_uvd_init        ! Initialization of U & V input data 
    8482 
    8583      rhd  (:,:,:  ) = 0._wp   ;   rhop (:,:,:  ) = 0._wp      ! set one for all to 0 at level jpk 
     
    103101         ub   (:,:,:) = 0._wp   ;   un   (:,:,:) = 0._wp 
    104102         vb   (:,:,:) = 0._wp   ;   vn   (:,:,:) = 0._wp   
    105          rotb (:,:,:) = 0._wp   ;   rotn (:,:,:) = 0._wp 
    106          hdivb(:,:,:) = 0._wp   ;   hdivn(:,:,:) = 0._wp 
     103                                    hdivn(:,:,:) = 0._wp 
    107104         ! 
    108105         IF( cp_cfg == 'eel' ) THEN 
     
    119116            ENDIF 
    120117            IF ( ln_uvd_init .AND. lk_c1d ) THEN ! read 3D U and V data at nit000 
    121                CALL wrk_alloc( jpi, jpj, jpk, 2, zuvd ) 
     118               CALL wrk_alloc( jpi,jpj,jpk,2,  zuvd ) 
    122119               CALL dta_uvd( nit000, zuvd ) 
    123120               ub(:,:,:) = zuvd(:,:,:,1) ;  un(:,:,:) = ub(:,:,:) 
    124121               vb(:,:,:) = zuvd(:,:,:,2) ;  vn(:,:,:) = vb(:,:,:) 
    125                CALL wrk_dealloc( jpi, jpj, jpk, 2, zuvd ) 
     122               CALL wrk_dealloc( jpi,jpj,jpk,2,  zuvd ) 
    126123            ENDIF 
    127124         ENDIF 
    128          ! 
    129          CALL eos( tsb, rhd, rhop, gdept_0(:,:,:) )        ! before potential and in situ densities 
    130 #if ! defined key_c1d 
    131          IF( ln_zps .AND. .NOT. ln_isfcav)                                 & 
    132             &            CALL zps_hde    ( nit000, jpts, tsb, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    133             &                                            rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
    134          IF( ln_zps .AND.       ln_isfcav)                                 & 
    135             &            CALL zps_hde_isf( nit000, jpts, tsb, gtsu, gtsv,  &    ! Partial steps for top cell (ISF) 
    136             &                                            rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
    137             &                                     gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
    138 #endif 
    139125         !    
    140126         ! - ML - sshn could be modified by istate_eel, so that initialization of fse3t_b is done here 
     
    142128            DO jk = 1, jpk 
    143129               fse3t_b(:,:,jk) = fse3t_n(:,:,jk) 
    144             ENDDO 
     130            END DO 
    145131         ENDIF 
    146132         !  
     
    163149      ! 
    164150      ! 
    165       un_b(:,:) = 0._wp ; vn_b(:,:) = 0._wp 
    166       ub_b(:,:) = 0._wp ; vb_b(:,:) = 0._wp 
     151      un_b(:,:) = 0._wp   ;  vn_b(:,:) = 0._wp 
     152      ub_b(:,:) = 0._wp   ;  vb_b(:,:) = 0._wp 
    167153      ! 
    168154      DO jk = 1, jpkm1 
     
    202188      !! References :  Philander ??? 
    203189      !!---------------------------------------------------------------------- 
    204       INTEGER  :: ji, jj, jk 
    205       REAL(wp) ::   zsal = 35.50 
     190      INTEGER  ::   ji, jj, jk 
     191      REAL(wp) ::   zsal = 35.50_wp 
    206192      !!---------------------------------------------------------------------- 
    207193      ! 
     
    233219      !!                and relative vorticity fields 
    234220      !!---------------------------------------------------------------------- 
    235       USE divcur     ! hor. divergence & rel. vorticity      (div_cur routine) 
     221      USE divhor     ! hor. divergence      (div_hor routine) 
    236222      USE iom 
    237223      ! 
     
    282268            tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
    283269            ! 
    284             ! set the dynamics: U,V, hdiv, rot (and ssh if necessary) 
     270            ! set the dynamics: U,V, hdiv (and ssh if necessary) 
    285271            ! ---------------- 
    286272            ! Start EEL5 configuration with barotropic geostrophic velocities  
     
    318304            ENDIF 
    319305            ! 
    320             CALL div_cur( nit000 )                  ! horizontal divergence and relative vorticity (curl) 
     306!!gm  Check  here call to div_hor should not be necessary 
     307!!gm         div_hor call runoffs  not sure they are defined at that level 
     308            CALL div_hor( nit000 )                  ! horizontal divergence and relative vorticity (curl) 
    321309            ! N.B. the vertical velocity will be computed from the horizontal divergence field 
    322310            ! in istate by a call to wzv routine 
     
    371359      !! 
    372360      !! ** Method  : - set temprature field 
    373       !!              - set salinity field 
     361      !!              - set salinity   field 
    374362      !!---------------------------------------------------------------------- 
    375363      INTEGER :: ji, jj, jk  ! dummy loop indices 
     
    458446      !!---------------------------------------------------------------------- 
    459447      USE dynspg          ! surface pressure gradient             (dyn_spg routine) 
    460       USE divcur          ! hor. divergence & rel. vorticity      (div_cur routine) 
     448      USE divhor          ! hor. divergence                       (div_hor routine) 
    461449      USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    462450      ! 
     
    467455      !!---------------------------------------------------------------------- 
    468456      ! 
    469       CALL wrk_alloc( jpi, jpj, jpk, zprn) 
     457      CALL wrk_alloc( jpi,jpj,jpk,  zprn) 
    470458      ! 
    471459      IF(lwp) WRITE(numout,*)  
     
    544532      indic = 0 
    545533      CALL dyn_spg( nit000, indic )       ! surface pressure gradient 
    546  
     534      ! 
    547535      ! the new velocity is ua*rdt 
    548  
     536      ! 
    549537      CALL lbc_lnk( ua, 'U', -1. ) 
    550538      CALL lbc_lnk( va, 'V', -1. ) 
     
    556544      un(:,:,:) = ub(:,:,:) 
    557545      vn(:,:,:) = vb(:,:,:) 
    558         
    559       ! Compute the divergence and curl 
    560  
    561       CALL div_cur( nit000 )            ! now horizontal divergence and curl 
    562  
    563       hdivb(:,:,:) = hdivn(:,:,:)       ! set the before to the now value 
    564       rotb (:,:,:) = rotn (:,:,:)       ! set the before to the now value 
    565       ! 
    566       CALL wrk_dealloc( jpi, jpj, jpk, zprn) 
     546      ! 
     547!!gm  Check  here call to div_hor should not be necessary 
     548!!gm         div_hor call runoffs  not sure they are defined at that level 
     549      CALL div_hor( nit000 )            ! now horizontal divergence 
     550      ! 
     551      CALL wrk_dealloc( jpi,jpj,jpk,   zprn) 
    567552      ! 
    568553   END SUBROUTINE istate_uvg 
Note: See TracChangeset for help on using the changeset viewer.