Changeset 8143


Ignore:
Timestamp:
2017-06-06T15:55:44+02:00 (3 years ago)
Author:
gm
Message:

#1880 (HPC-09) - step-7: top/bottom drag computed at T-points, zdfbfr.F90 replaced by zdfdrg.F90 + changes in namelist

Location:
branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM
Files:
1 deleted
44 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/AMM12/EXP00/namelist_cfg

    r7990 r8143  
    174174   filtide      = 'bdydta/amm12_bdytide_'         !  file name root of tidal forcing files 
    175175/ 
    176 !----------------------------------------------------------------------- 
    177 &nambfr        !   bottom friction 
    178 !----------------------------------------------------------------------- 
    179    nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    180                            !                              = 2 : nonlinear friction 
    181    rn_bfri2    =    2.5e-3 !  bottom drag coefficient (non linear case) 
    182    rn_bfeb2    =    0.0e0  !  bottom turbulent kinetic energy background  (m2/s2) 
    183    ln_loglayer =    .true. !  loglayer bottom friction (only effect when nn_bfr = 2) 
    184    rn_bfrz0    =    0.003  !  bottom roughness (only effect when ln_loglayer = .true.) 
     176 
     177!----------------------------------------------------------------------- 
     178&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     179!----------------------------------------------------------------------- 
     180   ln_NONE    = .false.    !  free-slip       : Cd = 0                  (F => fill namdrg_bot 
     181   ln_lin     = .false.    !      linear  drag: Cd = Cd0 Uc0                   &   namdrg_top) 
     182   ln_non_lin = .false.    !  non-linear  drag: Cd = Cd0 |U| 
     183   ln_loglayer= .true.     !  logarithmic drag: Cd = vkarmn/log(z/z0) |U| 
     184   ! 
     185   ln_drgimp  = .true.     !  implicit top/bottom friction flag 
     186/ 
     187!----------------------------------------------------------------------- 
     188&namdrg_bot        !   BOTTOM friction                                   
     189!----------------------------------------------------------------------- 
     190   rn_Cd0     =  2.5e-3   !  drag coefficient [-] 
     191   rn_Uc0     =  0.4      !  ref. velocity [m/s] (linear drag=Cd0*Uc0)  
     192   rn_Cdmax   =  0.1      !  drag value maximum [-] (logarithmic drag) 
     193   rn_ke0     =  0.0e0    !  background kinetic energy  [m2/s2] (non-linear cases) 
     194   rn_z0      =  0.003    !  roughness [m]  (ln_loglayer=T) 
     195   ln_boost   = .false.   !  =T regional boost of Cd0 ; =F constant 
     196      rn_boost=  50.         !  local boost factor  [-] 
    185197/ 
    186198!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/C1D_PAPA/EXP00/namelist_cfg

    r7990 r8143  
    183183/ 
    184184!----------------------------------------------------------------------- 
    185 &nambfr        !   bottom friction 
    186 !----------------------------------------------------------------------- 
    187    nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
     185&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     186!----------------------------------------------------------------------- 
     187   ln_non_lin = .true.     !  non-linear  drag: Cd = Cd0 |U| 
    188188/ 
    189189!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/GYRE_BFM/EXP00/namelist_cfg

    r7990 r8143  
    124124/ 
    125125!----------------------------------------------------------------------- 
    126 &nambfr        !   bottom friction 
    127 !----------------------------------------------------------------------- 
    128    nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
     126&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     127!----------------------------------------------------------------------- 
     128   ln_non_lin = .true.     !  non-linear  drag: Cd = Cd0 |U| 
    129129/ 
    130130!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/namelist_cfg

    r7990 r8143  
    7878/ 
    7979!----------------------------------------------------------------------- 
    80 &nambfr        !   bottom friction 
    81 !----------------------------------------------------------------------- 
    82    nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
     80&namdrg        !   top/bottom friction 
     81!----------------------------------------------------------------------- 
     82   ln_non_lin = .true.     !  non-linear  drag: Cd = Cd0 |U| 
    8383/ 
    8484!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/ORCA2_LIM3_PISCES/EXP00/1_namelist_cfg

    r7990 r8143  
    8686/ 
    8787!----------------------------------------------------------------------- 
    88 &nambfr        !   bottom friction 
    89 !----------------------------------------------------------------------- 
     88&namdrg        !   bottom friction 
     89!----------------------------------------------------------------------- 
     90   ln_lin = .true.         !     linear  drag: Cd = Cd0 Uc0 
    9091/ 
    9192!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/ORCA2_LIM3_PISCES/EXP00/namelist_cfg

    r7990 r8143  
    9999/ 
    100100!----------------------------------------------------------------------- 
    101 &nambfr        !   bottom friction 
    102 !----------------------------------------------------------------------- 
     101&namdrg        !   top/bottom friction 
     102!----------------------------------------------------------------------- 
     103   ln_lin = .true.         !     linear  drag: Cd = Cd0 Uc0 
    103104/ 
    104105!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/SHARED/namelist_ref

    r8093 r8143  
    88!!                                    namsbc_apr, namsbc_ssr, namsbc_alb, namsbc_wave) 
    99!!              4 - lateral boundary (namlbc, namagrif, nambdy, nambdy_tide) 
    10 !!              5 - bottom  boundary (nambfr, nambbc, nambbl) 
     10!!              5 - bottom  boundary (namdrg, namdrg_top, namdrg_bot, nambbc, nambbl) 
    1111!!              6 - Tracer           (nameos, namtra_adv, namtra_ldf, namtra_ldfeiv, namtra_dmp) 
    1212!!              7 - dynamics         (namdyn_adv, namdyn_vor, namdyn_hpg, namdyn_spg, namdyn_ldf) 
     
    603603!!                ***  top/Bottom boundary condition  ***             !! 
    604604!!====================================================================== 
    605 !!   nambfr        bottom friction                                      (default: NONE) 
    606 !!   namtfr        top    friction                                      (default: NONE) 
     605!!   namdrg        top/bottom drag coefficient                          (default: NONE) 
     606!!   namdrg_top    top    friction                                      (ln_isfcav=T) 
     607!!   namdrg_bot    bottom friction                                       
    607608!!   nambbc        bottom temperature boundary condition                (default: NO) 
    608609!!   nambbl        bottom boundary layer scheme                         (default: NO) 
     
    610611! 
    611612!----------------------------------------------------------------------- 
    612 &nambfr        !   bottom friction                                      (default: linear) 
    613 !----------------------------------------------------------------------- 
    614    nn_bfr      =    1      !  type of top/bottom drag: free slip (=0), linear drag (=1),  
    615    !                       !  nonlinear drag (=2), nonlinear with logarithmic formulation (=3) 
    616    ln_bfrimp   = .true.    !  implicit bottom friction (requires ln_zdfexp = .false. if true) 
    617    ln_loglayer = .false.   !  logarithmic formulation (non linear case only) 
    618                            ! 
    619    rn_bfri1    =    4.e-4  !  bottom drag coefficient (linear case) 
    620    rn_bfri2    =    1.e-3  !  bottom drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
    621    rn_bfri2_max=    1.e-1  !  max. bottom drag coefficient (non linear case and ln_loglayer=T) 
    622    rn_bfeb2    =    2.5e-3 !  bottom turbulent kinetic energy background  (m2/s2) 
    623    rn_bfrz0    =    3.e-3  !  bottom roughness [m] if ln_loglayer=T 
    624    ln_bfr2d    = .false.   !  horizontal variation of the bottom friction coef (read a 2D mask file ) 
    625       rn_bfrien   =    50.    !  local boost factor 
    626    ! 
    627    rn_tfri1    =    4.e-4  !  top drag coefficient (linear case) 
    628    rn_tfri2    =    2.5e-3 !  top drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
    629    rn_tfri2_max=    1.e-1  !  max. top drag coefficient (non linear case and ln_loglayer=T) 
    630    rn_tfeb2    =    0.0    !  top turbulent kinetic energy background  (m2/s2) 
    631    rn_tfrz0    =    3.e-3  !  top roughness [m] if ln_loglayer=T 
    632    ln_tfr2d    = .false.   !  horizontal variation of the top friction coef (read a 2D mask file ) 
    633       rn_tfrien   =   50.     !  local boost factor 
    634 / 
    635  
    636 !----------------------------------------------------------------------- 
    637 &namdrg            !   top/bottom drag coeeficient                      (default: NO selection) 
     613&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
    638614!----------------------------------------------------------------------- 
    639615   ln_NONE    = .false.    !  free-slip       : Cd = 0                  (F => fill namdrg_bot 
     
    644620   ln_drgimp  = .true.     !  implicit top/bottom friction flag 
    645621/ 
    646  
    647622!----------------------------------------------------------------------- 
    648623&namdrg_top        !   TOP friction                                     (ln_isfcav=T) 
     
    656631      rn_boost=  50.          !  local boost factor  [-] 
    657632/ 
    658  
    659633!----------------------------------------------------------------------- 
    660634&namdrg_bot        !   BOTTOM friction                                   
     
    668642      rn_boost=  50.         !  local boost factor  [-] 
    669643/ 
    670  
    671644!----------------------------------------------------------------------- 
    672645&nambbc        !   bottom temperature boundary condition                (default: NO) 
     
    980953   ln_mxl0     = .true.    !  surface mixing length scale = F(wind stress) (T) or not (F) 
    981954   rn_mxl0     =   0.04    !  surface  buoyancy lenght scale minimum value 
    982 !!gm   ln_drg      = .false.   !  top/bottom friction added as boundary condition of TKE 
     955   ln_drg      = .false.   !  top/bottom friction added as boundary condition of TKE 
    983956   ln_lc       = .true.    !  Langmuir cell parameterisation (Axell 2002) 
    984    rn_lc       =   0.15    !  coef. associated to Langmuir cells 
     957      rn_lc       =   0.15    !  coef. associated to Langmuir cells 
    985958   nn_etau     =   1       !  penetration of tke below the mixed layer (ML) due to NIWs 
    986                            !        = 0 no penetration 
    987                            !        = 1 add a tke source below the ML 
    988                            !        = 2 add a tke source just at the base of the ML 
    989                            !        = 3 as = 1 applied on HF part of the stress           (ln_cpl=T) 
    990    rn_efr      =   0.05    !  fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 
    991    nn_htau     =   1       !  type of exponential decrease of tke penetration below the ML 
    992                            !        = 0  constant 10 m length scale 
    993                            !        = 1  0.5m at the equator to 30m poleward of 40 degrees 
     959                              !        = 0 none ; = 1 add a tke source below the ML 
     960                              !        = 2 add a tke source just at the base of the ML 
     961                              !        = 3 as = 1 applied on HF part of the stress           (ln_cpl=T) 
     962      rn_efr      =   0.05    !  fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 
     963      nn_htau     =   1       !  type of exponential decrease of tke penetration below the ML 
     964                              !        = 0  constant 10 m length scale 
     965                              !        = 1  0.5m at the equator to 30m poleward of 40 degrees 
    994966/ 
    995967!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/ISOMIP/EXP00/namelist_cfg

    r8017 r8143  
    55&namusr_def    !   ISOMIP user defined namelist   
    66!----------------------------------------------------------------------- 
     7   ln_zco      = .false.      ! z-coordinate 
    78   ln_zps      = .true.       ! z-partial-step coordinate 
     9   ln_sco      = .false.      ! s-coordinate 
    810   rn_lam0     =   0.0        !  longitude of first raw and column T-point (jphgr_msh = 1) 
    911   rn_phi0     = -80.0        ! latitude  of first raw and column T-point (jphgr_msh = 1) 
     
    179181/ 
    180182!----------------------------------------------------------------------- 
    181 &nambfr        !   bottom friction 
    182 !----------------------------------------------------------------------- 
    183    nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    184                            !                              = 2 : nonlinear friction 
    185    rn_bfri1    =    4.e-4  !  bottom drag coefficient (linear case) 
    186    rn_bfri2    =    1.e-3  !  bottom drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
    187    rn_bfri2_max =   1.e-1  !  max. bottom drag coefficient (non linear case and ln_loglayer=T) 
    188    rn_bfeb2    =    2.5e-3 !  bottom turbulent kinetic energy background  (m2/s2) 
    189    rn_bfrz0    =    3.e-3  !  bottom roughness [m] if ln_loglayer=T 
    190    ln_bfr2d    = .false.   !  horizontal variation of the bottom friction coef (read a 2D mask file ) 
    191    rn_bfrien   =    50.    !  local multiplying factor of bfr (ln_bfr2d=T) 
    192    rn_tfri1    =    4.e-4  !  top drag coefficient (linear case) 
    193    rn_tfri2    =    2.5e-3 !  top drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
    194    rn_tfri2_max =   1.e-1  !  max. top drag coefficient (non linear case and ln_loglayer=T) 
    195    rn_tfeb2    =    0.0    !  top turbulent kinetic energy background  (m2/s2) 
    196    rn_tfrz0    =    3.e-3  !  top roughness [m] if ln_loglayer=T 
    197    ln_tfr2d    = .false.   !  horizontal variation of the top friction coef (read a 2D mask file ) 
    198    rn_tfrien   =    50.    !  local multiplying factor of tfr (ln_tfr2d=T) 
    199  
    200    ln_bfrimp   = .true.    !  implicit bottom friction (requires ln_zdfexp = .false. if true) 
    201    ln_loglayer = .false.   !  logarithmic formulation (non linear case) 
     183&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     184!----------------------------------------------------------------------- 
     185   ln_non_lin = .false.    !  non-linear  drag: Cd = Cd0 |U| 
     186/ 
     187!----------------------------------------------------------------------- 
     188&namdrg_top        !   TOP friction                                     (ln_isfcav=T) 
     189!----------------------------------------------------------------------- 
     190   rn_Cd0     =  2.5e-3    !  drag coefficient [-] 
     191   rn_Uc0     =  1.6       !  ref. velocity [m/s] (linear drag=Cd0*Uc0)  
     192   rn_Cdmax   =  0.1       !  drag value maximum [-] (logarithmic drag) 
     193   rn_ke0     =  0.0e-0    !  background kinetic energy  [m2/s2] (non-linear cases) 
     194   rn_z0      =  3.0e-3    !  roughness [m] (ln_loglayer=T) 
     195   ln_boost   = .false.    !  =T regional boost of Cd0 ; =F constant 
     196      rn_boost=  50.          !  local boost factor  [-] 
    202197/ 
    203198!----------------------------------------------------------------------- 
     
    215210&nameos        !   ocean physical parameters 
    216211!----------------------------------------------------------------------- 
    217    ln_teos10   = .false.         !  = Use TEOS-10 equation of state 
    218212   ln_eos80    = .true.          !  = Use EOS80 equation of state 
    219    !                             !  rd(T,S,Z)*rau0 = -a0*(1+.5*lambda*dT+mu*Z+nu*dS)*dT+b0*dS 
    220213/ 
    221214!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT2_flux_cen2_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT2_flux_ubs_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT2_vect_eenH_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT2_vect_een_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT2_vect_ene_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT2_vect_ens_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT4_flux_cen2_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT4_flux_ubs_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT4_vect_eenH_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT4_vect_een_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT4_vect_ene_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_FCT4_vect_ens_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/LOCK_EXCHANGE/EXP00/namelist_cfg

    r7954 r8143  
    6868/ 
    6969!----------------------------------------------------------------------- 
    70 &nambfr        !   bottom friction 
    71 !----------------------------------------------------------------------- 
    72    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    73                            !                              = 2 : nonlinear friction 
     70&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     71!----------------------------------------------------------------------- 
     72   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    7473/ 
    7574!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/OVERFLOW/EXP00/namelist_cfg

    r7954 r8143  
    6262/ 
    6363!----------------------------------------------------------------------- 
    64 &nambfr        !   bottom friction 
    65 !----------------------------------------------------------------------- 
    66    nn_bfr      =    0      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    67                            !                              = 2 : nonlinear friction 
     64&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     65!----------------------------------------------------------------------- 
     66   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    6867/ 
    6968!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/SAS_BIPER/EXP00/namelist_cfg

    r7990 r8143  
    8282/ 
    8383!----------------------------------------------------------------------- 
    84 &nambfr        !   bottom friction 
    85 !----------------------------------------------------------------------- 
     84&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     85!----------------------------------------------------------------------- 
     86   ln_NONE    = .false.    !  free-slip       : Cd = 0                   
    8687/ 
    8788!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/CONFIG/TEST_CASES/WAD/EXP00/namelist_cfg

    r7990 r8143  
    173173/ 
    174174!----------------------------------------------------------------------- 
    175 &nambfr        !   bottom friction 
    176 !----------------------------------------------------------------------- 
    177    nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
    178    !rn_bfri2    =    1.e-5  !  bottom drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
    179    !rn_bfri2_max =   1.e-4  !  max. bottom drag coefficient (non linear case and ln_loglayer=T) 
    180    rn_bfri2    =    1.e-5  !  bottom drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
    181    rn_bfri2_max =   1.e-4  !  max. bottom drag coefficient (non linear case and ln_loglayer=T) 
    182    !rn_bfeb2    =    2.5e-3 !  bottom turbulent kinetic energy background  (m2/s2) 
    183    !rn_bfrz0    =    3.e-3  !  bottom roughness [m] if ln_loglayer=T 
    184    ln_loglayer = .true.    !  logarithmic formulation (non linear case) 
     175&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
     176!----------------------------------------------------------------------- 
     177   ln_loglayer= .false.    !  logarithmic drag: Cd = vkarmn/log(z/z0) |U| 
     178/ 
     179!----------------------------------------------------------------------- 
     180&namdrg_bot        !   BOTTOM friction                                   
     181!----------------------------------------------------------------------- 
     182   rn_Cd0     =  1.e-4    !  drag coefficient [-] 
     183   rn_Uc0     =  0.1      !  ref. velocity [m/s] (linear drag=Cd0*Uc0)  
     184   rn_Cdmax   =  1.e-4    !  drag value maximum [-] (logarithmic drag) 
     185   rn_ke0     =  2.5e-3   !  background kinetic energy  [m2/s2] (non-linear cases) 
     186   rn_z0      =  3.e-3    !  roughness [m] (ln_loglayer=T) 
     187   ln_boost   = .false.   !  =T regional boost of Cd0 ; =F constant 
     188      rn_boost=  50.         !  local boost factor  [-] 
    185189/ 
    186190!----------------------------------------------------------------------- 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/CRS/crs.F90

    r8093 r8143  
    151151 
    152152      ! Vertical diffusion 
    153       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)  ::  avs_crs           !: salinity vertical diffusivity coeff. [m2/s] at w-point 
     153      REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)  ::  avt_crs           !: temperature vertical diffusivity coeff. [m2/s] at w-point 
     154      REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)  ::  avs_crs           !: salinity    vertical diffusivity coeff. [m2/s] at w-point 
    154155 
    155156      ! Mixing and Mixed Layer Depth 
     
    235236         &     fr_i_crs(jpi_crs,jpj_crs), sfx_crs(jpi_crs ,jpj_crs),  STAT=ierr(12)  ) 
    236237 
    237      ALLOCATE( tsn_crs(jpi_crs,jpj_crs,jpk,jpts), avs_crs(jpi_crs,jpj_crs,jpk), STAT=ierr(13) ) 
     238     ALLOCATE( tsn_crs(jpi_crs,jpj_crs,jpk,jpts), avt_crs(jpi_crs,jpj_crs,jpk),   & 
     239         &                                        avs_crs(jpi_crs,jpj_crs,jpk), STAT=ierr(13) ) 
    238240 
    239241      ALLOCATE( nmln_crs(jpi_crs,jpj_crs) , hmld_crs(jpi_crs,jpj_crs) , & 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/CRS/crsfld.F90

    r8093 r8143  
    210210      !  free memory 
    211211 
    212       !  avt, avs 
    213 !!gm BUG   TOP always uses avs !!! 
     212      !  avs 
    214213      SELECT CASE ( nn_crs_kz ) 
    215214         CASE ( 0 ) 
    216             CALL crs_dom_ope( avt, 'VOL', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     215            CALL crs_dom_ope( avt, 'VOL', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     216            CALL crs_dom_ope( avs, 'VOL', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
    217217         CASE ( 1 ) 
    218             CALL crs_dom_ope( avt, 'MAX', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     218            CALL crs_dom_ope( avt, 'MAX', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     219            CALL crs_dom_ope( avs, 'MAX', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
    219220         CASE ( 2 ) 
    220             CALL crs_dom_ope( avt, 'MIN', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     221            CALL crs_dom_ope( avt, 'MIN', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     222            CALL crs_dom_ope( avs, 'MIN', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
    221223      END SELECT 
    222224      ! 
    223       CALL iom_put( "avt", avs_crs )   !  Kz 
     225      CALL iom_put( "avt", avt_crs )   !  Kz on T 
     226      CALL iom_put( "avs", avs_crs )   !  Kz on S 
    224227       
    225228      !  sbc fields   
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r8093 r8143  
    6161   USE diurnal_bulk    ! diurnal warm layer 
    6262   USE cool_skin       ! Cool skin 
    63    USE wrk_nemo        ! working array 
    6463 
    6564   IMPLICIT NONE 
     
    183182         DO jj = 2, jpjm1 
    184183            DO ji = fs_2, fs_jpim1   ! vector opt. 
    185 !!gm old 
    186 !!gm BUG  missing x 0.5 
    187                zztmpx = (  bfrua(ji  ,jj) * un(ji  ,jj,mbku(ji  ,jj))  & 
    188                       &  + bfrua(ji-1,jj) * un(ji-1,jj,mbku(ji-1,jj))  )       
    189                zztmpy = (  bfrva(ji,  jj) * vn(ji,jj  ,mbkv(ji,jj  ))  & 
    190                       &  + bfrva(ji,jj-1) * vn(ji,jj-1,mbkv(ji,jj-1))  )  
    191                z2d(ji,jj) = rau0 * SQRT( zztmpx * zztmpx + zztmpy * zztmpy ) * tmask(ji,jj,1)  
    192 !!gm 
    193184               zztmp2 = (  ( rCdU_bot(ji+1,jj)+rCdU_bot(ji  ,jj) ) * un(ji  ,jj,mbku(ji  ,jj))  )**2   & 
    194185                  &   + (  ( rCdU_bot(ji  ,jj)+rCdU_bot(ji-1,jj) ) * un(ji-1,jj,mbku(ji-1,jj))  )**2   & 
     
    196187                  &   + (  ( rCdU_bot(ji,jj  )+rCdU_bot(ji,jj-1) ) * vn(ji,jj-1,mbkv(ji,jj-1))  )**2 
    197188               z2d(ji,jj) = zztmp * SQRT( zztmp2 ) * tmask(ji,jj,1)  
    198 !!gm new end 
    199189               ! 
    200             ENDDO 
    201          ENDDO 
     190            END DO 
     191         END DO 
    202192         CALL lbc_lnk( z2d, 'T', 1. ) 
    203193         CALL iom_put( "taubot", z2d )            
     
    449439      REAL(wp) ::   zsto, zout, zmax, zjulian                ! local scalars 
    450440      ! 
    451       REAL(wp), POINTER, DIMENSION(:,:)   :: zw2d       ! 2D workspace 
    452       REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d       ! 3D workspace 
     441      REAL(wp), DIMENSION(jpi,jpj)   :: zw2d       ! 2D workspace 
     442      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw3d       ! 3D workspace 
    453443      !!---------------------------------------------------------------------- 
    454444      !  
    455445      IF( nn_timing == 1 )   CALL timing_start('dia_wri') 
    456446      ! 
    457                              CALL wrk_alloc( jpi,jpj      , zw2d ) 
    458       IF( .NOT.ln_linssh )   CALL wrk_alloc( jpi,jpj,jpk  , zw3d ) 
    459       ! 
    460       ! Output the initial state and forcings 
    461       IF( ninist == 1 ) THEN                        
     447      IF( ninist == 1 ) THEN     !==  Output the initial state and forcings  ==! 
    462448         CALL dia_wri_state( 'output.init', kt ) 
    463449         ninist = 0 
     
    467453      ! ----------------- 
    468454 
    469       ! local variable for debugging 
    470       ll_print = .FALSE. 
     455      ll_print = .FALSE.                  ! local variable for debugging 
    471456      ll_print = ll_print .AND. lwp 
    472457 
     
    891876      ENDIF 
    892877      ! 
    893                              CALL wrk_dealloc( jpi , jpj        , zw2d ) 
    894       IF( .NOT.ln_linssh )   CALL wrk_dealloc( jpi , jpj , jpk  , zw3d ) 
    895       ! 
    896878      IF( nn_timing == 1 )   CALL timing_stop('dia_wri') 
    897879      ! 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynbfr.F90

    r8093 r8143  
    55   !!============================================================================== 
    66   !! History :  3.2  ! 2008-11  (A. C. Coward)  Original code 
    7    !!            3.4  ! 2011-09  (H. Liu) Make it consistent with semi-implicit 
    8    !!                            Bottom friction (ln_bfrimp = .true.)  
     7   !!            3.4  ! 2011-09  (H. Liu) Make it consistent with semi-implicit Bottom friction (ln_drgimp =T)  
     8   !!            4.0  ! 2017-05  (G. Madec)  drag coef. defined at t-point (zdfdrg.F90) 
    99   !!---------------------------------------------------------------------- 
    1010 
     
    1414   USE oce            ! ocean dynamics and tracers variables 
    1515   USE dom_oce        ! ocean space and time domain variables  
    16    USE zdf_oce        ! ocean vertical physics variables 
    17 !!gm new 
     16   USE zdf_oce        ! vertical physics: variables 
    1817   USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    19 !!gm old 
    20    USE zdfbfr         ! ocean bottom friction variables 
    21 !!gm 
    2218   USE trd_oce        ! trends: ocean variables 
    2319   USE trddyn         ! trend manager: dynamics 
     
    2622   USE prtctl         ! Print control 
    2723   USE timing         ! Timing 
    28    USE wrk_nemo       ! Memory Allocation 
    2924 
    3025   IMPLICIT NONE 
     
    3631#  include "vectopt_loop_substitute.h90" 
    3732   !!---------------------------------------------------------------------- 
    38    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     33   !! NEMO/OPA 4.0 , NEMO Consortium (2017) 
    3934   !! $Id$ 
    4035   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    4742      !! 
    4843      !! ** Purpose :   compute the bottom friction ocean dynamics physics. 
     44      !! 
     45      !!              only for explicit bottom friction form 
     46      !!              implicit bfr is implemented in dynzdf_imp 
    4947      !! 
    5048      !! ** Action  :   (ua,va)   momentum trend increased by bottom friction trend 
     
    6159      IF( nn_timing == 1 )  CALL timing_start('dyn_bfr') 
    6260      ! 
    63 !!gm issue: better to put the logical in step to control the call of zdf_bfr 
    64 !!          ==> change the logical from ln_bfrimp to ln_bfr_exp !! 
    65       IF( .NOT.ln_bfrimp) THEN     ! only for explicit bottom friction form 
    66                                     ! implicit bfr is implemented in dynzdf_imp 
    67  
    6861!!gm bug : time step is only rdt (not 2 rdt if euler start !) 
    6962         zm1_2dt = - 1._wp / ( 2._wp * rdt ) 
    7063 
    71          IF( l_trddyn ) THEN      ! trends: store the input trends 
    72             ALLOCATE( ztrdu(jpi,jpj,jpk) , ztrdv(jpi,jpj,jpk) ) 
    73             ztrdu(:,:,:) = ua(:,:,:) 
     64      IF( l_trddyn ) THEN      ! trends: store the input trends 
     65         ALLOCATE( ztrdu(jpi,jpj,jpk) , ztrdv(jpi,jpj,jpk) ) 
     66         ztrdu(:,:,:) = ua(:,:,:) 
    7467            ztrdv(:,:,:) = va(:,:,:) 
    75          ENDIF 
     68      ENDIF 
    7669 
    7770 
     71      DO jj = 2, jpjm1 
     72         DO ji = 2, jpim1 
     73            ikbu = mbku(ji,jj)          ! deepest wet ocean u- & v-levels 
     74            ikbv = mbkv(ji,jj) 
     75            ! 
     76            ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt) 
     77            zCdu = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / e3u_n(ji,jj,ikbu) 
     78            zCdv = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / e3v_n(ji,jj,ikbv) 
     79            ! 
     80            ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * ub(ji,jj,ikbu) 
     81            va(ji,jj,ikbv) = va(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * vb(ji,jj,ikbv) 
     82         END DO 
     83      END DO 
     84      ! 
     85      IF( ln_isfcav ) THEN        ! ocean cavities 
    7886         DO jj = 2, jpjm1 
    7987            DO ji = 2, jpim1 
    80                ikbu = mbku(ji,jj)          ! deepest wet ocean u- & v-levels 
    81                ikbv = mbkv(ji,jj) 
     88               ikbu = miku(ji,jj)          ! first wet ocean u- & v-levels 
     89               ikbv = mikv(ji,jj) 
    8290               ! 
    8391               ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt) 
    84 !!gm old 
    85                ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + MAX(  bfrua(ji,jj) / e3u_n(ji,jj,ikbu) , zm1_2dt  ) * ub(ji,jj,ikbu) 
    86                va(ji,jj,ikbv) = va(ji,jj,ikbv) + MAX(  bfrva(ji,jj) / e3v_n(ji,jj,ikbv) , zm1_2dt  ) * vb(ji,jj,ikbv) 
    87 !!gm new 
    88 !               zCdu = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / e3u_n(ji,jj,ikbu) 
    89 !               zCdv = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / e3v_n(ji,jj,ikbv) 
    90 !               ! 
    91 !               ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * ub(ji,jj,ikbu) 
    92 !               va(ji,jj,ikbv) = va(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * vb(ji,jj,ikbv) 
    93 !!gm 
    94             END DO 
     92               zCdu = 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / e3u_n(ji,jj,ikbu)    ! NB: Cdtop masked 
     93               zCdv = 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / e3v_n(ji,jj,ikbv) 
     94               ! 
     95               ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * ub(ji,jj,ikbu) 
     96               va(ji,jj,ikbv) = va(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * vb(ji,jj,ikbv) 
     97           END DO 
    9598         END DO 
    96          ! 
    97          IF( ln_isfcav ) THEN        ! ocean cavities 
    98             DO jj = 2, jpjm1 
    99                DO ji = 2, jpim1 
    100                   ikbu = miku(ji,jj)          ! first wet ocean u- & v-levels 
    101                   ikbv = mikv(ji,jj) 
    102                   ! 
    103                   ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt) 
    104 !!gm old 
    105                  ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + MAX(  tfrua(ji,jj) / e3u_n(ji,jj,ikbu) , zm1_2dt  ) * ub(ji,jj,ikbu) & 
    106                     &             * (1.-umask(ji,jj,1)) 
    107                  va(ji,jj,ikbv) = va(ji,jj,ikbv) + MAX(  tfrva(ji,jj) / e3v_n(ji,jj,ikbv) , zm1_2dt  ) * vb(ji,jj,ikbv) & 
    108                     &             * (1.-vmask(ji,jj,1)) 
    109 !!gm new 
    110 !                  zCdu = 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / e3u_n(ji,jj,ikbu)    ! NB: Cdtop masked 
    111 !                  zCdv = 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / e3v_n(ji,jj,ikbv) 
    112 !                  ! 
    113 !                  ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * ub(ji,jj,ikbu) 
    114 !                  va(ji,jj,ikbv) = va(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * vb(ji,jj,ikbv) 
    115 !!gm 
    116               END DO 
    117            END DO 
    118          END IF 
    119         ! 
    120         IF( l_trddyn ) THEN      ! trends: send trends to trddyn for further diagnostics 
    121            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    122            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    123            CALL trd_dyn( ztrdu(:,:,:), ztrdv(:,:,:), jpdyn_bfr, kt ) 
    124            DEALLOCATE( ztrdu, ztrdv ) 
    125         ENDIF 
    126         !                                          ! print mean trends (used for debugging) 
    127         IF(ln_ctl)   CALL prt_ctl( tab3d_1=ua, clinfo1=' bfr  - Ua: ', mask1=umask,               & 
    128            &                       tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    129         ! 
    130       ENDIF     ! end explicit bottom friction 
     99      ENDIF 
     100      ! 
     101      IF( l_trddyn ) THEN      ! trends: send trends to trddyn for further diagnostics 
     102         ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
     103         ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
     104         CALL trd_dyn( ztrdu(:,:,:), ztrdv(:,:,:), jpdyn_bfr, kt ) 
     105         DEALLOCATE( ztrdu, ztrdv ) 
     106      ENDIF 
     107      !                                          ! print mean trends (used for debugging) 
     108      IF(ln_ctl)   CALL prt_ctl( tab3d_1=ua, clinfo1=' bfr  - Ua: ', mask1=umask,               & 
     109         &                       tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    131110      ! 
    132111      IF( nn_timing == 1 )  CALL timing_stop('dyn_bfr') 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r8093 r8143  
    1616   !!             3.7  ! 2015-11  (J. Chanut) free surface simplification 
    1717   !!              -   ! 2016-12  (G. Madec, E. Clementi) update for Stoke-Drift divergence 
     18   !!             4.0  ! 2017-05  (G. Madec)  drag coef. defined at t-point (zdfdrg.F90) 
    1819   !!--------------------------------------------------------------------- 
    1920 
     
    2728   USE dom_oce         ! ocean space and time domain 
    2829   USE sbc_oce         ! surface boundary condition: ocean 
    29    USE zdf_oce         ! Bottom friction coefts 
    30 !!gm new 
    31    USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    32 !!gm 
     30   USE zdf_oce         ! vertical physics: variables 
     31   USE zdfdrg          ! vertical physics: top/bottom drag coef. 
    3332   USE sbcisf          ! ice shelf variable (fwfisf) 
    3433   USE sbcapr          ! surface boundary condition: atmospheric pressure 
     
    4342   USE updtide         ! tide potential 
    4443   USE sbcwave         ! surface wave 
     44   USE diatmb          ! Top,middle,bottom output 
     45#if defined key_agrif 
     46   USE agrif_opa_interp ! agrif 
     47#endif 
     48#if defined key_asminc    
     49   USE asminc          ! Assimilation increment 
     50#endif 
    4551   ! 
    4652   USE in_out_manager  ! I/O manager 
     
    5056   USE iom             ! IOM library 
    5157   USE restart         ! only for lrst_oce 
    52    USE wrk_nemo        ! Memory Allocation 
    5358   USE timing          ! Timing     
    54    USE diatmb          ! Top,middle,bottom output 
    55 #if defined key_agrif 
    56    USE agrif_opa_interp ! agrif 
    57 #endif 
    58 #if defined key_asminc    
    59    USE asminc          ! Assimilation increment 
    60 #endif 
    61  
    6259 
    6360   IMPLICIT NONE 
     
    6966   PUBLIC ts_rst            !    "      "     "    " 
    7067 
    71    INTEGER, SAVE :: icycle  ! Number of barotropic sub-steps for each internal step nn_baro <= 2.5 nn_baro 
    72    REAL(wp),SAVE :: rdtbt   ! Barotropic time step 
    73  
    74    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::   wgtbtp1, wgtbtp2   !: 1st & 2nd weights used in time filtering of barotropic fields 
    75  
    76    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  zwz          !: ff_f/h at F points 
    77    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  ftnw, ftne   !: triad of coriolis parameter 
    78    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  ftsw, ftse   !: (only used with een vorticity scheme) 
    79  
    8068   !! Time filtered arrays at baroclinic time step: 
    81    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   un_adv , vn_adv     !: Advection vel. at "now" barocl. step 
     69   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   un_adv , vn_adv   !: Advection vel. at "now" barocl. step 
     70 
     71   INTEGER , SAVE ::   icycle   ! Number of barotropic sub-steps for each internal step nn_baro <= 2.5 nn_baro 
     72   REAL(wp), SAVE ::   rdtbt    ! Barotropic time step 
     73   ! 
     74   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)   ::   wgtbtp1, wgtbtp2   ! 1st & 2nd weights used in time filtering of barotropic fields 
     75   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  zwz                 ! ff_f/h at F points 
     76   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  ftnw, ftne          ! triad of coriolis parameter 
     77   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  ftsw, ftse          ! (only used with een vorticity scheme) 
     78 
     79   REAL(wp) ::   r1_12 = 1._wp / 12._wp   ! local ratios 
     80   REAL(wp) ::   r1_8  = 0.125_wp         ! 
     81   REAL(wp) ::   r1_4  = 0.25_wp          ! 
     82   REAL(wp) ::   r1_2  = 0.5_wp           ! 
    8283 
    8384   !! * Substitutions 
    8485#  include "vectopt_loop_substitute.h90" 
    8586   !!---------------------------------------------------------------------- 
    86    !! NEMO/OPA 3.5 , NEMO Consortium (2013) 
     87   !! NEMO/OPA 4.0 , NEMO Consortium (2017) 
    8788   !! $Id$ 
    8889   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    140141      INTEGER, INTENT(in)  ::   kt   ! ocean time-step index 
    141142      ! 
    142       LOGICAL  ::   ll_fw_start        ! if true, forward integration  
    143       LOGICAL  ::   ll_init             ! if true, special startup of 2d equations 
    144       LOGICAL  ::   ll_tmp1, ll_tmp2            ! local logical variables used in W/D 
    145143      INTEGER  ::   ji, jj, jk, jn        ! dummy loop indices 
    146       INTEGER  ::   ikbu, ikbv, noffset      ! local integers 
    147       INTEGER  ::   iktu, iktv               ! local integers 
    148       REAL(wp) ::   zmdi 
    149       REAL(wp) ::   zraur, z1_2dt_b, z2dt_bf    ! local scalars 
    150       REAL(wp) ::   zx1, zy1, zx2, zy2          !   -      - 
    151       REAL(wp) ::   z1_12, z1_8, z1_4, z1_2     !   -      - 
    152       REAL(wp) ::   zu_spg, zv_spg              !   -      - 
    153       REAL(wp) ::   zhura, zhvra                !   -      - 
    154       REAL(wp) ::   za0, za1, za2, za3          !   -      - 
    155       REAL(wp) ::   zztmp                       !   -      - 
    156       REAL(wp), DIMENSION(jpi,jpj) :: zsshp2_e 
    157       REAL(wp), DIMENSION(jpi,jpj) :: zu_trd, zv_trd, zu_frc, zv_frc, zssh_frc 
    158       REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zhdiv 
    159       REAL(wp), DIMENSION(jpi,jpj) :: zhup2_e, zhvp2_e, zhust_e, zhvst_e 
    160       REAL(wp), DIMENSION(jpi,jpj) :: zsshu_a, zsshv_a 
    161       REAL(wp), DIMENSION(jpi,jpj) :: zhf 
     144      LOGICAL  ::   ll_fw_start           ! =T : forward integration  
     145      LOGICAL  ::   ll_init               ! =T : special startup of 2d equations 
     146      LOGICAL  ::   ll_tmp1, ll_tmp2      ! local logical variables used in W/D 
     147      INTEGER  ::   ikbu, iktu, noffset   ! local integers 
     148      INTEGER  ::   ikbv, iktv            !   -      - 
     149      REAL(wp) ::   z1_2dt_b, z2dt_bf        ! local scalars 
     150      REAL(wp) ::   zx1, zx2, zu_spg, zhura  !   -      - 
     151      REAL(wp) ::   zy1, zy2, zv_spg, zhvra  !   -      - 
     152      REAL(wp) ::   za0, za1, za2, za3       !   -      - 
     153      REAL(wp) ::   zmdi, zztmp              !   -      - 
     154      REAL(wp), DIMENSION(jpi,jpj) :: zsshp2_e, zhf 
     155      REAL(wp), DIMENSION(jpi,jpj) :: zwx, zu_trd, zu_frc, zssh_frc 
     156      REAL(wp), DIMENSION(jpi,jpj) :: zwy, zv_trd, zv_frc, zhdiv 
     157      REAL(wp), DIMENSION(jpi,jpj) :: zsshu_a, zhup2_e, zhust_e 
     158      REAL(wp), DIMENSION(jpi,jpj) :: zsshv_a, zhvp2_e, zhvst_e 
    162159      REAL(wp), DIMENSION(jpi,jpj) :: zCdU_u, zCdU_v   ! top/bottom stress at u- & v-points 
    163160      ! 
     
    170167      ! 
    171168      zmdi=1.e+20                               !  missing data indicator for masking 
    172       !                                         !* Local constant initialization 
    173       z1_12 = 1._wp / 12._wp  
    174       z1_8  = 0.125_wp                                    
    175       z1_4  = 0.25_wp 
    176       z1_2  = 0.5_wp      
    177       zraur = 1._wp / rau0 
     169      ! 
    178170      !                                            ! reciprocal of baroclinic time step  
    179171      IF( kt == nit000 .AND. neuler == 0 ) THEN   ;   z2dt_bf =          rdt 
     
    210202      ENDIF 
    211203      ! 
    212 !!gm old/new 
    213204      IF( ln_isfcav ) THEN    ! top+bottom friction (ocean cavities) 
    214          zCdU_u(:,:) = bfrua(:,:) + tfrua(:,:) 
    215          zCdU_v(:,:) = bfrva(:,:) + tfrva(:,:) 
     205         DO jj = 2, jpjm1 
     206            DO ji = fs_2, fs_jpim1   ! vector opt. 
     207               zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 
     208               zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 
     209            END DO 
     210         END DO 
    216211      ELSE                    ! bottom friction only 
    217          zCdU_u(:,:) = bfrua(:,:) 
    218          zCdU_v(:,:) = bfrva(:,:) 
    219       ENDIF 
    220 !!gm new 
    221 !      IF( ln_isfcav ) THEN    ! top+bottom friction (ocean cavities) 
    222 !         DO jj = 2, jpjm1 
    223 !            DO ji = fs_2, fs_jpim1   ! vector opt. 
    224 !               zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 
    225 !               zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 
    226 !            END DO 
    227 !         END DO 
    228 !      ELSE                    ! bottom friction only 
    229 !         DO jj = 2, jpjm1 
    230 !            DO ji = fs_2, fs_jpim1   ! vector opt. 
    231 !               zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 
    232 !               zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) 
    233 !            END DO 
    234 !         END DO 
    235 !      ENDIF 
    236 !!gm       
     212         DO jj = 2, jpjm1 
     213            DO ji = fs_2, fs_jpim1   ! vector opt. 
     214               zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 
     215               zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) 
     216            END DO 
     217         END DO 
     218      ENDIF 
    237219      ! 
    238220      ! Set arrays to remove/compute coriolis trend. 
     
    287269!!gm  
    288270!!             
    289               IF ( .not. ln_sco ) THEN 
     271              IF( .NOT.ln_sco ) THEN 
    290272   
    291273   !!gm  agree the JC comment  : this should be done in a much clear way 
     
    338320      IF (.NOT.ln_bt_fw .AND.( neuler==0 .AND. kt==nit000+1 ) ) THEN 
    339321         ll_fw_start=.FALSE. 
    340          CALL ts_wgt(ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2) 
     322         CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 
    341323      ENDIF 
    342324                           
     
    387369               zx2 = ( zwx(ji  ,jj) + zwx(ji  ,jj+1) ) * r1_e2v(ji,jj) 
    388370               ! energy conserving formulation for planetary vorticity term 
    389                zu_trd(ji,jj) = z1_4 * ( zwz(ji  ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 
    390                zv_trd(ji,jj) =-z1_4 * ( zwz(ji-1,jj  ) * zx1 + zwz(ji,jj) * zx2 ) 
     371               zu_trd(ji,jj) =   r1_4 * ( zwz(ji  ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 
     372               zv_trd(ji,jj) = - r1_4 * ( zwz(ji-1,jj  ) * zx1 + zwz(ji,jj) * zx2 ) 
    391373            END DO 
    392374         END DO 
     
    395377         DO jj = 2, jpjm1 
    396378            DO ji = fs_2, fs_jpim1   ! vector opt. 
    397                zy1 =   z1_8 * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1) & 
     379               zy1 =   r1_8 * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1) & 
    398380                 &            + zwy(ji  ,jj  ) + zwy(ji+1,jj  ) ) * r1_e1u(ji,jj) 
    399                zx1 = - z1_8 * ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1) & 
     381               zx1 = - r1_8 * ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1) & 
    400382                 &            + zwx(ji  ,jj  ) + zwx(ji  ,jj+1) ) * r1_e2v(ji,jj) 
    401383               zu_trd(ji,jj)  = zy1 * ( zwz(ji  ,jj-1) + zwz(ji,jj) ) 
     
    407389         DO jj = 2, jpjm1 
    408390            DO ji = fs_2, fs_jpim1   ! vector opt. 
    409                zu_trd(ji,jj) = + z1_12 * r1_e1u(ji,jj) * (  ftne(ji,jj  ) * zwy(ji  ,jj  ) & 
     391               zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * (  ftne(ji,jj  ) * zwy(ji  ,jj  ) & 
    410392                &                                         + ftnw(ji+1,jj) * zwy(ji+1,jj  ) & 
    411393                &                                         + ftse(ji,jj  ) * zwy(ji  ,jj-1) & 
    412394                &                                         + ftsw(ji+1,jj) * zwy(ji+1,jj-1) ) 
    413                zv_trd(ji,jj) = - z1_12 * r1_e2v(ji,jj) * (  ftsw(ji,jj+1) * zwx(ji-1,jj+1) & 
     395               zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * (  ftsw(ji,jj+1) * zwx(ji-1,jj+1) & 
    414396                &                                         + ftse(ji,jj+1) * zwx(ji  ,jj+1) & 
    415397                &                                         + ftnw(ji,jj  ) * zwx(ji-1,jj  ) & 
     
    423405      !                                   ! ---------------------------------------------------- 
    424406      IF( .NOT.ln_linssh ) THEN                 ! Variable volume : remove surface pressure gradient 
    425         IF( ln_wd ) THEN                        ! Calculating and applying W/D gravity filters 
    426            DO jj = 2, jpjm1 
    427               DO ji = 2, jpim1  
    428                 ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
    429                      &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) .AND.            & 
    430                      &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 
     407         IF( ln_wd ) THEN                        ! Calculating and applying W/D gravity filters 
     408            DO jj = 2, jpjm1 
     409               DO ji = 2, jpim1  
     410                  ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     411                     &      MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) .AND.            & 
     412                     &      MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 
    431413                     &                                                         > rn_wdmin1 + rn_wdmin2 
    432                 ll_tmp2 = ( ABS( sshn(ji+1,jj)             -   sshn(ji  ,jj))  > 1.E-12 ).AND.( & 
    433                      &    MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
    434                      &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 ) 
    435     
    436                 IF(ll_tmp1) THEN 
    437                   zcpx(ji,jj) = 1.0_wp 
    438                 ELSE IF(ll_tmp2) THEN 
    439                   ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
    440                   zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 
    441                      &             / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    442                 ELSE 
    443                   zcpx(ji,jj) = 0._wp 
    444                 END IF 
    445           
    446                 ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     414                  ll_tmp2 = ( ABS( sshn(ji+1,jj)             -   sshn(ji  ,jj))  > 1.E-12 ).AND.( & 
     415                     &      MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     416                     &      MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 ) 
     417                     ! 
     418                  IF(ll_tmp1) THEN 
     419                     zcpx(ji,jj) = 1.0_wp 
     420                  ELSE IF(ll_tmp2) THEN   ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     421                     zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     422                        &             / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
     423                  ELSE 
     424                     zcpx(ji,jj) = 0._wp 
     425                  ENDIF 
     426                  ! 
     427                  ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
    447428                     &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) .AND.            & 
    448429                     &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 
    449430                     &                                                         > rn_wdmin1 + rn_wdmin2 
    450                 ll_tmp2 = ( ABS( sshn(ji,jj)               -   sshn(ji,jj+1))  > 1.E-12 ).AND.( & 
     431                  ll_tmp2 = ( ABS( sshn(ji,jj)               -   sshn(ji,jj+1))  > 1.E-12 ).AND.( & 
    451432                     &    MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
    452433                     &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 ) 
    453     
    454                 IF(ll_tmp1) THEN 
    455                   zcpy(ji,jj) = 1.0_wp 
    456                 ELSE IF(ll_tmp2) THEN 
    457                   ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
    458                   zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 
    459                               &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    460                 ELSE 
    461                   zcpy(ji,jj) = 0._wp 
    462                 END IF 
    463               END DO 
    464            END DO 
    465   
    466            DO jj = 2, jpjm1 
    467               DO ji = 2, jpim1 
    468                  zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj  ) - sshn(ji  ,jj ) )   & 
    469                         &                        * r1_e1u(ji,jj) * zcpx(ji,jj) 
    470                  zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji  ,jj+1) - sshn(ji  ,jj ) )   & 
    471                         &                        * r1_e2v(ji,jj) * zcpy(ji,jj) 
    472               END DO 
    473            END DO 
    474  
     434                     ! 
     435                  IF(ll_tmp1) THEN 
     436                     zcpy(ji,jj) = 1.0_wp 
     437                  ELSE IF(ll_tmp2) THEN 
     438                     ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     439                     zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     440                        &             / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
     441                  ELSE 
     442                     zcpy(ji,jj) = 0._wp 
     443                  ENDIF 
     444               END DO 
     445            END DO 
     446            ! 
     447            DO jj = 2, jpjm1 
     448               DO ji = 2, jpim1 
     449                  zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj  ) - sshn(ji  ,jj ) )   & 
     450                     &                            * r1_e1u(ji,jj) * zcpx(ji,jj) 
     451                  zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji  ,jj+1) - sshn(ji  ,jj ) )   & 
     452                     &                            * r1_e2v(ji,jj) * zcpy(ji,jj) 
     453               END DO 
     454            END DO 
     455            ! 
    475456         ELSE 
    476  
    477            DO jj = 2, jpjm1 
    478               DO ji = fs_2, fs_jpim1   ! vector opt. 
    479                  zu_trd(ji,jj) = zu_trd(ji,jj) - grav * (  sshn(ji+1,jj  ) - sshn(ji  ,jj  )  ) * r1_e1u(ji,jj) 
    480                  zv_trd(ji,jj) = zv_trd(ji,jj) - grav * (  sshn(ji  ,jj+1) - sshn(ji  ,jj  )  ) * r1_e2v(ji,jj)  
    481               END DO 
    482            END DO 
    483         ENDIF 
    484  
    485       ENDIF 
    486  
     457            ! 
     458            DO jj = 2, jpjm1 
     459               DO ji = fs_2, fs_jpim1   ! vector opt. 
     460                  zu_trd(ji,jj) = zu_trd(ji,jj) - grav * (  sshn(ji+1,jj  ) - sshn(ji  ,jj  )  ) * r1_e1u(ji,jj) 
     461                  zv_trd(ji,jj) = zv_trd(ji,jj) - grav * (  sshn(ji  ,jj+1) - sshn(ji  ,jj  )  ) * r1_e2v(ji,jj)  
     462               END DO 
     463            END DO 
     464         ENDIF 
     465         ! 
     466      ENDIF 
     467      ! 
    487468      DO jj = 2, jpjm1                          ! Remove coriolis term (and possibly spg) from barotropic trend 
    488469         DO ji = fs_2, fs_jpim1 
     
    492473      END DO  
    493474      ! 
    494       !                 ! Add bottom stress contribution from baroclinic velocities:       
    495       IF (ln_bt_fw) THEN 
     475      !                 ! Add BOTTOM stress contribution from baroclinic velocities:       
     476      IF( ln_bt_fw ) THEN 
    496477         DO jj = 2, jpjm1                           
    497478            DO ji = fs_2, fs_jpim1   ! vector opt. 
     
    518499         DO jj = 2, jpjm1                           
    519500            DO ji = fs_2, fs_jpim1   ! vector opt. 
    520 !!gm old 
    521                zu_frc(ji,jj) = zu_frc(ji,jj) + MAX( r1_hu_n(ji,jj) * bfrua(ji,jj) , zztmp ) * zwx(ji,jj) 
    522                zv_frc(ji,jj) = zv_frc(ji,jj) + MAX( r1_hv_n(ji,jj) * bfrva(ji,jj) , zztmp ) * zwy(ji,jj) 
    523 !!gm new 
    524 !               zu_frc(ji,jj) = zu_frc(ji,jj) + MAX( r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) * zwx(ji,jj) 
    525 !               zv_frc(ji,jj) = zv_frc(ji,jj) + MAX( r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) , zztmp ) * zwy(ji,jj) 
    526 !!gm 
     501               zu_frc(ji,jj) = zu_frc(ji,jj) + MAX( r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) * zwx(ji,jj) 
     502               zv_frc(ji,jj) = zv_frc(ji,jj) + MAX( r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) , zztmp ) * zwy(ji,jj) 
    527503            END DO 
    528504         END DO 
     
    530506         DO jj = 2, jpjm1                           
    531507            DO ji = fs_2, fs_jpim1   ! vector opt. 
    532 !!gm old 
    533                zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * bfrua(ji,jj) * zwx(ji,jj) 
    534                zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * bfrva(ji,jj) * zwy(ji,jj) 
    535 !!gm new 
    536 !               zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zwx(ji,jj) 
    537 !               zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zwy(ji,jj) 
    538 !!gm 
     508               zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zwx(ji,jj) 
     509               zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zwy(ji,jj) 
    539510            END DO 
    540511         END DO 
    541512      END IF 
    542513      ! 
    543       !                                         ! Add top stress contribution from baroclinic velocities:       
    544       IF( ln_bt_fw ) THEN 
     514      IF( ln_isfcav ) THEN       ! Add TOP stress contribution from baroclinic velocities:       
     515         IF( ln_bt_fw ) THEN 
     516            DO jj = 2, jpjm1 
     517               DO ji = fs_2, fs_jpim1   ! vector opt. 
     518                  iktu = miku(ji,jj) 
     519                  iktv = mikv(ji,jj) 
     520                  zwx(ji,jj) = un(ji,jj,iktu) - un_b(ji,jj) ! NOW top baroclinic velocities 
     521                  zwy(ji,jj) = vn(ji,jj,iktv) - vn_b(ji,jj) 
     522               END DO 
     523            END DO 
     524         ELSE 
     525            DO jj = 2, jpjm1 
     526               DO ji = fs_2, fs_jpim1   ! vector opt. 
     527                  iktu = miku(ji,jj) 
     528                  iktv = mikv(ji,jj) 
     529                  zwx(ji,jj) = ub(ji,jj,iktu) - ub_b(ji,jj) ! BEFORE top baroclinic velocities 
     530                  zwy(ji,jj) = vb(ji,jj,iktv) - vb_b(ji,jj) 
     531               END DO 
     532            END DO 
     533         ENDIF 
     534         ! 
     535         ! Note that the "unclipped" top friction parameter is used even with explicit drag 
     536         DO jj = 2, jpjm1               
     537            DO ji = fs_2, fs_jpim1   ! vector opt. 
     538               zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * zwx(ji,jj) 
     539               zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * zwy(ji,jj) 
     540            END DO 
     541         END DO 
     542      ENDIF 
     543      !        
     544      IF( ln_bt_fw ) THEN                        ! Add wind forcing 
    545545         DO jj = 2, jpjm1 
    546546            DO ji = fs_2, fs_jpim1   ! vector opt. 
    547                iktu = miku(ji,jj) 
    548                iktv = mikv(ji,jj) 
    549                zwx(ji,jj) = un(ji,jj,iktu) - un_b(ji,jj) ! NOW top baroclinic velocities 
    550                zwy(ji,jj) = vn(ji,jj,iktv) - vn_b(ji,jj) 
     547               zu_frc(ji,jj) =  zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu_n(ji,jj) 
     548               zv_frc(ji,jj) =  zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv_n(ji,jj) 
    551549            END DO 
    552550         END DO 
    553551      ELSE 
     552         zztmp = r1_rau0 * r1_2 
    554553         DO jj = 2, jpjm1 
    555554            DO ji = fs_2, fs_jpim1   ! vector opt. 
    556                iktu = miku(ji,jj) 
    557                iktv = mikv(ji,jj) 
    558                zwx(ji,jj) = ub(ji,jj,iktu) - ub_b(ji,jj) ! BEFORE top baroclinic velocities 
    559                zwy(ji,jj) = vb(ji,jj,iktv) - vb_b(ji,jj) 
    560             END DO 
    561          END DO 
    562       ENDIF 
    563       ! 
    564       ! Note that the "unclipped" top friction parameter is used even with explicit drag 
    565       DO jj = 2, jpjm1               
    566          DO ji = fs_2, fs_jpim1   ! vector opt. 
    567             zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * tfrua(ji,jj) * zwx(ji,jj) 
    568             zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * tfrva(ji,jj) * zwy(ji,jj) 
    569          END DO 
    570       END DO 
    571       !        
    572       IF (ln_bt_fw) THEN                        ! Add wind forcing 
    573          DO jj = 2, jpjm1 
    574             DO ji = fs_2, fs_jpim1   ! vector opt. 
    575                zu_frc(ji,jj) =  zu_frc(ji,jj) + zraur * utau(ji,jj) * r1_hu_n(ji,jj) 
    576                zv_frc(ji,jj) =  zv_frc(ji,jj) + zraur * vtau(ji,jj) * r1_hv_n(ji,jj) 
    577             END DO 
    578          END DO 
    579       ELSE 
    580          DO jj = 2, jpjm1 
    581             DO ji = fs_2, fs_jpim1   ! vector opt. 
    582                zu_frc(ji,jj) =  zu_frc(ji,jj) + zraur * z1_2 * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 
    583                zv_frc(ji,jj) =  zv_frc(ji,jj) + zraur * z1_2 * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 
     555               zu_frc(ji,jj) =  zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 
     556               zv_frc(ji,jj) =  zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 
    584557            END DO 
    585558         END DO 
    586559      ENDIF   
    587560      ! 
    588       IF ( ln_apr_dyn ) THEN                    ! Add atm pressure forcing 
    589          IF (ln_bt_fw) THEN 
     561      IF( ln_apr_dyn ) THEN                     ! Add atm pressure forcing 
     562         IF( ln_bt_fw ) THEN 
    590563            DO jj = 2, jpjm1               
    591564               DO ji = fs_2, fs_jpim1   ! vector opt. 
     
    597570            END DO 
    598571         ELSE 
     572            zztmp = grav * r1_2 
    599573            DO jj = 2, jpjm1               
    600574               DO ji = fs_2, fs_jpim1   ! vector opt. 
    601                   zu_spg =  grav * z1_2 * (  ssh_ib (ji+1,jj  ) - ssh_ib (ji,jj)    & 
    602                       &                    + ssh_ibb(ji+1,jj  ) - ssh_ibb(ji,jj)  ) * r1_e1u(ji,jj) 
    603                   zv_spg =  grav * z1_2 * (  ssh_ib (ji  ,jj+1) - ssh_ib (ji,jj)    & 
    604                       &                    + ssh_ibb(ji  ,jj+1) - ssh_ibb(ji,jj)  ) * r1_e2v(ji,jj) 
     575                  zu_spg = zztmp * (  ssh_ib (ji+1,jj  ) - ssh_ib (ji,jj)    & 
     576                      &             + ssh_ibb(ji+1,jj  ) - ssh_ibb(ji,jj)  ) * r1_e1u(ji,jj) 
     577                  zv_spg = zztmp * (  ssh_ib (ji  ,jj+1) - ssh_ib (ji,jj)    & 
     578                      &             + ssh_ibb(ji  ,jj+1) - ssh_ibb(ji,jj)  ) * r1_e2v(ji,jj) 
    605579                  zu_frc(ji,jj) = zu_frc(ji,jj) + zu_spg 
    606580                  zv_frc(ji,jj) = zv_frc(ji,jj) + zv_spg 
     
    613587      !                                         ! Surface net water flux and rivers 
    614588      IF (ln_bt_fw) THEN 
    615          zssh_frc(:,:) = zraur * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) 
     589         zssh_frc(:,:) = r1_rau0 * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) 
    616590      ELSE 
    617          zssh_frc(:,:) = zraur * z1_2 * (  emp(:,:) + emp_b(:,:) - rnf(:,:) - rnf_b(:,:)   & 
    618                 &                        + fwfisf(:,:) + fwfisf_b(:,:)                     ) 
     591         zztmp = r1_rau0 * r1_2 
     592         zssh_frc(:,:) = zztmp * (  emp(:,:) + emp_b(:,:) - rnf(:,:) - rnf_b(:,:)   & 
     593                &                 + fwfisf(:,:) + fwfisf_b(:,:)                     ) 
    619594      ENDIF 
    620595      ! 
     
    712687            DO jj = 2, jpjm1                                    ! Sea Surface Height at u- & v-points 
    713688               DO ji = 2, fs_jpim1   ! Vector opt. 
    714                   zwx(ji,jj) = z1_2 * ssumask(ji,jj)  * r1_e1e2u(ji,jj)     & 
     689                  zwx(ji,jj) = r1_2 * ssumask(ji,jj)  * r1_e1e2u(ji,jj)     & 
    715690                     &              * ( e1e2t(ji  ,jj) * zsshp2_e(ji  ,jj)  & 
    716691                     &              +   e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) 
    717                   zwy(ji,jj) = z1_2 * ssvmask(ji,jj)  * r1_e1e2v(ji,jj)     & 
     692                  zwy(ji,jj) = r1_2 * ssvmask(ji,jj)  * r1_e1e2v(ji,jj)     & 
    718693                     &              * ( e1e2t(ji,jj  ) * zsshp2_e(ji,jj  )  & 
    719694                     &              +   e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) 
     
    789764            DO jj = 2, jpjm1 
    790765               DO ji = 2, jpim1      ! NO Vector Opt. 
    791                   zsshu_a(ji,jj) = z1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj)    & 
     766                  zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj)    & 
    792767                     &              * ( e1e2t(ji  ,jj  )  * ssha_e(ji  ,jj  ) & 
    793768                     &              +   e1e2t(ji+1,jj  )  * ssha_e(ji+1,jj  ) ) 
    794                   zsshv_a(ji,jj) = z1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj)    & 
     769                  zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj)    & 
    795770                     &              * ( e1e2t(ji  ,jj  )  * ssha_e(ji  ,jj  ) & 
    796771                     &              +   e1e2t(ji  ,jj+1)  * ssha_e(ji  ,jj+1) ) 
     
    868843            DO jj = 2, jpjm1                             
    869844               DO ji = 2, jpim1 
    870                   zx1 = z1_2 * ssumask(ji  ,jj) *  r1_e1e2u(ji  ,jj)    & 
     845                  zx1 = r1_2 * ssumask(ji  ,jj) *  r1_e1e2u(ji  ,jj)    & 
    871846                     &      * ( e1e2t(ji  ,jj  ) * zsshp2_e(ji  ,jj)    & 
    872847                     &      +   e1e2t(ji+1,jj  ) * zsshp2_e(ji+1,jj  ) ) 
    873                   zy1 = z1_2 * ssvmask(ji  ,jj) *  r1_e1e2v(ji  ,jj  )  & 
     848                  zy1 = r1_2 * ssvmask(ji  ,jj) *  r1_e1e2v(ji  ,jj  )  & 
    874849                     &       * ( e1e2t(ji ,jj  ) * zsshp2_e(ji  ,jj  )  & 
    875850                     &       +   e1e2t(ji ,jj+1) * zsshp2_e(ji  ,jj+1) ) 
     
    895870                  zx1 = ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1) ) * r1_e2v(ji,jj) 
    896871                  zx2 = ( zwx(ji  ,jj  ) + zwx(ji  ,jj+1) ) * r1_e2v(ji,jj) 
    897                   zu_trd(ji,jj) = z1_4 * ( zwz(ji  ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 
    898                   zv_trd(ji,jj) =-z1_4 * ( zwz(ji-1,jj  ) * zx1 + zwz(ji,jj) * zx2 ) 
     872                  zu_trd(ji,jj) = r1_4 * ( zwz(ji  ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 
     873                  zv_trd(ji,jj) =-r1_4 * ( zwz(ji-1,jj  ) * zx1 + zwz(ji,jj) * zx2 ) 
    899874               END DO 
    900875            END DO 
     
    903878            DO jj = 2, jpjm1 
    904879               DO ji = fs_2, fs_jpim1   ! vector opt. 
    905                   zy1 =   z1_8 * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1) & 
     880                  zy1 =   r1_8 * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1) & 
    906881                   &             + zwy(ji  ,jj  ) + zwy(ji+1,jj  ) ) * r1_e1u(ji,jj) 
    907                   zx1 = - z1_8 * ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1) & 
     882                  zx1 = - r1_8 * ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1) & 
    908883                   &             + zwx(ji  ,jj  ) + zwx(ji  ,jj+1) ) * r1_e2v(ji,jj) 
    909884                  zu_trd(ji,jj)  = zy1 * ( zwz(ji  ,jj-1) + zwz(ji,jj) ) 
     
    915890            DO jj = 2, jpjm1 
    916891               DO ji = fs_2, fs_jpim1   ! vector opt. 
    917                   zu_trd(ji,jj) = + z1_12 * r1_e1u(ji,jj) * (  ftne(ji,jj  ) * zwy(ji  ,jj  ) & 
     892                  zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * (  ftne(ji,jj  ) * zwy(ji  ,jj  ) & 
    918893                     &                                       + ftnw(ji+1,jj) * zwy(ji+1,jj  ) & 
    919894                     &                                       + ftse(ji,jj  ) * zwy(ji  ,jj-1) &  
    920895                     &                                       + ftsw(ji+1,jj) * zwy(ji+1,jj-1) ) 
    921                   zv_trd(ji,jj) = - z1_12 * r1_e2v(ji,jj) * (  ftsw(ji,jj+1) * zwx(ji-1,jj+1) &  
     896                  zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * (  ftsw(ji,jj+1) * zwx(ji-1,jj+1) &  
    922897                     &                                       + ftse(ji,jj+1) * zwx(ji  ,jj+1) & 
    923898                     &                                       + ftnw(ji,jj  ) * zwx(ji-1,jj  ) &  
     
    10821057         vn_adv(:,:) = zwy(:,:) * r1_hv_n(:,:) 
    10831058      ELSE 
    1084          un_adv(:,:) = z1_2 * ( ub2_b(:,:) + zwx(:,:) ) * r1_hu_n(:,:) 
    1085          vn_adv(:,:) = z1_2 * ( vb2_b(:,:) + zwy(:,:) ) * r1_hv_n(:,:) 
     1059         un_adv(:,:) = r1_2 * ( ub2_b(:,:) + zwx(:,:) ) * r1_hu_n(:,:) 
     1060         vn_adv(:,:) = r1_2 * ( vb2_b(:,:) + zwy(:,:) ) * r1_hv_n(:,:) 
    10861061      END IF 
    10871062 
     
    11011076         DO jj = 1, jpjm1 
    11021077            DO ji = 1, jpim1      ! NO Vector Opt. 
    1103                zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1)  * r1_e1e2u(ji,jj) & 
     1078               zsshu_a(ji,jj) = r1_2 * umask(ji,jj,1)  * r1_e1e2u(ji,jj) & 
    11041079                  &              * ( e1e2t(ji  ,jj) * ssha(ji  ,jj)    & 
    11051080                  &              +   e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 
    1106                zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1)  * r1_e1e2v(ji,jj) & 
     1081               zsshv_a(ji,jj) = r1_2 * vmask(ji,jj,1)  * r1_e1e2v(ji,jj) & 
    11071082                  &              * ( e1e2t(ji,jj  ) * ssha(ji,jj  )    & 
    11081083                  &              +   e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 
     
    12991274      INTEGER  ::   ji ,jj              ! dummy loop indices 
    13001275      REAL(wp) ::   zxr2, zyr2, zcmax   ! local scalar 
    1301       REAL(wp), POINTER, DIMENSION(:,:) ::   zcu 
     1276      REAL(wp), DIMENSION(jpi,jpj) ::   zcu 
    13021277      !!---------------------------------------------------------------------- 
    13031278      ! 
    13041279      ! Max courant number for ext. grav. waves 
    1305       ! 
    1306       CALL wrk_alloc( jpi,jpj,   zcu ) 
    13071280      ! 
    13081281      DO jj = 1, jpj 
     
    13711344      ENDIF 
    13721345      ! 
    1373       CALL wrk_dealloc( jpi,jpj,   zcu ) 
    1374       ! 
    13751346   END SUBROUTINE dyn_spg_ts_init 
    13761347 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90

    r8093 r8143  
    99   !!            3.3  !  2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps 
    1010   !!            3.4  !  2012-01  (H. Liu) Semi-implicit bottom friction 
     11   !!            4.0  !  2017-05  (G. Madec)  drag coef. defined at t-point (zdfdrg.F90) 
    1112   !!---------------------------------------------------------------------- 
    1213 
     
    2223   USE dynadv   , ONLY: ln_dynadv_vec ! Momentum advection form 
    2324   USE zdf_oce        ! ocean vertical physics 
    24 !!gm new 
    2525   USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    26 !!gm old 
    27    USE zdfbfr         ! Bottom friction setup 
    28 !!gm 
    2926   ! 
    3027   USE in_out_manager ! I/O manager 
    3128   USE lib_mpp        ! MPP library 
    32    USE wrk_nemo       ! Memory Allocation 
    3329   USE timing         ! Timing 
    3430 
     
    6561      !!      with the following surface/top/bottom boundary condition: 
    6662      !!      surface: wind stress input (averaged over kt-1/2 & kt+1/2) 
    67       !!      top & bottom : top stress (iceshelf-ocean) & bottom stress (cf zdfbfr.F) 
     63      !!      top & bottom : top stress (iceshelf-ocean) & bottom stress (cf zdfdrg.F90) 
    6864      !! 
    6965      !! ** Action :   (ua,va) after velocity  
     
    7672      REAL(wp) ::   zzwi, ze3ua   ! local scalars 
    7773      REAL(wp) ::   zzws, ze3va   !   -      - 
    78       REAL(wp), POINTER, DIMENSION(:,:,:) ::  zwi, zwd, zws 
     74      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zwi, zwd, zws 
    7975      !!---------------------------------------------------------------------- 
    8076      ! 
    8177      IF( nn_timing == 1 )  CALL timing_start('dyn_zdf_imp') 
    82       ! 
    83       CALL wrk_alloc( jpi,jpj,jpk, zwi, zwd, zws )  
    8478      ! 
    8579      IF( kt == nit000 ) THEN 
     
    115109      ! column vector of the tri-diagonal matrix equation 
    116110      ! 
    117       IF( ln_bfrimp ) THEN 
     111      IF( ln_drgimp ) THEN 
    118112         DO jj = 2, jpjm1 
    119113            DO ji = 2, jpim1 
    120114               ikbu = mbku(ji,jj)       ! ocean bottom level at u- and v-points  
    121115               ikbv = mbkv(ji,jj)       ! (deepest ocean u- and v-points) 
    122 !!gm old 
    123                avmu(ji,jj,ikbu+1) = -bfrua(ji,jj) * e3uw_n(ji,jj,ikbu+1) 
    124                avmv(ji,jj,ikbv+1) = -bfrva(ji,jj) * e3vw_n(ji,jj,ikbv+1) 
    125 !!gm new 
    126 !               avmu(ji,jj,ikbu+1) = -0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * e3uw_n(ji,jj,ikbu+1) 
    127 !               avmv(ji,jj,ikbv+1) = -0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * e3vw_n(ji,jj,ikbv+1) 
    128 !!gm 
     116               avmu(ji,jj,ikbu+1) = -0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * e3uw_n(ji,jj,ikbu+1) 
     117               avmv(ji,jj,ikbv+1) = -0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * e3vw_n(ji,jj,ikbv+1) 
    129118            END DO 
    130119         END DO 
     
    134123                  ikbu = miku(ji,jj)       ! ocean top level at u- and v-points  
    135124                  ikbv = mikv(ji,jj)       ! (first wet ocean u- and v-points) 
    136 !!gm old 
    137                   IF( ikbu >= 2 )   avmu(ji,jj,ikbu) = -tfrua(ji,jj) * e3uw_n(ji,jj,ikbu) 
    138                   IF( ikbv >= 2 )   avmv(ji,jj,ikbv) = -tfrva(ji,jj) * e3vw_n(ji,jj,ikbv) 
    139 !!gm new 
    140125                  ! top Cd is masked (=0 outside cavities) no need of test on mik>=2 
    141 !                  avmu(ji,jj,ikbu) = -0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * e3uw_n(ji,jj,ikbu) 
    142 !                  avmv(ji,jj,ikbv) = -0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * e3vw_n(ji,jj,ikbv) 
    143 !!gm 
     126                  avmu(ji,jj,ikbu) = -0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * e3uw_n(ji,jj,ikbu) 
     127                  avmv(ji,jj,ikbv) = -0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * e3vw_n(ji,jj,ikbv) 
    144128               END DO 
    145129            END DO 
     
    152136      !            not lead to the effective stress seen over the whole barotropic loop.  
    153137      ! G. Madec : in linear free surface, e3u_a = e3u_n = e3u_0, so systematic use of e3u_a 
    154       IF( ln_bfrimp .AND. ln_dynspg_ts ) THEN 
     138      IF( ln_drgimp .AND. ln_dynspg_ts ) THEN 
    155139         DO jk = 1, jpkm1        ! remove barotropic velocities 
    156140            ua(:,:,jk) = ( ua(:,:,jk) - ua_b(:,:) ) * umask(:,:,jk) 
     
    163147               ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,ikbu) + r_vvl * e3u_a(ji,jj,ikbu) 
    164148               ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikbv) + r_vvl * e3v_a(ji,jj,ikbv) 
    165 !!gm old 
    166                ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * bfrua(ji,jj) * ua_b(ji,jj) / ze3ua 
    167                va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * bfrva(ji,jj) * va_b(ji,jj) / ze3va 
    168 !!gm new 
    169 !               ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * ua_b(ji,jj) / ze3ua 
    170 !               va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * va_b(ji,jj) / ze3va 
    171 !!gm 
     149               ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * ua_b(ji,jj) / ze3ua 
     150               va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * va_b(ji,jj) / ze3va 
    172151            END DO 
    173152         END DO 
     
    179158                  ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,ikbu) + r_vvl * e3u_a(ji,jj,ikbu) 
    180159                  ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikbv) + r_vvl * e3v_a(ji,jj,ikbv) 
    181 !!gm old 
    182                   ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * tfrua(ji,jj) * ua_b(ji,jj) / ze3ua 
    183                   va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * tfrva(ji,jj) * va_b(ji,jj) / ze3va 
    184 !!gm new 
    185160                  ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * ua_b(ji,jj) / ze3ua 
    186161                  va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * va_b(ji,jj) / ze3va 
    187 !!gm 
    188162               END DO 
    189163            END DO 
     
    342316         END DO 
    343317      END DO 
    344        
    345       ! J. Chanut: Lines below are useless ? 
    346       !! restore bottom layer avmu(v)  
    347       !!gm  I almost sure it is !!!! 
    348       IF( ln_bfrimp ) THEN 
    349         DO jj = 2, jpjm1 
    350            DO ji = 2, jpim1 
    351               ikbu = mbku(ji,jj)         ! ocean bottom level at u- and v-points  
    352               ikbv = mbkv(ji,jj)         ! (deepest ocean u- and v-points) 
    353               avmu(ji,jj,ikbu+1) = 0._wp 
    354               avmv(ji,jj,ikbv+1) = 0._wp 
    355            END DO 
    356         END DO 
    357         IF (ln_isfcav) THEN 
    358            DO jj = 2, jpjm1 
    359               DO ji = 2, jpim1 
    360                  ikbu = miku(ji,jj)         ! ocean top level at u- and v-points  
    361                  ikbv = mikv(ji,jj)         ! (first wet ocean u- and v-points) 
    362                  IF( ikbu > 1 )   avmu(ji,jj,ikbu) = 0._wp 
    363                  IF( ikbv > 1 )   avmv(ji,jj,ikbv) = 0._wp 
    364               END DO 
    365            END DO 
    366         ENDIF 
    367       ENDIF 
    368       ! 
    369       CALL wrk_dealloc( jpi,jpj,jpk,   zwi, zwd, zws)  
    370318      ! 
    371319      IF( nn_timing == 1 )   CALL timing_stop('dyn_zdf_imp') 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90

    r7816 r8143  
    55   !!                   shelf 
    66   !!====================================================================== 
    7    !! History :  3.2   !  2011-02  (C.Harris  ) Original code isf cav 
    8    !!            X.X   !  2006-02  (C. Wang   ) Original code bg03 
    9    !!            3.4   !  2013-03  (P. Mathiot) Merging + parametrization 
     7   !! History :  3.2  !  2011-02  (C.Harris  ) Original code isf cav 
     8   !!            X.X  !  2006-02  (C. Wang   ) Original code bg03 
     9   !!            3.4  !  2013-03  (P. Mathiot) Merging + parametrization 
    1010   !!---------------------------------------------------------------------- 
    1111 
    1212   !!---------------------------------------------------------------------- 
    13    !!   sbc_isf        : update sbc under ice shelf 
     13   !!   sbc_isf       : update sbc under ice shelf 
    1414   !!---------------------------------------------------------------------- 
    15    USE oce             ! ocean dynamics and tracers 
    16    USE dom_oce         ! ocean space and time domain 
    17    USE phycst          ! physical constants 
    18    USE eosbn2          ! equation of state 
    19    USE sbc_oce         ! surface boundary condition: ocean fields 
    20    USE zdfbfr          ! 
     15   USE oce            ! ocean dynamics and tracers 
     16   USE dom_oce        ! ocean space and time domain 
     17   USE phycst         ! physical constants 
     18   USE eosbn2         ! equation of state 
     19   USE sbc_oce        ! surface boundary condition: ocean fields 
     20   USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    2121   ! 
    22    USE in_out_manager  ! I/O manager 
    23    USE iom             ! I/O manager library 
    24    USE fldread         ! read input field at current time step 
    25    USE lbclnk          ! 
    26    USE wrk_nemo        ! Memory allocation 
    27    USE timing          ! Timing 
    28    USE lib_fortran     ! glob_sum 
     22   USE in_out_manager ! I/O manager 
     23   USE iom            ! I/O manager library 
     24   USE fldread        ! read input field at current time step 
     25   USE lbclnk         ! 
     26   USE wrk_nemo       ! Memory allocation 
     27   USE timing         ! Timing 
     28   USE lib_fortran    ! glob_sum 
    2929 
    3030   IMPLICIT NONE 
     
    7777CONTAINS 
    7878  
    79   SUBROUTINE sbc_isf(kt) 
     79  SUBROUTINE sbc_isf( kt ) 
    8080      !!--------------------------------------------------------------------- 
    8181      !!                  ***  ROUTINE sbc_isf  *** 
     
    9494      INTEGER               :: ji, jj, jk           ! loop index 
    9595      INTEGER               :: ikt, ikb             ! loop index 
    96       REAL(wp), DIMENSION (:,:), POINTER :: zt_frz, zdep ! freezing temperature (zt_frz) at depth (zdep)  
     96      REAL(wp), DIMENSION(jpi,jpj) ::   zt_frz, zdep  ! freezing temperature (zt_frz) at depth (zdep)  
    9797      REAL(wp), DIMENSION(:,:,:), POINTER :: zfwfisf3d, zqhcisf3d, zqlatisf3d 
    9898      REAL(wp), DIMENSION(:,:  ), POINTER :: zqhcisf2d 
     
    100100      ! 
    101101      IF( MOD( kt-1, nn_fsbc) == 0 ) THEN 
    102          ! allocation 
    103          CALL wrk_alloc( jpi,jpj, zt_frz, zdep  ) 
    104102 
    105103         ! compute salt and heat flux 
     
    204202            CALL wrk_dealloc( jpi,jpj,     zqhcisf2d                        ) 
    205203          END IF 
    206           ! deallocation 
    207           CALL wrk_dealloc( jpi,jpj, zt_frz, zdep  ) 
    208204          ! 
    209205        END IF 
     
    254250  END FUNCTION 
    255251 
     252 
    256253  SUBROUTINE sbc_isf_init 
    257254      !!--------------------------------------------------------------------- 
     
    289286 
    290287      IF ( lwp ) WRITE(numout,*) 
    291       IF ( lwp ) WRITE(numout,*) 'sbc_isf: heat flux of the ice shelf' 
    292       IF ( lwp ) WRITE(numout,*) '~~~~~~~~~' 
    293       IF ( lwp ) WRITE(numout,*) 'sbcisf :'  
    294       IF ( lwp ) WRITE(numout,*) '~~~~~~~~' 
     288      IF ( lwp ) WRITE(numout,*) 'sbc_isf_init : heat flux of the ice shelf' 
     289      IF ( lwp ) WRITE(numout,*) '~~~~~~~~~~~' 
    295290      IF ( lwp ) WRITE(numout,*) '        nn_isf      = ', nn_isf 
    296291      IF ( lwp ) WRITE(numout,*) '        nn_isfblk   = ', nn_isfblk 
     
    299294      IF ( lwp ) WRITE(numout,*) '        rn_gammat0  = ', rn_gammat0   
    300295      IF ( lwp ) WRITE(numout,*) '        rn_gammas0  = ', rn_gammas0   
    301       IF ( lwp ) WRITE(numout,*) '        rn_tfri2    = ', rn_tfri2  
     296      IF ( lwp ) WRITE(numout,*) '        rn_Cd0      = ', r_Cdmin_top  
    302297      ! 
    303298      ! Allocate public variable 
     
    305300      ! 
    306301      ! initialisation 
    307       qisf(:,:)        = 0._wp  ; fwfisf  (:,:) = 0._wp 
    308       risf_tsc(:,:,:)  = 0._wp  ; fwfisf_b(:,:) = 0._wp 
     302      qisf    (:,:)    = 0._wp   ;  fwfisf  (:,:) = 0._wp 
     303      risf_tsc(:,:,:)  = 0._wp   ;  fwfisf_b(:,:) = 0._wp 
    309304      ! 
    310305      ! define isf tbl tickness, top and bottom indice 
     
    312307      CASE ( 1 )  
    313308         rhisf_tbl(:,:) = rn_hisf_tbl 
    314          misfkt(:,:)    = mikt(:,:)         ! same indice for bg03 et cav => used in isfdiv 
     309         misfkt   (:,:) = mikt(:,:)         ! same indice for bg03 et cav => used in isfdiv 
    315310 
    316311      CASE ( 2 , 3 ) 
     
    346341            DO jj = 1, jpj 
    347342                ik = 2 
     343!!gm potential bug: use gdepw_0 not _n 
    348344                DO WHILE ( ik <= mbkt(ji,jj) .AND. gdepw_n(ji,jj,ik) < rzisf_tbl(ji,jj) ) ;  ik = ik + 1 ;  END DO 
    349345                misfkt(ji,jj) = ik-1 
     
    354350         ! as in nn_isf == 1 
    355351         rhisf_tbl(:,:) = rn_hisf_tbl 
    356          misfkt(:,:)    = mikt(:,:)         ! same indice for bg03 et cav => used in isfdiv 
     352         misfkt   (:,:) = mikt(:,:)         ! same indice for bg03 et cav => used in isfdiv 
    357353          
    358354         ! load variable used in fldread (use for temporal interpolation of isf fwf forcing) 
     
    377373            ! determine the deepest level influenced by the boundary layer 
    378374            DO jk = ikt+1, mbkt(ji,jj) 
    379                IF ( (SUM(e3t_n(ji,jj,ikt:jk-1)) < rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) ) ikb = jk 
     375               IF( (SUM(e3t_n(ji,jj,ikt:jk-1)) < rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) )  ikb = jk 
    380376            END DO 
    381377            rhisf_tbl(ji,jj) = MIN(rhisf_tbl(ji,jj), SUM(e3t_n(ji,jj,ikt:ikb))) ! limit the tbl to water thickness. 
     
    390386  END SUBROUTINE sbc_isf_init 
    391387 
     388 
    392389  SUBROUTINE sbc_isf_bg03(kt) 
    393390      !!--------------------------------------------------------------------- 
     
    402399      !!         interaction for climate models", Ocean Modelling 5(2003) 157-170. 
    403400      !!         (hereafter BG) 
    404       !! History : 
    405       !!         06-02  (C. Wang) Original code 
     401      !! History :  06-02  (C. Wang) Original code 
    406402      !!---------------------------------------------------------------------- 
    407403      INTEGER, INTENT ( in ) :: kt 
     
    415411      !!---------------------------------------------------------------------- 
    416412 
    417       IF ( nn_timing == 1 ) CALL timing_start('sbc_isf_bg03') 
     413      IF( nn_timing == 1 )  CALL timing_start('sbc_isf_bg03') 
    418414      ! 
    419415      DO ji = 1, jpi 
     
    441437               !add to salinity trend 
    442438            ELSE 
    443                qisf(ji,jj) = 0._wp ; fwfisf(ji,jj) = 0._wp 
     439               qisf(ji,jj) = 0._wp   ;  fwfisf(ji,jj) = 0._wp 
    444440            END IF 
    445441         END DO 
    446442      END DO 
    447443      ! 
    448       IF( nn_timing == 1 )  CALL timing_stop('sbc_isf_bg03') 
     444      IF( nn_timing == 1 )   CALL timing_stop('sbc_isf_bg03') 
    449445      ! 
    450446  END SUBROUTINE sbc_isf_bg03 
     447 
    451448 
    452449  SUBROUTINE sbc_isf_cav( kt ) 
     
    463460      !!                emp, emps  : update freshwater flux below ice shelf 
    464461      !!--------------------------------------------------------------------- 
    465       INTEGER, INTENT(in)          ::   kt         ! ocean time step 
     462      INTEGER, INTENT(in) ::   kt   ! ocean time step 
    466463      ! 
    467464      INTEGER  ::   ji, jj     ! dummy loop indices 
    468465      INTEGER  ::   nit 
     466      LOGICAL  ::   lit 
    469467      REAL(wp) ::   zlamb1, zlamb2, zlamb3 
    470468      REAL(wp) ::   zeps1,zeps2,zeps3,zeps4,zeps6,zeps7 
     
    472470      REAL(wp) ::   zeps = 1.e-20_wp         
    473471      REAL(wp) ::   zerr 
    474       REAL(wp), DIMENSION(:,:), POINTER ::   zfrz 
    475       REAL(wp), DIMENSION(:,:), POINTER ::   zgammat, zgammas  
    476       REAL(wp), DIMENSION(:,:), POINTER ::   zfwflx, zhtflx, zhtflx_b 
    477       LOGICAL  ::   lit 
     472      REAL(wp), DIMENSION(jpi,jpj) ::   zfrz 
     473      REAL(wp), DIMENSION(jpi,jpj) ::   zgammat, zgammas  
     474      REAL(wp), DIMENSION(jpi,jpj) ::   zfwflx, zhtflx, zhtflx_b 
    478475      !!--------------------------------------------------------------------- 
    479476      ! coeficient for linearisation of potential tfreez 
     
    484481      IF( nn_timing == 1 )  CALL timing_start('sbc_isf_cav') 
    485482      ! 
    486       CALL wrk_alloc( jpi,jpj, zfrz  , zgammat, zgammas  ) 
    487       CALL wrk_alloc( jpi,jpj, zfwflx, zhtflx , zhtflx_b ) 
    488  
    489483      ! initialisation 
    490484      zgammat(:,:) = rn_gammat0 ; zgammas (:,:) = rn_gammas0 
     
    578572      CALL iom_put('isfgammas', zgammas) 
    579573      !  
    580       CALL wrk_dealloc( jpi,jpj, zfrz  , zgammat, zgammas  ) 
    581       CALL wrk_dealloc( jpi,jpj, zfwflx, zhtflx , zhtflx_b ) 
    582       ! 
    583574      IF( nn_timing == 1 )  CALL timing_stop('sbc_isf_cav') 
    584575      ! 
     
    600591      INTEGER  :: ikt                         
    601592      INTEGER  :: ji, jj                     ! loop index 
    602       REAL(wp), DIMENSION(:,:), POINTER :: zustar           ! U, V at T point and friction velocity 
    603593      REAL(wp) :: zdku, zdkv                 ! U, V shear  
    604594      REAL(wp) :: zPr, zSc, zRc              ! Prandtl, Scmidth and Richardson number  
     
    614604      REAL(wp), PARAMETER :: znu   = 1.95e-6_wp ! kinamatic viscosity of sea water (m2.s-1) 
    615605      REAL(wp), DIMENSION(2) :: zts, zab 
     606      REAL(wp), DIMENSION(jpi,jpj) :: zustar   ! U, V at T point and friction velocity 
    616607      !!--------------------------------------------------------------------- 
    617       CALL wrk_alloc( jpi,jpj, zustar ) 
    618608      ! 
    619609      SELECT CASE ( nn_gammablk ) 
     
    626616         !! Jenkins et al., 2010, JPO, p2298-2312 
    627617         !! Adopted by Asay-Davis et al. (2015) 
    628  
    629          !! compute ustar (eq. 24) 
    630          zustar(:,:) = SQRT( rn_tfri2 * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + rn_tfeb2) ) 
     618!!gm  I don't understand the u* expression in those papers... (see for example zdfglf module) 
     619!!    for me ustar= Cd0 * |U|  not  (Cd0)^1/2 * |U| ....  which is what you can find in Jenkins et al. 
     620 
     621         !! compute ustar (eq. 24)           !! NB: here r_Cdmin_top = rn_Cd0 read in namdrg_top namelist) 
     622         zustar(:,:) = SQRT( r_Cdmin_top * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + r_ke0_top) ) 
    631623 
    632624         !! Compute gammats 
     
    638630         !! as MOL depends of flux and flux depends of MOL, best will be iteration (TO DO) 
    639631         !! compute ustar 
    640          zustar(:,:) = SQRT( rn_tfri2 * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + rn_tfeb2) ) 
     632         zustar(:,:) = SQRT( r_Cdmin_top * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + r_ke0_top) ) 
    641633 
    642634         !! compute Pr and Sc number (can be improved) 
     
    649641 
    650642         !! compute gamma 
    651          DO ji=2,jpi 
    652             DO jj=2,jpj 
     643         DO ji = 2, jpi 
     644            DO jj = 2, jpj 
    653645               ikt = mikt(ji,jj) 
    654646 
    655                IF (zustar(ji,jj) == 0._wp) THEN           ! only for kt = 1 I think 
     647               IF( zustar(ji,jj) == 0._wp ) THEN           ! only for kt = 1 I think 
    656648                  pgt = rn_gammat0 
    657649                  pgs = rn_gammas0 
    658650               ELSE 
    659651                  !! compute Rc number (as done in zdfric.F90) 
     652!!gm better to do it like in the new zdfric.F90   i.e. avm weighted Ri computation 
     653!!gm moreover, use Max(rn2,0) to take care of static instabilities.... 
    660654                  zcoef = 0.5_wp / e3w_n(ji,jj,ikt) 
    661655                  !                                            ! shear of horizontal velocity 
     
    703697         CALL lbc_lnk(pgs(:,:),'T',1.) 
    704698      END SELECT 
    705       CALL wrk_dealloc( jpi,jpj, zustar ) 
    706699      ! 
    707700   END SUBROUTINE sbc_isf_gammats 
    708701 
     702 
    709703   SUBROUTINE sbc_isf_tbl( pvarin, pvarout, cd_ptin ) 
    710704      !!---------------------------------------------------------------------- 
     
    714708      !! 
    715709      !!---------------------------------------------------------------------- 
    716       REAL(wp), DIMENSION(:,:,:), INTENT( in  ) :: pvarin 
    717       REAL(wp), DIMENSION(:,:)  , INTENT( out ) :: pvarout 
    718       CHARACTER(len=1),           INTENT( in  ) :: cd_ptin ! point of variable in/out 
    719       ! 
    720       REAL(wp) :: ze3, zhk 
     710      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) :: pvarin 
     711      REAL(wp), DIMENSION(:,:)  , INTENT(  out) :: pvarout 
     712      CHARACTER(len=1),           INTENT(in   ) :: cd_ptin ! point of variable in/out 
     713      ! 
     714      INTEGER ::   ji, jj, jk                ! loop index 
     715      INTEGER ::   ikt, ikb                    ! top and bottom index of the tbl 
     716      REAL(wp) ::   ze3, zhk 
    721717      REAL(wp), DIMENSION(:,:), POINTER :: zhisf_tbl ! thickness of the tbl 
    722  
    723       INTEGER :: ji, jj, jk                  ! loop index 
    724       INTEGER :: ikt, ikb                    ! top and bottom index of the tbl 
    725718      !!---------------------------------------------------------------------- 
    726719      ! allocation 
     
    736729               ikt = miku(ji,jj) ; ikb = miku(ji,jj) 
    737730               ! thickness of boundary layer at least the top level thickness 
    738                zhisf_tbl(ji,jj) = MAX(rhisf_tbl_0(ji,jj), e3u_n(ji,jj,ikt)) 
     731               zhisf_tbl(ji,jj) = MAX( rhisf_tbl_0(ji,jj) , e3u_n(ji,jj,ikt) ) 
    739732 
    740733               ! determine the deepest level influenced by the boundary layer 
     
    755748            END DO 
    756749         END DO 
    757          DO jj = 2,jpj 
    758             DO ji = 2,jpi 
     750         DO jj = 2, jpj 
     751            DO ji = 2, jpi 
     752!!gm a wet-point only average should be used here !!! 
    759753               pvarout(ji,jj) = 0.5_wp * (pvarout(ji,jj) + pvarout(ji-1,jj)) 
    760754            END DO 
     
    786780            END DO 
    787781         END DO 
    788          DO jj = 2,jpj 
    789             DO ji = 2,jpi 
     782         DO jj = 2, jpj 
     783            DO ji = 2, jpi 
     784!!gm a wet-point only average should be used here !!! 
    790785               pvarout(ji,jj) = 0.5_wp * (pvarout(ji,jj) + pvarout(ji,jj-1)) 
    791786            END DO 
     
    882877      ! 
    883878   END SUBROUTINE sbc_isf_div 
     879 
    884880   !!====================================================================== 
    885881END MODULE sbcisf 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/TRD/trd_oce.F90

    r7646 r8143  
    7272   INTEGER, PUBLIC, PARAMETER ::   jpdyn_atf  = 10     !: Asselin time filter 
    7373   INTEGER, PUBLIC, PARAMETER ::   jpdyn_tau  = 11     !: surface stress 
    74    INTEGER, PUBLIC, PARAMETER ::   jpdyn_bfri = 12     !: implicit bottom friction (ln_bfrimp=.TRUE.) 
     74   INTEGER, PUBLIC, PARAMETER ::   jpdyn_bfri = 12     !: implicit bottom friction (ln_drgimp=.TRUE.) 
    7575   INTEGER, PUBLIC, PARAMETER ::   jpdyn_ken  = 13     !: use for calculation of KE 
    7676   ! 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/TRD/trddyn.F90

    r6140 r8143  
    1515   USE oce            ! ocean dynamics and tracers variables 
    1616   USE dom_oce        ! ocean space and time domain variables 
    17    USE zdf_oce        ! ocean vertical physics variables 
     17   USE phycst         ! physical constants 
     18   USE sbc_oce        ! surface boundary condition: ocean 
     19   USE zdf_oce        ! ocean vertical physics: variables 
     20   USE zdfdrg         ! ocean vertical physics: bottom friction 
    1821   USE trd_oce        ! trends: ocean variables 
    19    USE zdfbfr         ! bottom friction 
    20    USE sbc_oce        ! surface boundary condition: ocean 
    21    USE phycst         ! physical constants 
    2222   USE trdken         ! trends: Kinetic ENergy  
    2323   USE trdglo         ! trends: global domain averaged 
    2424   USE trdvor         ! trends: vertical averaged vorticity  
    2525   USE trdmxl         ! trends: mixed layer averaged  
     26   ! 
    2627   USE in_out_manager ! I/O manager 
    2728   USE lbclnk         ! lateral boundary condition  
    2829   USE iom            ! I/O manager library 
    2930   USE lib_mpp        ! MPP library 
    30    USE wrk_nemo       ! Memory allocation 
    3131 
    3232   IMPLICIT NONE 
    3333   PRIVATE 
    3434 
    35    PUBLIC trd_dyn        ! called by all dynXX modules 
     35   PUBLIC trd_dyn        ! called by all dynXXX modules 
    3636 
    3737   !! * Substitutions 
    3838#  include "vectopt_loop_substitute.h90" 
    3939   !!---------------------------------------------------------------------- 
    40    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     40   !! NEMO/OPA 4.0 , NEMO Consortium (2017) 
    4141   !! $Id$ 
    4242   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    103103      INTEGER ::   ji, jj, jk   ! dummy loop indices 
    104104      INTEGER ::   ikbu, ikbv   ! local integers 
    105       REAL(wp), POINTER, DIMENSION(:,:)   ::   z2dx, z2dy   ! 2D workspace  
    106       REAL(wp), POINTER, DIMENSION(:,:,:) ::   z3dx, z3dy   ! 3D workspace  
     105      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   z2dx, z2dy   ! 2D workspace  
     106      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   z3dx, z3dy   ! 3D workspace  
    107107      !!---------------------------------------------------------------------- 
    108108      ! 
     
    118118      CASE( jpdyn_keg )   ;   CALL iom_put( "utrd_keg", putrd )    ! Kinetic Energy gradient (or had) 
    119119                              CALL iom_put( "vtrd_keg", pvtrd ) 
    120                               CALL wrk_alloc( jpi, jpj, jpk, z3dx, z3dy ) 
     120                              ALLOCATE( z3dx(jpi,jpj,jpk) , z3dy(jpi,jpj,jpk) ) 
    121121                              z3dx(:,:,:) = 0._wp                  ! U.dxU & V.dyV (approximation) 
    122122                              z3dy(:,:,:) = 0._wp 
     
    133133                              CALL iom_put( "utrd_udx", z3dx  ) 
    134134                              CALL iom_put( "vtrd_vdy", z3dy  ) 
    135                               CALL wrk_dealloc( jpi, jpj, jpk, z3dx, z3dy ) 
    136       CASE( jpdyn_zad )   ;   CALL iom_put( "utrd_zad", putrd )    ! vertical   advection 
     135                              DEALLOCATE( z3dx , z3dy ) 
     136      CASE( jpdyn_zad )   ;   CALL iom_put( "utrd_zad", putrd )    ! vertical advection 
    137137                              CALL iom_put( "vtrd_zad", pvtrd ) 
    138       CASE( jpdyn_ldf )   ;   CALL iom_put( "utrd_ldf", putrd )    ! lateral diffusion 
     138      CASE( jpdyn_ldf )   ;   CALL iom_put( "utrd_ldf", putrd )    ! lateral  diffusion 
    139139                              CALL iom_put( "vtrd_ldf", pvtrd ) 
    140140      CASE( jpdyn_zdf )   ;   CALL iom_put( "utrd_zdf", putrd )    ! vertical diffusion  
    141141                              CALL iom_put( "vtrd_zdf", pvtrd ) 
     142                              ! 
    142143                              !                                    ! wind stress trends 
    143                               CALL wrk_alloc( jpi, jpj, z2dx, z2dy ) 
     144                              ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) ) 
    144145                              z2dx(:,:) = ( utau_b(:,:) + utau(:,:) ) / ( e3u_n(:,:,1) * rau0 ) 
    145146                              z2dy(:,:) = ( vtau_b(:,:) + vtau(:,:) ) / ( e3v_n(:,:,1) * rau0 ) 
    146147                              CALL iom_put( "utrd_tau", z2dx ) 
    147148                              CALL iom_put( "vtrd_tau", z2dy ) 
    148                               CALL wrk_dealloc( jpi, jpj, z2dx, z2dy ) 
    149       CASE( jpdyn_bfr )       ! called if ln_bfrimp=T 
    150                               CALL iom_put( "utrd_bfr", putrd )    ! bottom friction (explicit case) 
    151                               CALL iom_put( "vtrd_bfr", pvtrd ) 
    152       CASE( jpdyn_atf )   ;   CALL iom_put( "utrd_atf", putrd )        ! asselin filter trends  
    153                               CALL iom_put( "vtrd_atf", pvtrd ) 
    154       CASE( jpdyn_bfri )  ;   IF( ln_bfrimp ) THEN                     ! bottom friction (implicit case) 
    155                                  CALL wrk_alloc( jpi, jpj, jpk, z3dx, z3dy ) 
     149                              DEALLOCATE( z2dx , z2dy ) 
     150                              !                                    ! bottom stress tends (implicit case) 
     151                              IF( ln_drgimp ) THEN 
     152                                 ALLOCATE( z3dx(jpi,jpj,jpk) , z3dy(jpi,jpj,jpk) ) 
    156153                                 z3dx(:,:,:) = 0._wp   ;   z3dy(:,:,:) = 0._wp  ! after velocity known (now filed at this stage) 
    157154                                 DO jk = 1, jpkm1 
     
    160157                                          ikbu = mbku(ji,jj)          ! deepest ocean u- & v-levels 
    161158                                          ikbv = mbkv(ji,jj) 
    162                                           z3dx(ji,jj,jk) = bfrua(ji,jj) * un(ji,jj,ikbu) / e3u_n(ji,jj,ikbu) 
    163                                           z3dy(ji,jj,jk) = bfrva(ji,jj) * vn(ji,jj,ikbv) / e3v_n(ji,jj,ikbv) 
     159                                          z3dx(ji,jj,jk) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) )*un(ji,jj,ikbu)/e3u_n(ji,jj,ikbu) 
     160                                          z3dy(ji,jj,jk) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) )*vn(ji,jj,ikbv)/e3v_n(ji,jj,ikbv) 
    164161                                       END DO 
    165162                                    END DO 
    166163                                 END DO 
    167                                  CALL lbc_lnk( z3dx, 'U', -1. ) ; CALL lbc_lnk( z3dy, 'V', -1. ) 
    168                                  CALL iom_put( "utrd_bfri", z3dx ) 
    169                                  CALL iom_put( "vtrd_bfri", z3dy ) 
    170                                  CALL wrk_dealloc( jpi, jpj, jpk, z3dx, z3dy ) 
    171                               ENDIF 
     164                                 CALL lbc_lnk( z3dx, 'U', -1. )   ;   CALL lbc_lnk( z3dy, 'V', -1. ) 
     165                                 CALL iom_put( "utrd_bfr", z3dx ) 
     166                                 CALL iom_put( "vtrd_bfr", z3dy ) 
     167                                 DEALLOCATE( z3dx , z3dy ) 
     168                              ENDIF 
     169      CASE( jpdyn_bfr )       ! called if ln_drgimp=F 
     170                              CALL iom_put( "utrd_bfr", putrd )    ! bottom friction (explicit case) 
     171                              CALL iom_put( "vtrd_bfr", pvtrd ) 
     172      CASE( jpdyn_atf )   ;   CALL iom_put( "utrd_atf", putrd )        ! asselin filter trends  
     173                              CALL iom_put( "vtrd_atf", pvtrd ) 
    172174      END SELECT 
    173175      ! 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/TRD/trdglo.F90

    r7931 r8143  
    99 
    1010   !!---------------------------------------------------------------------- 
    11    !!   trd_glo      : domain averaged budget of trends (including kinetic energy and T^2 trends) 
    12    !!   glo_dyn_wri  : print dynamic trends in ocean.output file 
    13    !!   glo_tra_wri  : print global T & T^2 trends in ocean.output file 
    14    !!   trd_glo_init : initialization step 
     11   !!   trd_glo       : domain averaged budget of trends (including kinetic energy and T^2 trends) 
     12   !!   glo_dyn_wri   : print dynamic trends in ocean.output file 
     13   !!   glo_tra_wri   : print global T & T^2 trends in ocean.output file 
     14   !!   trd_glo_init  : initialization step 
    1515   !!---------------------------------------------------------------------- 
    16    USE oce             ! ocean dynamics and tracers variables 
    17    USE dom_oce         ! ocean space and time domain variables 
    18    USE sbc_oce         ! surface boundary condition: ocean 
    19    USE trd_oce         ! trends: ocean variables 
    20    USE phycst          ! physical constants 
    21    USE ldftra          ! lateral diffusion: eddy diffusivity & EIV coeff. 
    22    USE ldfdyn          ! ocean dynamics: lateral physics 
    23    USE zdf_oce         ! ocean vertical physics 
    24    USE zdfbfr          ! bottom friction 
    25    USE zdfddm          ! ocean vertical physics: double diffusion 
    26    USE eosbn2          ! equation of state 
    27    USE phycst          ! physical constants 
     16   USE oce            ! ocean dynamics and tracers variables 
     17   USE dom_oce        ! ocean space and time domain variables 
     18   USE sbc_oce        ! surface boundary condition: ocean 
     19   USE trd_oce        ! trends: ocean variables 
     20   USE phycst         ! physical constants 
     21   USE ldftra         ! lateral diffusion: eddy diffusivity & EIV coeff. 
     22   USE ldfdyn         ! ocean dynamics: lateral physics 
     23   USE zdf_oce        ! ocean vertical physics 
     24   USE zdfdrg         ! ocean vertical physics: bottom friction 
     25   USE zdfddm         ! ocean vertical physics: double diffusion 
     26   USE eosbn2         ! equation of state 
     27   USE phycst         ! physical constants 
    2828   ! 
    29    USE lib_mpp         ! distibuted memory computing library 
    30    USE in_out_manager  ! I/O manager 
    31    USE iom             ! I/O manager library 
    32    USE wrk_nemo        ! Memory allocation 
     29   USE lib_mpp        ! distibuted memory computing library 
     30   USE in_out_manager ! I/O manager 
     31   USE iom            ! I/O manager library 
    3332 
    3433   IMPLICIT NONE 
     
    7776      INTEGER ::   ikbu, ikbv      ! local integers 
    7877      REAL(wp)::   zvm, zvt, zvs, z1_2rau0   ! local scalars 
    79       REAL(wp), POINTER, DIMENSION(:,:)  :: ztswu, ztswv, z2dx, z2dy   ! 2D workspace  
    80       !!---------------------------------------------------------------------- 
    81  
    82       CALL wrk_alloc( jpi, jpj, ztswu, ztswv, z2dx, z2dy ) 
    83  
     78      REAL(wp), DIMENSION(jpi,jpj)  :: ztswu, ztswv, z2dx, z2dy   ! 2D workspace  
     79      !!---------------------------------------------------------------------- 
     80      ! 
    8481      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN 
    8582         ! 
     
    123120               DO jj = 1, jpjm1 
    124121                  DO ji = 1, jpim1 
    125                      zvt = ptrdx(ji,jj,jk) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)   & 
    126                         &                  * e1u    (ji  ,jj  ) * e2u    (ji,jj) * e3u_n(ji,jj,jk) 
    127                      zvs = ptrdy(ji,jj,jk) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)   & 
    128                         &                  * e1v    (ji  ,jj  ) * e2v    (ji,jj) * e3u_n(ji,jj,jk) 
     122                     zvt = ptrdx(ji,jj,jk) * tmask_i(ji+1,jj) * tmask_i(ji,jj) * umask(ji,jj,jk)   & 
     123                        &                                     * e1e2u  (ji,jj) * e3u_n(ji,jj,jk) 
     124                     zvs = ptrdy(ji,jj,jk) * tmask_i(ji,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)   & 
     125                        &                                     * e1e2v  (ji,jj) * e3u_n(ji,jj,jk) 
    129126                     umo(ktrd) = umo(ktrd) + zvt 
    130127                     vmo(ktrd) = vmo(ktrd) + zvs 
     
    138135               DO jj = 1, jpjm1 
    139136                  DO ji = 1, jpim1 
    140                      zvt = ( utau_b(ji,jj) + utau(ji,jj) ) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)   & 
    141                         &                       * z1_2rau0 * e1u    (ji  ,jj  ) * e2u    (ji,jj) 
    142                      zvs = ( vtau_b(ji,jj) + vtau(ji,jj) ) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)   & 
    143                         &                       * z1_2rau0 * e1v    (ji  ,jj  ) * e2v    (ji,jj) * e3u_n(ji,jj,jk) 
     137                     zvt = ( utau_b(ji,jj) + utau(ji,jj) ) * tmask_i(ji+1,jj) * tmask_i(ji,jj) * umask(ji,jj,jk)   & 
     138                        &                                                     * z1_2rau0       * e1e2u(ji,jj) 
     139                     zvs = ( vtau_b(ji,jj) + vtau(ji,jj) ) * tmask_i(ji,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)   & 
     140                        &                                                     * z1_2rau0       * e1e2v(ji,jj) 
    144141                     umo(jpdyn_tau) = umo(jpdyn_tau) + zvt 
    145142                     vmo(jpdyn_tau) = vmo(jpdyn_tau) + zvs 
     
    151148            IF( ktrd == jpdyn_atf ) THEN     ! last trend (asselin time filter) 
    152149               ! 
    153                IF( ln_bfrimp ) THEN                   ! implicit bfr case: compute separately the bottom friction  
     150               IF( ln_drgimp ) THEN                   ! implicit drag case: compute separately the bottom friction  
    154151                  z1_2rau0 = 0.5_wp / rau0 
    155152                  DO jj = 1, jpjm1 
     
    157154                        ikbu = mbku(ji,jj)                  ! deepest ocean u- & v-levels 
    158155                        ikbv = mbkv(ji,jj) 
    159                         zvt = bfrua(ji,jj) * un(ji,jj,ikbu) * e1u(ji,jj) * e2v(ji,jj) 
    160                         zvs = bfrva(ji,jj) * vn(ji,jj,ikbv) * e1v(ji,jj) * e2v(ji,jj) 
     156                        zvt = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * un(ji,jj,ikbu) * e1e2u(ji,jj) 
     157                        zvs = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vn(ji,jj,ikbv) * e1e2v(ji,jj) 
    161158                        umo(jpdyn_bfri) = umo(jpdyn_bfri) + zvt 
    162159                        vmo(jpdyn_bfri) = vmo(jpdyn_bfri) + zvs 
     
    165162                  END DO 
    166163               ENDIF 
     164!!gm top drag case is missing  
    167165               !  
    168166               CALL glo_dyn_wri( kt )                 ! print the results in ocean.output 
     
    178176      ENDIF 
    179177      ! 
    180       CALL wrk_dealloc( jpi, jpj, ztswu, ztswv, z2dx, z2dy ) 
    181       ! 
    182178   END SUBROUTINE trd_glo 
    183179 
     
    193189      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    194190      REAL(wp) ::   zcof         ! local scalar 
    195       REAL(wp), POINTER, DIMENSION(:,:,:)  ::  zkx, zky, zkz, zkepe   
    196       !!---------------------------------------------------------------------- 
    197  
    198       CALL wrk_alloc( jpi, jpj, jpk, zkx, zky, zkz, zkepe ) 
     191      REAL(wp), DIMENSION(jpi,jpj,jpk)  ::  zkx, zky, zkz, zkepe   
     192      !!---------------------------------------------------------------------- 
    199193 
    200194      ! I. Momentum trends 
     
    283277            &      + vmo(jpdyn_bfr) + vmo(jpdyn_atf) ) / tvolv 
    284278            WRITE (numout,9513) umo(jpdyn_tau) / tvolu, vmo(jpdyn_tau) / tvolv 
    285             IF( ln_bfrimp )   WRITE (numout,9514) umo(jpdyn_bfri) / tvolu, vmo(jpdyn_bfri) / tvolv 
     279            IF( ln_drgimp )   WRITE (numout,9514) umo(jpdyn_bfri) / tvolu, vmo(jpdyn_bfri) / tvolv 
    286280         ENDIF 
    287281 
     
    322316            &      + hke(jpdyn_bfr) + hke(jpdyn_atf) ) / tvolt 
    323317            WRITE (numout,9533) hke(jpdyn_tau) / tvolt 
    324             IF( ln_bfrimp )   WRITE (numout,9534) hke(jpdyn_bfri) / tvolt 
     318            IF( ln_drgimp )   WRITE (numout,9534) hke(jpdyn_bfri) / tvolt 
    325319         ENDIF 
    326320 
     
    372366      ENDIF 
    373367      ! 
    374       CALL wrk_dealloc( jpi, jpj, jpk, zkx, zky, zkz, zkepe ) 
    375       ! 
    376368   END SUBROUTINE glo_dyn_wri 
    377369 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90

    r7646 r8143  
    1313   USE oce            ! ocean dynamics and tracers variables 
    1414   USE dom_oce        ! ocean space and time domain variables 
     15   USE phycst         ! physical constants 
    1516   USE sbc_oce        ! surface boundary condition: ocean 
    1617   USE zdf_oce        ! ocean vertical physics variables 
     18   USE zdfdrg         ! ocean vertical physics: bottom friction 
     19!!gm   USE dynhpg          ! hydrostatic pressure gradient    
     20   USE ldftra         ! ocean active tracers lateral physics 
    1721   USE trd_oce        ! trends: ocean variables 
    18 !!gm   USE dynhpg          ! hydrostatic pressure gradient    
    19    USE zdfbfr         ! bottom friction 
    20    USE ldftra         ! ocean active tracers lateral physics 
    21    USE phycst         ! physical constants 
    2222   USE trdvor         ! ocean vorticity trends  
    2323   USE trdglo         ! trends:global domain averaged 
     
    2727   USE iom            ! I/O manager library 
    2828   USE lib_mpp        ! MPP library 
    29    USE wrk_nemo       ! Memory allocation 
    3029   USE ldfslp         ! Isopycnal slopes 
    3130 
     
    7473      !!          diagnose separately the KE trend associated with wind stress 
    7574      !!              - bottom friction case (jpdyn_bfr): 
    76       !!          explicit case (ln_bfrimp=F): bottom trend put in the 1st level  
     75      !!          explicit case (ln_drgimp=F): bottom trend put in the 1st level  
    7776      !!                                       of putrd, pvtrd 
    7877      ! 
     
    8685      INTEGER ::   ikbu  , ikbv     ! local integers 
    8786      INTEGER ::   ikbum1, ikbvm1   !   -       - 
    88       REAL(wp), POINTER, DIMENSION(:,:)   ::   z2dx, z2dy, zke2d   ! 2D workspace  
    89       REAL(wp), POINTER, DIMENSION(:,:,:) ::   zke                 ! 3D workspace  
    90       !!---------------------------------------------------------------------- 
    91       ! 
    92       CALL wrk_alloc( jpi, jpj, jpk, zke ) 
     87      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::   z2dx, z2dy, zke2d   ! 2D workspace  
     88      REAL(wp), DIMENSION(jpi,jpj,jpk)      ::   zke                 ! 3D workspace  
     89      !!---------------------------------------------------------------------- 
    9390      ! 
    9491      CALL lbc_lnk( putrd, 'U', -1. )   ;   CALL lbc_lnk( pvtrd, 'V', -1. )      ! lateral boundary conditions 
     
    125122         CASE( jpdyn_zdf )   ;   CALL iom_put( "ketrd_zdf"   , zke )    ! vertical diffusion  
    126123         !                   !                                          ! wind stress trends 
    127                                  CALL wrk_alloc( jpi, jpj, z2dx, z2dy, zke2d ) 
     124                                 ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) , zke2d(jpi,jpj) ) 
    128125                           z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1e2u(:,:) * umask(:,:,1) 
    129126                           z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1e2v(:,:) * vmask(:,:,1) 
     
    136133                           END DO 
    137134                                 CALL iom_put( "ketrd_tau"   , zke2d )  !  
    138                                  CALL wrk_dealloc( jpi, jpj     , z2dx, z2dy, zke2d ) 
     135                                 DEALLOCATE( z2dx , z2dy , zke2d ) 
    139136         CASE( jpdyn_bfr )   ;   CALL iom_put( "ketrd_bfr"   , zke )    ! bottom friction (explicit case)  
    140137!!gm TO BE DONE properly 
    141 !!gm only valid if ln_bfrimp=F otherwise the bottom stress as to be recomputed at the end of the computation.... 
    142 !         IF(.NOT. ln_bfrimp) THEN 
     138!!gm only valid if ln_drgimp=F otherwise the bottom stress as to be recomputed at the end of the computation.... 
     139!         IF(.NOT. ln_drgimp) THEN 
    143140!            DO jj = 1, jpj    !    
    144141!               DO ji = 1, jpi 
     
    163160!! reflechir a une possible sauvegarde du "vrai" un,vn pour le calcul de atf.... 
    164161! 
    165 !         IF( ln_bfrimp ) THEN                                          ! bottom friction (implicit case) 
     162!         IF( ln_drgimp ) THEN                                          ! bottom friction (implicit case) 
    166163!            DO jj = 1, jpj                                                  ! after velocity known (now filed at this stage) 
    167164!               DO ji = 1, jpi 
     
    192189      END SELECT 
    193190      ! 
    194       CALL wrk_dealloc( jpi, jpj, jpk, zke ) 
    195       ! 
    196191   END SUBROUTINE trd_ken 
    197192 
     
    207202      !! ** Work only for full steps and partial steps (ln_hpg_zco or ln_hpg_zps) 
    208203      !!----------------------------------------------------------------------  
    209       INTEGER, INTENT(in) ::   kt    ! ocean time-step index 
    210       !! 
    211       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(  out) ::   pconv 
    212       ! 
    213       INTEGER  ::   ji, jj, jk                       ! dummy loop indices 
    214       INTEGER  ::   iku, ikv                         ! temporary integers 
    215       REAL(wp) ::   zcoef                            ! temporary scalars 
    216       REAL(wp), POINTER, DIMENSION(:,:,:) ::  zconv  ! temporary conv on W-grid 
    217       !!---------------------------------------------------------------------- 
    218       ! 
    219       CALL wrk_alloc( jpi,jpj,jpk, zconv ) 
     204      INTEGER                   , INTENT(in   ) ::   kt      ! ocean time-step index 
     205      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pconv   !  
     206      ! 
     207      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     208      INTEGER  ::   iku, ikv     ! local integers 
     209      REAL(wp) ::   zcoef        ! local scalars 
     210      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zconv  ! 3D workspace 
     211      !!---------------------------------------------------------------------- 
    220212      ! 
    221213      ! Local constant initialization  
     
    240232      END DO 
    241233      ! 
    242       CALL wrk_dealloc( jpi,jpj,jpk, zconv )       
    243       ! 
    244234   END SUBROUTINE ken_p2k 
    245235 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90

    r8093 r8143  
    44   !! Ocean physics : define vertical mixing variables 
    55   !!===================================================================== 
    6    !! history :  1.0  !  2002-06  (G. Madec) Original code 
    7    !!            3.2  !  2009-07  (G.Madec) addition of avm 
     6   !! history :  1.0  !  2002-06  (G. Madec)  Original code 
     7   !!            3.2  !  2009-07  (G. Madec)  addition of avm 
     8   !!            4.0  !  2017-05  (G. Madec)  drag coef. defined at t-point (zdfdrg.F90) 
    89   !!---------------------------------------------------------------------- 
    910   USE par_oce        ! ocean parameters 
     
    5455   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   avtb_2d        !: horizontal shape of background Kz profile 
    5556   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:)     ::   avmb , avtb    !: background profile of avm and avt 
    56 !!gm 
    57    REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   bfrua, bfrva   !: bottom friction coefficients 
    58    REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   tfrua, tfrva   !: top    friction coefficients 
    59 !!gm 
    6057 
    6158   !!---------------------------------------------------------------------- 
     
    7168      !!---------------------------------------------------------------------- 
    7269      ! 
    73       ALLOCATE( avm  (jpi,jpj,jpk) , avt  (jpi,jpj,jpk) , avs(jpi,jpj,jpk) ,        & 
    74          &      avm_k(jpi,jpj,jpk) , avt_k(jpi,jpj,jpk) , en (jpi,jpj,jpk) ,        &  
    75          &      avmb(jpk) , bfrua(jpi,jpj) , tfrua(jpi, jpj) ,                      & 
    76          &      avtb(jpk) , bfrva(jpi,jpj) , tfrva(jpi, jpj) , avtb_2d(jpi,jpj) ,   & 
    77          &      avmu(jpi,jpj,jpk), avmv(jpi,jpj,jpk) ,  STAT = zdf_oce_alloc ) 
     70      ALLOCATE( avm  (jpi,jpj,jpk) , avt  (jpi,jpj,jpk) , avs(jpi,jpj,jpk) ,   & 
     71         &      avm_k(jpi,jpj,jpk) , avt_k(jpi,jpj,jpk) , en (jpi,jpj,jpk) ,   &  
     72         &      avmb(jpk) , avtb(jpk) ,  avtb_2d(jpi,jpj) ,                    & 
     73         &      avmu(jpi,jpj,jpk), avmv(jpi,jpj,jpk)    , STAT = zdf_oce_alloc ) 
    7874         ! 
    7975      IF( zdf_oce_alloc /= 0 )   CALL ctl_warn('zdf_oce_alloc: failed to allocate arrays') 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfdrg.F90

    r8093 r8143  
    1010   !!            3.4  ! 2011-11  (H. Liu) implementation of semi-implicit bottom friction option 
    1111   !!                 ! 2012-06  (H. Liu) implementation of Log Layer bottom friction option 
     12   !!            4.0  ! 2017-05  (G. Madec) zdfbfr becomes zdfdrg + variable names change 
     13   !!                                     + drag defined at t-point + new user interface + top drag (ocean cavities) 
    1214   !!---------------------------------------------------------------------- 
    1315 
     
    2830   USE prtctl         ! Print control 
    2931   USE timing         ! Timing 
    30    USE wrk_nemo       ! Memory Allocation 
    3132 
    3233   IMPLICIT NONE 
     
    137138      ENDIF 
    138139      ! 
    139 !!gm to be moved at the end of zdfphy 
    140       CALL lbc_lnk( pCdU, 'T', 1. )    ! Lateral boundary condition 
    141 !!gm end 
    142       ! 
    143140      IF(ln_ctl)   CALL prt_ctl( tab2d_1=pCdU, clinfo1=' Cd*U ') 
    144141      ! 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r8093 r8143  
    88   !!            3.3  !  2010-10  (C. Bricaud)  Add in the reference 
    99   !!            4.0  !  2017-04  (G. Madec)  remove CPP keys & avm at t-point only  
     10   !!             -   !  2017-05  (G. Madec)  add top friction as boundary condition 
    1011   !!---------------------------------------------------------------------- 
    1112 
     
    1920   USE domvvl         ! ocean space and time domain : variable volume layer 
    2021   USE zdf_oce        ! ocean vertical physics 
    21 !!gm old 
    22    USE zdfbfr  , ONLY : rn_tfrz0, rn_bfrz0   ! top/bottom roughness 
    23 !!gm new 
    2422   USE zdfdrg  , ONLY : r_z0_top , r_z0_bot   ! top/bottom roughness 
    2523   USE zdfdrg  , ONLY : rCdU_top , rCdU_bot   ! top/bottom friction 
    26 !!gm 
    2724   USE sbc_oce        ! surface boundary condition: ocean 
    2825   USE phycst         ! physical constants 
     
    3229   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    3330   USE lib_mpp        ! MPP manager 
    34    USE wrk_nemo       ! work arrays 
    3531   USE prtctl         ! Print control 
    3632   USE in_out_manager ! I/O manager 
     
    149145      REAL(wp) ::   prod, buoy, diss, zdiss, sm         !   -      - 
    150146      REAL(wp) ::   gh, gm, shr, dif, zsqen, zavt, zavm !   -      - 
     147      REAL(wp) ::   zmsku, zmskv                        !   -      - 
    151148      REAL(wp), DIMENSION(jpi,jpj)     ::   zdep 
    152149      REAL(wp), DIMENSION(jpi,jpj)     ::   zkar 
     
    158155      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zwall_psi   ! Wall function use in the wb case (ln_sigpsi) 
    159156      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   psi         ! psi at time now 
    160       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z_elem_a    ! element of the first  matrix diagonal 
    161       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z_elem_b    ! element of the second matrix diagonal 
    162       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z_elem_c    ! element of the third  matrix diagonal 
     157      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zd_lw, zd_up, zdiag   ! lower, upper  and diagonal of the matrix 
    163158      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zstt, zstm  ! stability function on tracer and momentum 
    164159      !!-------------------------------------------------------------------- 
     
    171166      ustar2_top (:,:) = 0._wp   ;   zwall_psi(:,:,:) = 0._wp 
    172167      ustar2_bot (:,:) = 0._wp 
    173        
    174  
    175168 
    176169      ! Compute surface, top and bottom friction at T-points 
     
    181174            ustar2_surf(ji,jj) = r1_rau0 * taum(ji,jj) * tmask(ji,jj,1) 
    182175            !    
    183 !!gm old 
    184             ! bottom friction (explicit before friction)         
    185             ! Note that we chose here not to bound the friction as in dynbfr)    
    186             ztx2 = (  bfrua(ji,jj)  * ub(ji,jj,mbku(ji,jj)) + bfrua(ji-1,jj) * ub(ji-1,jj,mbku(ji-1,jj))  )   &          
    187                & * ( 1._wp - 0.5_wp * umask(ji,jj,1) * umask(ji-1,jj,1)  )       
    188             zty2 = (  bfrva(ji,jj)  * vb(ji,jj,mbkv(ji,jj)) + bfrva(ji,jj-1) * vb(ji,jj-1,mbkv(ji,jj-1))  )   &          
    189                & * ( 1._wp - 0.5_wp * vmask(ji,jj,1) * vmask(ji,jj-1,1)  )       
    190             ustar2_bot(ji,jj) = SQRT( ztx2 * ztx2 + zty2 * zty2 ) * tmask(ji,jj,1)          
    191          END DO          
    192       END DO     
    193 !!gm new 
    194176!!gm Rq we may add here r_ke0(_top/_bot) ?  ==>> think about that... 
    195 !          ! bottom friction (explicit before friction) 
    196 !          zmsku = ( 2._wp - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) 
    197 !          zmskv = ( 2._wp - vmask(ji,jj-1,mbkt(ji,jj)) * vmask(ji,jj,mbkt(ji,jj)) )     ! (CAUTION: CdU<0) 
    198 !          ustar2_bot(ji,jj) = - rCdU_bot(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mbkt(ji,jj))+ub(ji-1,jj,mbkt(ji,jj)) ) )**2 
    199 !             &                                         + ( zmskv*( vb(ji,jj,mbkt(ji,jj))+vb(ji,jj-1,mbkt(ji,jj)) ) )**2  ) 
    200 !       END DO 
    201 !    END DO 
    202 !    IF( ln_isfcav ) THEN       !top friction 
    203 !       DO jj = 2, jpjm1 
    204 !          DO ji = fs_2, fs_jpim1   ! vector opt. 
    205 !             zmsku = ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) 
    206 !             zmskv = ( 2. - vmask(ji,jj-1,mikt(ji,jj)) * vmask(ji,jj,mikt(ji,jj)) )     ! (CAUTION: CdU<0) 
    207 !             ustar2_top(ji,jj) = - rCdU_top(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mikt(ji,jj))+ub(ji-1,jj,mikt(ji,jj)) ) )**2 
    208 !                &                                         + ( zmskv*( vb(ji,jj,mikt(ji,jj))+vb(ji,jj-1,mikt(ji,jj)) ) )**2  ) 
    209 !          END DO 
    210 !       END DO 
    211 !    ENDIF 
    212 !!gm 
     177          ! bottom friction (explicit before friction) 
     178          zmsku = ( 2._wp - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) 
     179          zmskv = ( 2._wp - vmask(ji,jj-1,mbkt(ji,jj)) * vmask(ji,jj,mbkt(ji,jj)) )     ! (CAUTION: CdU<0) 
     180          ustar2_bot(ji,jj) = - rCdU_bot(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mbkt(ji,jj))+ub(ji-1,jj,mbkt(ji,jj)) ) )**2  & 
     181             &                                         + ( zmskv*( vb(ji,jj,mbkt(ji,jj))+vb(ji,jj-1,mbkt(ji,jj)) ) )**2  ) 
     182         END DO 
     183      END DO 
     184      IF( ln_isfcav ) THEN       !top friction 
     185         DO jj = 2, jpjm1 
     186            DO ji = fs_2, fs_jpim1   ! vector opt. 
     187               zmsku = ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) 
     188               zmskv = ( 2. - vmask(ji,jj-1,mikt(ji,jj)) * vmask(ji,jj,mikt(ji,jj)) )     ! (CAUTION: CdU<0) 
     189               ustar2_top(ji,jj) = - rCdU_top(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mikt(ji,jj))+ub(ji-1,jj,mikt(ji,jj)) ) )**2  & 
     190                  &                                         + ( zmskv*( vb(ji,jj,mikt(ji,jj))+vb(ji,jj-1,mikt(ji,jj)) ) )**2  ) 
     191            END DO 
     192         END DO 
     193      ENDIF 
     194    
    213195      SELECT CASE ( nn_z0_met )      !==  Set surface roughness length  ==! 
    214196      CASE ( 0 )                          ! Constant roughness           
     
    261243      ! The surface boundary condition are set after 
    262244      ! The bottom boundary condition are also set after. In standard e(bottom)=0. 
    263       ! z_elem_b : diagonal z_elem_c : upper diagonal z_elem_a : lower diagonal 
     245      ! zdiag : diagonal zd_up : upper diagonal zd_lw : lower diagonal 
    264246      ! Warning : after this step, en : right hand side of the matrix 
    265247 
     
    296278               zcof = rfact_tke * tmask(ji,jj,jk) 
    297279               !                                               ! lower diagonal 
    298                z_elem_a(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) ) / ( e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk) ) 
     280               zd_lw(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) ) / ( e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk) ) 
    299281               !                                               ! upper diagonal 
    300                z_elem_c(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) ) / ( e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk) ) 
     282               zd_up(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) ) / ( e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk) ) 
    301283               !                                               ! diagonal 
    302                z_elem_b(ji,jj,jk) = 1._wp - z_elem_a(ji,jj,jk) - z_elem_c(ji,jj,jk)  + rdt * zdiss * wmask(ji,jj,jk)  
     284               zdiag(ji,jj,jk) = 1._wp - zd_lw(ji,jj,jk) - zd_up(ji,jj,jk)  + rdt * zdiss * wmask(ji,jj,jk)  
    303285               !                                               ! right hand side in en 
    304286               en(ji,jj,jk) = en(ji,jj,jk) + rdt * zesh2 * wmask(ji,jj,jk) 
     
    307289      END DO 
    308290      ! 
    309       z_elem_b(:,:,jpk) = 1._wp 
     291      zdiag(:,:,jpk) = 1._wp 
    310292      ! 
    311293      ! Set surface condition on zwall_psi (1 at the bottom) 
     
    318300      SELECT CASE ( nn_bc_surf ) 
    319301      ! 
    320       CASE ( 0 )             ! Dirichlet case 
     302      CASE ( 0 )             ! Dirichlet boundary condition (set e at k=1 & 2)  
    321303      ! First level 
    322       en(:,:,1) = rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3 
    323       en(:,:,1) = MAX(en(:,:,1), rn_emin)  
    324       z_elem_a(:,:,1) = en(:,:,1) 
    325       z_elem_c(:,:,1) = 0._wp 
    326       z_elem_b(:,:,1) = 1._wp 
     304      en   (:,:,1) = MAX(  rn_emin , rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3  ) 
     305      zd_lw(:,:,1) = en(:,:,1) 
     306      zd_up(:,:,1) = 0._wp 
     307      zdiag(:,:,1) = 1._wp 
    327308      !  
    328309      ! One level below 
    329       en(:,:,2) = rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+gdepw_n(:,:,2)) & 
    330          &               / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 
    331       en(:,:,2) = MAX(en(:,:,2), rn_emin ) 
    332       z_elem_a(:,:,2) = 0._wp  
    333       z_elem_c(:,:,2) = 0._wp 
    334       z_elem_b(:,:,2) = 1._wp 
    335       ! 
    336       ! 
    337       CASE ( 1 )             ! Neumann boundary condition on d(e)/dz 
     310      en   (:,:,2) =  MAX(  rc02r * ustar2_surf(:,:) * (  1._wp + rsbc_tke1 * ((zhsro(:,:)+gdepw_n(:,:,2))   & 
     311         &                 / zhsro(:,:) )**(1.5_wp*ra_sf)  )**(2._wp/3._wp)                      , rn_emin   ) 
     312      zd_lw(:,:,2) = 0._wp  
     313      zd_up(:,:,2) = 0._wp 
     314      zdiag(:,:,2) = 1._wp 
     315      ! 
     316      ! 
     317      CASE ( 1 )             ! Neumann boundary condition (set d(e)/dz) 
    338318      ! 
    339319      ! Dirichlet conditions at k=1 
    340       en(:,:,1)       = rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3 
    341       en(:,:,1)       = MAX(en(:,:,1), rn_emin)       
    342       z_elem_a(:,:,1) = en(:,:,1) 
    343       z_elem_c(:,:,1) = 0._wp 
    344       z_elem_b(:,:,1) = 1._wp 
     320      en   (:,:,1) = MAX(  rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3 , rn_emin  ) 
     321      zd_lw(:,:,1) = en(:,:,1) 
     322      zd_up(:,:,1) = 0._wp 
     323      zdiag(:,:,1) = 1._wp 
    345324      ! 
    346325      ! at k=2, set de/dz=Fw 
    347326      !cbr 
    348       z_elem_b(:,:,2) = z_elem_b(:,:,2) +  z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 
    349       z_elem_a(:,:,2) = 0._wp 
    350       zkar(:,:)       = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans*gdept_n(:,:,1)/zhsro(:,:)) )) 
    351       zflxs(:,:)      = rsbc_tke2 * ustar2_surf(:,:)**1.5_wp * zkar(:,:) & 
    352           &                       * ((zhsro(:,:)+gdept_n(:,:,1)) / zhsro(:,:) )**(1.5_wp*ra_sf) 
    353  
    354       en(:,:,2) = en(:,:,2) + zflxs(:,:)/e3w_n(:,:,2) 
     327      zdiag(:,:,2) = zdiag(:,:,2) +  zd_lw(:,:,2) ! Remove zd_lw from zdiag 
     328      zd_lw(:,:,2) = 0._wp 
     329      zkar (:,:)   = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans*gdept_n(:,:,1)/zhsro(:,:)) )) 
     330      zflxs(:,:)   = rsbc_tke2 * ustar2_surf(:,:)**1.5_wp * zkar(:,:) & 
     331          &                    * (  ( zhsro(:,:)+gdept_n(:,:,1) ) / zhsro(:,:) )**(1.5_wp*ra_sf) 
     332!!gm why not   :                        * ( 1._wp + gdept_n(:,:,1) / zhsro(:,:) )**(1.5_wp*ra_sf) 
     333      en(:,:,2) = en(:,:,2) + zflxs(:,:) / e3w_n(:,:,2) 
    355334      ! 
    356335      ! 
     
    377356               ! Dirichlet condition applied at:  
    378357               !     Bottom level (ibot)      &      Just above it (ibotm1)    
    379                z_elem_a(ji,jj,ibot) = 0._wp   ;   z_elem_a(ji,jj,ibotm1) = 0._wp 
    380                z_elem_c(ji,jj,ibot) = 0._wp   ;   z_elem_c(ji,jj,ibotm1) = 0._wp 
    381                z_elem_b(ji,jj,ibot) = 1._wp   ;   z_elem_b(ji,jj,ibotm1) = 1._wp 
    382                en      (ji,jj,ibot) = z_en    ;   en      (ji,jj,ibotm1) = z_en 
    383             END DO 
    384          END DO 
    385 !!gm new 
    386 !         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
    387 !            DO jj = 2, jpjm1 
    388 !               DO ji = fs_2, fs_jpim1   ! vector opt. 
    389 !                  itop   = mikt(ji,jj)       ! k   top w-point 
    390 !                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
    391 !                  !                                                ! mask at the ocean surface points 
    392 !                  z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 
    393 !                  ! 
    394 !                  ! Dirichlet condition applied at:  
    395 !                  !     top level (itop)         &      Just below it (itopp1)   
    396 !                  z_elem_a(ji,jj,itop) = 0._wp   ;   z_elem_a(ji,jj,ipotm1) = 0._wp 
    397 !                  z_elem_c(ji,jj,itop) = 0._wp   ;   z_elem_c(ji,jj,itopp1) = 0._wp 
    398 !                  z_elem_b(ji,jj,itop) = 1._wp   ;   z_elem_b(ji,jj,itopp1) = 1._wp 
    399 !                  en      (ji,jj,itop) = zen     ;   en      (ji,jj,itopp1) = z_en 
    400 !               END DO 
    401 !            END DO 
    402 !         ENDIF 
    403 !!gm 
     358               zd_lw(ji,jj,ibot) = 0._wp   ;   zd_lw(ji,jj,ibotm1) = 0._wp 
     359               zd_up(ji,jj,ibot) = 0._wp   ;   zd_up(ji,jj,ibotm1) = 0._wp 
     360               zdiag(ji,jj,ibot) = 1._wp   ;   zdiag(ji,jj,ibotm1) = 1._wp 
     361               en   (ji,jj,ibot) = z_en    ;   en   (ji,jj,ibotm1) = z_en 
     362            END DO 
     363         END DO 
     364         ! 
     365         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
     366            DO jj = 2, jpjm1 
     367               DO ji = fs_2, fs_jpim1   ! vector opt. 
     368                  itop   = mikt(ji,jj)       ! k   top w-point 
     369                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     370                  !                                                ! mask at the ocean surface points 
     371                  z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 
     372                  ! 
     373 !!gm TO BE VERIFIED !!! 
     374                  ! Dirichlet condition applied at:  
     375                  !     top level (itop)         &      Just below it (itopp1)    
     376                  zd_lw(ji,jj,itop) = 0._wp   ;   zd_lw(ji,jj,itopp1) = 0._wp 
     377                  zd_up(ji,jj,itop) = 0._wp   ;   zd_up(ji,jj,itopp1) = 0._wp 
     378                  zdiag(ji,jj,itop) = 1._wp   ;   zdiag(ji,jj,itopp1) = 1._wp 
     379                  en   (ji,jj,itop) = z_en    ;   en   (ji,jj,itopp1) = z_en 
     380               END DO 
     381            END DO 
     382         ENDIF 
    404383         ! 
    405384      CASE ( 1 )             ! Neumman boundary condition 
     
    415394               !     Bottom level (ibot)      &      Just above it (ibotm1)    
    416395               !         Dirichlet            !         Neumann 
    417                z_elem_a(ji,jj,ibot) = 0._wp   !   ! Remove z_elem_c from z_elem_b 
    418                z_elem_b(ji,jj,ibot) = 1._wp   ;   z_elem_b(ji,jj,ibotm1) = z_elem_b(ji,jj,ibotm1) + z_elem_c(ji,jj,ibotm1) 
    419                z_elem_c(ji,jj,ibot) = 0._wp   ;   z_elem_c(ji,jj,ibotm1) = 0._wp 
    420             END DO 
    421          END DO 
    422 !!gm new 
    423 !         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
    424 !            DO jj = 2, jpjm1 
    425 !               DO ji = fs_2, fs_jpim1   ! vector opt. 
    426 !                  itop   = mikt(ji,jj)       ! k   top w-point 
    427 !                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
    428 !                  !                                                ! mask at the ocean surface points 
    429 !                  z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 
    430 !                  ! 
    431 !                  ! Bottom level Dirichlet condition: 
    432 !                  !     Bottom level (ibot)      &      Just above it (ibotm1)    
    433 !                  !         Dirichlet            !         Neumann 
    434 !                  z_elem_a(ji,jj,itop) = 0._wp   !   ! Remove z_elem_c from z_elem_b 
    435 !                  z_elem_b(ji,jj,itop) = 1._wp   ;   z_elem_b(ji,jj,itopp1) = z_elem_b(ji,jj,itopp1) + z_elem_c(ji,jj,itopp1) 
    436 !                  z_elem_c(ji,jj,itop) = 0._wp   ;   z_elem_c(ji,jj,itopp1) = 0._wp 
    437 !               END DO 
    438 !            END DO 
    439 !         ENDIF 
    440 !!gm 
     396               zd_lw(ji,jj,ibot) = 0._wp   !   ! Remove zd_up from zdiag 
     397               zdiag(ji,jj,ibot) = 1._wp   ;   zdiag(ji,jj,ibotm1) = zdiag(ji,jj,ibotm1) + zd_up(ji,jj,ibotm1) 
     398               zd_up(ji,jj,ibot) = 0._wp   ;   zd_up(ji,jj,ibotm1) = 0._wp 
     399            END DO 
     400         END DO 
     401         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
     402            DO jj = 2, jpjm1 
     403               DO ji = fs_2, fs_jpim1   ! vector opt. 
     404                  itop   = mikt(ji,jj)       ! k   top w-point 
     405                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     406                  !                                                ! mask at the ocean surface points 
     407                  z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 
     408                  ! 
     409                  ! Bottom level Dirichlet condition: 
     410                  !     Bottom level (ibot)      &      Just above it (ibotm1)    
     411                  !         Dirichlet            !         Neumann 
     412                  zd_lw(ji,jj,itop) = 0._wp   !   ! Remove zd_up from zdiag 
     413                  zdiag(ji,jj,itop) = 1._wp   ;   zdiag(ji,jj,itopp1) = zdiag(ji,jj,itopp1) + zd_up(ji,jj,itopp1) 
     414                  zd_up(ji,jj,itop) = 0._wp   ;   zd_up(ji,jj,itopp1) = 0._wp 
     415               END DO 
     416            END DO 
     417         ENDIF 
    441418         ! 
    442419      END SELECT 
     
    448425         DO jj = 2, jpjm1 
    449426            DO ji = fs_2, fs_jpim1    ! vector opt. 
    450                z_elem_b(ji,jj,jk) = z_elem_b(ji,jj,jk) - z_elem_a(ji,jj,jk) * z_elem_c(ji,jj,jk-1) / z_elem_b(ji,jj,jk-1) 
     427               zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 
    451428            END DO 
    452429         END DO 
     
    455432         DO jj = 2, jpjm1 
    456433            DO ji = fs_2, fs_jpim1    ! vector opt. 
    457                z_elem_a(ji,jj,jk) = en(ji,jj,jk) - z_elem_a(ji,jj,jk) / z_elem_b(ji,jj,jk-1) * z_elem_a(ji,jj,jk-1) 
     434               zd_lw(ji,jj,jk) = en(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) * zd_lw(ji,jj,jk-1) 
    458435            END DO 
    459436         END DO 
     
    462439         DO jj = 2, jpjm1 
    463440            DO ji = fs_2, fs_jpim1    ! vector opt. 
    464                en(ji,jj,jk) = ( z_elem_a(ji,jj,jk) - z_elem_c(ji,jj,jk) * en(ji,jj,jk+1) ) / z_elem_b(ji,jj,jk) 
     441               en(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * en(ji,jj,jk+1) ) / zdiag(ji,jj,jk) 
    465442            END DO 
    466443         END DO 
     
    519496      ! Resolution of a tridiagonal linear system by a "methode de chasse" 
    520497      ! computation from level 2 to jpkm1  (e(1) already computed and e(jpk)=0 ). 
    521       ! z_elem_b : diagonal z_elem_c : upper diagonal z_elem_a : lower diagonal 
     498      ! zdiag : diagonal zd_up : upper diagonal zd_lw : lower diagonal 
    522499      ! Warning : after this step, en : right hand side of the matrix 
    523500 
     
    543520               diss = rpsi2 * zratio * zwall(ji,jj,jk) * eps(ji,jj,jk) 
    544521               ! 
    545                zdir = 0.5_wp + SIGN( 0.5_wp, prod + buoy )   ! zdir =1(=0) if shear(ji,jj,jk)+buoy >0(<0) 
     522               zdir = 0.5_wp + SIGN( 0.5_wp, prod + buoy )     ! zdir =1(=0) if shear(ji,jj,jk)+buoy >0(<0) 
    546523               ! 
    547524               zesh2 = zdir * ( prod + buoy )          + (1._wp - zdir ) * prod                        ! production term 
     
    551528               zcof = rfact_psi * zwall_psi(ji,jj,jk) * tmask(ji,jj,jk) 
    552529               !                                               ! lower diagonal 
    553                z_elem_a(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) ) / ( e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk) ) 
     530               zd_lw(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) ) / ( e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk) ) 
    554531               !                                               ! upper diagonal 
    555                z_elem_c(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) ) / ( e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk) ) 
     532               zd_up(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) ) / ( e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk) ) 
    556533               !                                               ! diagonal 
    557                z_elem_b(ji,jj,jk) = 1._wp - z_elem_a(ji,jj,jk) - z_elem_c(ji,jj,jk) + rdt * zdiss * wmask(ji,jj,jk) 
     534               zdiag(ji,jj,jk) = 1._wp - zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) + rdt * zdiss * wmask(ji,jj,jk) 
    558535               !                                               ! right hand side in psi 
    559536               psi(ji,jj,jk) = psi(ji,jj,jk) + rdt * zesh2 * wmask(ji,jj,jk) 
     
    562539      END DO 
    563540      ! 
    564       z_elem_b(:,:,jpk) = 1._wp 
     541      zdiag(:,:,jpk) = 1._wp 
    565542 
    566543      ! Surface boundary condition on psi 
     
    574551         zdep    (:,:)   = zhsro(:,:) * rl_sf ! Cosmetic 
    575552         psi     (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 
    576          z_elem_a(:,:,1) = psi(:,:,1) 
    577          z_elem_c(:,:,1) = 0._wp 
    578          z_elem_b(:,:,1) = 1._wp 
     553         zd_lw(:,:,1) = psi(:,:,1) 
     554         zd_up(:,:,1) = 0._wp 
     555         zdiag(:,:,1) = 1._wp 
    579556         ! 
    580557         ! One level below 
     
    582559         zdep    (:,:)   = (zhsro(:,:) + gdepw_n(:,:,2)) * zkar(:,:) 
    583560         psi     (:,:,2) = rc0**rpp * en(:,:,2)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 
    584          z_elem_a(:,:,2) = 0._wp 
    585          z_elem_c(:,:,2) = 0._wp 
    586          z_elem_b(:,:,2) = 1._wp 
     561         zd_lw(:,:,2) = 0._wp 
     562         zd_up(:,:,2) = 0._wp 
     563         zdiag(:,:,2) = 1._wp 
    587564         !  
    588565      CASE ( 1 )             ! Neumann boundary condition on d(psi)/dz 
     
    591568         zdep    (:,:)   = zhsro(:,:) * rl_sf 
    592569         psi     (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 
    593          z_elem_a(:,:,1) = psi(:,:,1) 
    594          z_elem_c(:,:,1) = 0._wp 
    595          z_elem_b(:,:,1) = 1._wp 
     570         zd_lw(:,:,1) = psi(:,:,1) 
     571         zd_up(:,:,1) = 0._wp 
     572         zdiag(:,:,1) = 1._wp 
    596573         ! 
    597574         ! Neumann condition at k=2 
    598          z_elem_b(:,:,2) = z_elem_b(:,:,2) +  z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 
    599          z_elem_a(:,:,2) = 0._wp 
     575         zdiag(:,:,2) = zdiag(:,:,2) +  zd_lw(:,:,2) ! Remove zd_lw from zdiag 
     576         zd_lw(:,:,2) = 0._wp 
    600577         ! 
    601578         ! Set psi vertical flux at the surface: 
     
    613590      ! -------------------------------- 
    614591      ! 
    615       SELECT CASE ( nn_bc_bot ) 
     592!!gm should be done for ISF (top boundary cond.) 
     593!!gm so, totally new staff needed      ===>>> think about that ! 
     594! 
     595      SELECT CASE ( nn_bc_bot )     ! bottom boundary 
    616596      ! 
    617597      CASE ( 0 )             ! Dirichlet  
    618          !                      ! en(ibot) = u*^2 / Co2 and hmxl_n(ibot) = vkarmn * rn_bfrz0 
     598         !                      ! en(ibot) = u*^2 / Co2 and hmxl_n(ibot) = vkarmn * r_z0_bot 
    619599         !                      ! Balance between the production and the dissipation terms 
    620600         DO jj = 2, jpjm1 
     
    622602               ibot   = mbkt(ji,jj) + 1      ! k   bottom level of w-point 
    623603               ibotm1 = mbkt(ji,jj)          ! k-1 bottom level of w-point but >=1 
    624                zdep(ji,jj) = vkarmn * rn_bfrz0 
     604               zdep(ji,jj) = vkarmn * r_z0_bot 
    625605               psi (ji,jj,ibot) = rc0**rpp * en(ji,jj,ibot)**rmm * zdep(ji,jj)**rnn 
    626                z_elem_a(ji,jj,ibot) = 0._wp 
    627                z_elem_c(ji,jj,ibot) = 0._wp 
    628                z_elem_b(ji,jj,ibot) = 1._wp 
     606               zd_lw(ji,jj,ibot) = 0._wp 
     607               zd_up(ji,jj,ibot) = 0._wp 
     608               zdiag(ji,jj,ibot) = 1._wp 
    629609               ! 
    630610               ! Just above last level, Dirichlet condition again (GOTM like) 
    631                zdep(ji,jj) = vkarmn * ( rn_bfrz0 + e3t_n(ji,jj,ibotm1) ) 
     611               zdep(ji,jj) = vkarmn * ( r_z0_bot + e3t_n(ji,jj,ibotm1) ) 
    632612               psi (ji,jj,ibotm1) = rc0**rpp * en(ji,jj,ibot  )**rmm * zdep(ji,jj)**rnn 
    633                z_elem_a(ji,jj,ibotm1) = 0._wp 
    634                z_elem_c(ji,jj,ibotm1) = 0._wp 
    635                z_elem_b(ji,jj,ibotm1) = 1._wp 
     613               zd_lw(ji,jj,ibotm1) = 0._wp 
     614               zd_up(ji,jj,ibotm1) = 0._wp 
     615               zdiag(ji,jj,ibotm1) = 1._wp 
    636616            END DO 
    637617         END DO 
     
    645625               ! 
    646626               ! Bottom level Dirichlet condition: 
    647                zdep(ji,jj) = vkarmn * rn_bfrz0 
     627               zdep(ji,jj) = vkarmn * r_z0_bot 
    648628               psi (ji,jj,ibot) = rc0**rpp * en(ji,jj,ibot)**rmm * zdep(ji,jj)**rnn 
    649629               ! 
    650                z_elem_a(ji,jj,ibot) = 0._wp 
    651                z_elem_c(ji,jj,ibot) = 0._wp 
    652                z_elem_b(ji,jj,ibot) = 1._wp 
     630               zd_lw(ji,jj,ibot) = 0._wp 
     631               zd_up(ji,jj,ibot) = 0._wp 
     632               zdiag(ji,jj,ibot) = 1._wp 
    653633               ! 
    654634               ! Just above last level: Neumann condition with flux injection 
    655                z_elem_b(ji,jj,ibotm1) = z_elem_b(ji,jj,ibotm1) + z_elem_c(ji,jj,ibotm1) ! Remove z_elem_c from z_elem_b 
    656                z_elem_c(ji,jj,ibotm1) = 0. 
     635               zdiag(ji,jj,ibotm1) = zdiag(ji,jj,ibotm1) + zd_up(ji,jj,ibotm1) ! Remove zd_up from zdiag 
     636               zd_up(ji,jj,ibotm1) = 0. 
    657637               ! 
    658638               ! Set psi vertical flux at the bottom: 
    659                zdep(ji,jj) = rn_bfrz0 + 0.5_wp*e3t_n(ji,jj,ibotm1) 
     639               zdep(ji,jj) = r_z0_bot + 0.5_wp*e3t_n(ji,jj,ibotm1) 
    660640               zflxb = rsbc_psi2 * ( p_avm(ji,jj,ibot) + p_avm(ji,jj,ibotm1) )   & 
    661641                  &  * (0.5_wp*(en(ji,jj,ibot)+en(ji,jj,ibotm1)))**rmm * zdep(ji,jj)**(rnn-1._wp) 
     
    672652         DO jj = 2, jpjm1 
    673653            DO ji = fs_2, fs_jpim1    ! vector opt. 
    674                z_elem_b(ji,jj,jk) = z_elem_b(ji,jj,jk) - z_elem_a(ji,jj,jk) * z_elem_c(ji,jj,jk-1) / z_elem_b(ji,jj,jk-1) 
     654               zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 
    675655            END DO 
    676656         END DO 
     
    679659         DO jj = 2, jpjm1 
    680660            DO ji = fs_2, fs_jpim1    ! vector opt. 
    681                z_elem_a(ji,jj,jk) = psi(ji,jj,jk) - z_elem_a(ji,jj,jk) / z_elem_b(ji,jj,jk-1) * z_elem_a(ji,jj,jk-1) 
     661               zd_lw(ji,jj,jk) = psi(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) * zd_lw(ji,jj,jk-1) 
    682662            END DO 
    683663         END DO 
     
    686666         DO jj = 2, jpjm1 
    687667            DO ji = fs_2, fs_jpim1    ! vector opt. 
    688                psi(ji,jj,jk) = ( z_elem_a(ji,jj,jk) - z_elem_c(ji,jj,jk) * psi(ji,jj,jk+1) ) / z_elem_b(ji,jj,jk) 
     668               psi(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * psi(ji,jj,jk+1) ) / zdiag(ji,jj,jk) 
    689669            END DO 
    690670         END DO 
     
    814794      zstm(:,:,1) = zstm(:,:,2) 
    815795 
    816 !!gm should be done for ISF (top boundary cond.) 
    817796      DO jj = 2, jpjm1 
    818797         DO ji = fs_2, fs_jpim1   ! vector opt. 
     
    820799         END DO 
    821800      END DO 
     801!!gm should be done for ISF (top boundary cond.) 
     802!!gm so, totally new staff needed!!gm 
    822803 
    823804      ! Compute diffusivities/viscosities 
     
    904885         WRITE(numout,*) '      Type of closure                               nn_clos        = ', nn_clos 
    905886         WRITE(numout,*) '      Surface roughness (m)                         rn_hsro        = ', rn_hsro 
    906 !!gm old 
    907          WRITE(numout,*) '      top    roughness (m) (nambfr namelist)        rn_tfrz0       = ', rn_tfrz0 
    908          WRITE(numout,*) '      Bottom roughness (m) (nambfr namelist)        rn_bfrz0       = ', rn_bfrz0 
    909 !!gm new 
    910 !         WRITE(numout,*) '   Namelist namdrg_top/_bot used values:' 
    911 !         WRITE(numout,*) '      top    ocean cavity roughness (m)             rn_z0(_top)   = ', r_z0_top 
    912 !         WRITE(numout,*) '      Bottom seafloor     roughness (m)             rn_z0(_bot)   = ', r_z0_bot 
    913 !!gm 
     887         WRITE(numout,*) 
     888         WRITE(numout,*) '   Namelist namdrg_top/_bot:   used values:' 
     889         WRITE(numout,*) '      top    ocean cavity roughness (m)             rn_z0(_top)   = ', r_z0_top 
     890         WRITE(numout,*) '      Bottom seafloor     roughness (m)             rn_z0(_bot)   = ', r_z0_bot 
    914891         WRITE(numout,*) 
    915892      ENDIF 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90

    r7753 r8143  
    1111   !!   zdf_mxl      : Compute the turbocline and mixed layer depths. 
    1212   !!---------------------------------------------------------------------- 
    13    USE oce             ! ocean dynamics and tracers variables 
    14    USE dom_oce         ! ocean space and time domain variables 
    15    USE trc_oce, ONLY: l_offline         ! ocean space and time domain variables 
    16    USE zdf_oce         ! ocean vertical physics 
    17    USE in_out_manager  ! I/O manager 
    18    USE prtctl          ! Print control 
    19    USE phycst          ! physical constants 
    20    USE iom             ! I/O library 
    21    USE lib_mpp         ! MPP library 
    22    USE wrk_nemo        ! work arrays 
    23    USE timing          ! Timing 
     13   USE oce            ! ocean dynamics and tracers variables 
     14   USE dom_oce        ! ocean space and time domain variables 
     15   USE trc_oce  , ONLY: l_offline         ! ocean space and time domain variables 
     16   USE zdf_oce        ! ocean vertical physics 
     17   USE in_out_manager ! I/O manager 
     18   USE prtctl         ! Print control 
     19   USE phycst         ! physical constants 
     20   USE iom            ! I/O library 
     21   USE lib_mpp        ! MPP library 
     22   USE timing         ! Timing 
    2423 
    2524   IMPLICIT NONE 
    2625   PRIVATE 
    2726 
    28    PUBLIC   zdf_mxl       ! called by step.F90 
     27   PUBLIC   zdf_mxl   ! called by zdfphy.F90 
    2928 
    30    INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   nmln    !: number of level in the mixed layer (used by TOP) 
    31    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmld    !: mixing layer depth (turbocline)      [m] 
    32    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlp    !: mixed layer depth  (rho=rho0+zdcrit) [m] 
    33    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlpt   !: depth of the last T-point inside the mixed layer [m] 
     29   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   nmln    !: number of level in the mixed layer (used by LDF, ZDF, TRD, TOP) 
     30   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmld    !: mixing layer depth (turbocline)      [m]   (used by TOP) 
     31   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlp    !: mixed layer depth  (rho=rho0+zdcrit) [m]   (used by LDF) 
     32   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlpt   !: depth of the last T-point inside the mixed layer [m] (used by LDF) 
    3433 
    3534   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    !: density criterion for mixed layer depth 
     
    3736 
    3837   !!---------------------------------------------------------------------- 
    39    !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
     38   !! NEMO/OPA 4.0 , NEMO Consortium (2017) 
    4039   !! $Id$  
    4140   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    8079      INTEGER  ::   iikn, iiki, ikt ! local integer 
    8180      REAL(wp) ::   zN2_c           ! local scalar 
    82       INTEGER, POINTER, DIMENSION(:,:) ::   imld   ! 2D workspace 
     81      INTEGER, DIMENSION(jpi,jpj) ::   imld   ! 2D workspace 
    8382      !!---------------------------------------------------------------------- 
    8483      ! 
    8584      IF( nn_timing == 1 )  CALL timing_start('zdf_mxl') 
    8685      ! 
    87       CALL wrk_alloc( jpi,jpj, imld ) 
    88  
    8986      IF( kt == nit000 ) THEN 
    9087         IF(lwp) WRITE(numout,*) 
     
    144141      IF(ln_ctl)   CALL prt_ctl( tab2d_1=REAL(nmln,wp), clinfo1=' nmln : ', tab2d_2=hmlp, clinfo2=' hmlp : ', ovlap=1 ) 
    145142      ! 
    146       CALL wrk_dealloc( jpi,jpj, imld ) 
    147       ! 
    148143      IF( nn_timing == 1 )  CALL timing_stop('zdf_mxl') 
    149144      ! 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfphy.F90

    r8093 r8143  
    1414   USE zdf_oce        ! vertical physics: shared variables          
    1515   USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    16 !!gm old 
    17    USE zdfbfr         ! vertical physics: bottom friction 
    18 !!gm 
    1916   USE zdfsh2         ! vertical physics: shear production term of TKE 
    2017   USE zdfric         ! vertical physics: RIChardson dependent vertical mixing    
     
    195192      !                          !== top/bottom friction  ==! 
    196193      CALL zdf_drg_init 
    197 !!gm old 
    198       CALL zdf_bfr_init 
    199 !!gm 
    200194      ! 
    201195      !                          !== time-stepping  ==! 
     
    224218      !! --------------------------------------------------------------------- 
    225219      ! 
    226 !!gm old 
    227       CALL zdf_bfr( kt )                        !* bottom friction (if quadratic) 
    228 !!gm 
    229       ! 
    230 !      IF( l_zdfdrg ) THEN     !==  update top/bottom drag  ==!   (non-linear cases) 
    231 !         ! 
    232 !         !                       !* bottom drag 
    233 !         CALL zdf_drg( kt, mbkt    , r_Cdmin_bot, r_Cdmax_bot,   &   ! <<== in  
    234 !            &              r_z0_bot,   r_ke0_bot,    rCd0_bot,   & 
    235 !            &                                        rCdU_bot  )     ! ==>> out : bottom drag [m/s] 
    236 !         IF( ln_isfcav ) THEN    !* top drag   (ocean cavities) 
    237 !            CALL zdf_drg( kt, mikt    , r_Cdmin_top, r_Cdmax_top,   &   ! <<== in  
    238 !               &              r_z0_top,   r_ke0_top,    rCd0_top,   & 
    239 !               &                                        rCdU_top  )     ! ==>> out : bottom drag [m/s] 
    240 !         ENDIF 
    241 !      ENDIF 
     220      IF( l_zdfdrg ) THEN     !==  update top/bottom drag  ==!   (non-linear cases) 
     221         ! 
     222         !                       !* bottom drag 
     223         CALL zdf_drg( kt, mbkt    , r_Cdmin_bot, r_Cdmax_bot,   &   ! <<== in  
     224            &              r_z0_bot,   r_ke0_bot,    rCd0_bot,   & 
     225            &                                        rCdU_bot  )     ! ==>> out : bottom drag [m/s] 
     226         IF( ln_isfcav ) THEN    !* top drag   (ocean cavities) 
     227            CALL zdf_drg( kt, mikt    , r_Cdmin_top, r_Cdmax_top,   &   ! <<== in  
     228               &              r_z0_top,   r_ke0_top,    rCd0_top,   & 
     229               &                                        rCdU_top  )     ! ==>> out : bottom drag [m/s] 
     230         ENDIF 
     231      ENDIF 
    242232      ! 
    243233      !                       !==  Kz from chosen turbulent closure  ==!   (avm_k, avt_k) 
     
    290280      CALL lbc_lnk( avs  , 'W', 1. )                  !!gm  To be tested 
    291281      ! 
    292  
    293  
     282      IF( l_zdfdrg ) THEN     ! drag  have been updated (non-linear cases) 
     283         IF( ln_isfcav )   CALL lbc_lnk( rCdU_top, 'T', 1. )   ! top    drag 
     284                           CALL lbc_lnk( rCdU_bot, 'T', 1. )   ! bottom drag 
     285      ENDIF 
     286      ! 
    294287      CALL zdf_mxl( kt )                        !* mixed layer depth, and level 
    295288 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfric.F90

    r8093 r8143  
    6161   SUBROUTINE zdf_ric_init 
    6262      !!---------------------------------------------------------------------- 
    63       !!                 ***  ROUTINE zdfbfr_init  *** 
     63      !!                 ***  ROUTINE zdf_ric_init  *** 
    6464      !!                     
    6565      !! ** Purpose :   Initialization of the vertical eddy diffusivity and 
     
    9090      IF(lwp) THEN                   ! Control print 
    9191         WRITE(numout,*) 
    92          WRITE(numout,*) 'zdf_ric : Ri depend vertical mixing scheme' 
    93          WRITE(numout,*) '~~~~~~~' 
     92         WRITE(numout,*) 'zdf_ric_init : Ri depend vertical mixing scheme' 
     93         WRITE(numout,*) '~~~~~~~~~~~~' 
    9494         WRITE(numout,*) '   Namelist namzdf_ric : set Kz=F(Ri) parameters' 
    9595         WRITE(numout,*) '      maximum vertical viscosity        rn_avmri  = ', rn_avmri 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r8093 r8143  
    2828   !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability 
    2929   !!            4.0  !  2017-04  (G. Madec)  remove CPP ddm key & avm at t-point only  
     30   !!             -   !  2017-05  (G. Madec)  add top/bottom friction as boundary condition (ln_drg) 
    3031   !!---------------------------------------------------------------------- 
    3132 
     
    4344   USE sbc_oce        ! surface boundary condition: ocean 
    4445   USE zdf_oce        ! vertical physics: ocean variables 
    45 !!gm new 
    4646   USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    47 !!gm 
    4847   USE zdfmxl         ! vertical physics: mixed layer 
    4948#if defined key_agrif 
     
    5756   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    5857   USE prtctl         ! Print control 
    59    USE wrk_nemo       ! work arrays 
    6058   USE timing         ! Timing 
    6159   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     
    7977   REAL(wp) ::   rn_emin0  ! surface minimum value of tke   [m2/s2] 
    8078   REAL(wp) ::   rn_bshear ! background shear (>0) currently a numerical threshold (do not change it) 
     79   LOGICAL  ::   ln_drg    ! top/bottom friction forcing flag  
    8180   INTEGER  ::   nn_etau   ! type of depth penetration of surface tke (=0/1/2/3) 
    82    INTEGER  ::   nn_htau   ! type of tke profile of penetration (=0/1) 
    83    REAL(wp) ::   rn_efr    ! fraction of TKE surface value which penetrates in the ocean 
     81   INTEGER  ::      nn_htau   ! type of tke profile of penetration (=0/1) 
     82   REAL(wp) ::      rn_efr    ! fraction of TKE surface value which penetrates in the ocean 
    8483   LOGICAL  ::   ln_lc     ! Langmuir cells (LC) as a source term of TKE or not 
    85    REAL(wp) ::   rn_lc     ! coef to compute vertical velocity of Langmuir cells 
     84   REAL(wp) ::      rn_lc     ! coef to compute vertical velocity of Langmuir cells 
    8685 
    8786   REAL(wp) ::   ri_cri    ! critic Richardson number (deduced from rn_ediff and rn_ediss values) 
     
    204203      ! 
    205204      INTEGER ::   ji, jj, jk              ! dummy loop arguments 
    206 !!bfr      REAL(wp) ::   zebot, zmshu, zmskv      ! local scalars 
    207       REAL(wp) ::   zrhoa  = 1.22            ! Air density kg/m3 
    208       REAL(wp) ::   zcdrag = 1.5e-3          ! drag coefficient 
    209       REAL(wp) ::   zbbrau, zri              ! local scalars 
    210       REAL(wp) ::   zfact1, zfact2, zfact3   !   -         - 
    211       REAL(wp) ::   ztx2  , zty2  , zcof     !   -         - 
    212       REAL(wp) ::   ztau  , zdif             !   -         - 
    213       REAL(wp) ::   zus   , zwlc  , zind     !   -         - 
    214       REAL(wp) ::   zzd_up, zzd_lw           !   -         - 
     205      REAL(wp) ::   zetop, zebot, zmsku, zmskv ! local scalars 
     206      REAL(wp) ::   zrhoa  = 1.22              ! Air density kg/m3 
     207      REAL(wp) ::   zcdrag = 1.5e-3            ! drag coefficient 
     208      REAL(wp) ::   zbbrau, zri                ! local scalars 
     209      REAL(wp) ::   zfact1, zfact2, zfact3     !   -         - 
     210      REAL(wp) ::   ztx2  , zty2  , zcof       !   -         - 
     211      REAL(wp) ::   ztau  , zdif               !   -         - 
     212      REAL(wp) ::   zus   , zwlc  , zind       !   -         - 
     213      REAL(wp) ::   zzd_up, zzd_lw             !   -         - 
    215214      INTEGER , DIMENSION(jpi,jpj)     ::   imlc 
    216215      REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc 
     
    227226      ! 
    228227      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    229       !                     !  Surface boundary condition on tke 
    230       !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     228      !                     !  Surface/top/bottom boundary condition on tke 
     229      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     230       
     231      DO jj = 2, jpjm1            ! en(1)   = rn_ebb taum / rau0  (min value rn_emin0) 
     232         DO ji = fs_2, fs_jpim1   ! vector opt. 
     233            en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1) 
     234         END DO 
     235      END DO 
    231236      IF ( ln_isfcav ) THEN 
    232237         DO jj = 2, jpjm1            ! en(mikt(ji,jj))   = rn_emin 
     
    235240            END DO 
    236241         END DO 
    237       END IF 
    238       DO jj = 2, jpjm1            ! en(1)   = rn_ebb taum / rau0  (min value rn_emin0) 
    239          DO ji = fs_2, fs_jpim1   ! vector opt. 
    240             en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1) 
    241          END DO 
    242       END DO 
     242      ENDIF 
    243243       
    244 !!bfr   - start commented area 
    245244      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    246245      !                     !  Bottom boundary condition on tke 
    247246      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    248247      ! 
    249       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    250       ! Tests to date have found the bottom boundary condition on tke to have very little effect. 
    251       ! The condition is coded here for completion but commented out until there is proof that the 
    252       ! computational cost is justified 
    253       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    254       !                     en(bot)   = (rn_ebb0/rau0)*0.5*sqrt(u_botfr^2+v_botfr^2) (min value rn_emin) 
    255 !!gm old 
    256 !!    DO jj = 2, jpjm1 
    257 !!       DO ji = fs_2, fs_jpim1   ! vector opt. 
    258 !!          ztx2 = bfrua(ji-1,jj) * ub(ji-1,jj,mbku(ji-1,jj)) + & 
    259 !!                 bfrua(ji  ,jj) * ub(ji  ,jj,mbku(ji  ,jj) ) 
    260 !!          zty2 = bfrva(ji,jj  ) * vb(ji,jj  ,mbkv(ji,jj  )) + & 
    261 !!                 bfrva(ji,jj-1) * vb(ji,jj-1,mbkv(ji,jj-1) ) 
    262 !!          zebot = 0.001875_wp * SQRT( ztx2 * ztx2 + zty2 * zty2 )   !  where 0.001875 = (rn_ebb0/rau0) * 0.5 = 3.75*0.5/1000. 
    263 !!          en(ji,jj,mbkt(ji,jj)+1) = MAX( zebot, rn_emin ) * ssmask(ji,jj) 
    264 !!       END DO 
    265 !!    END DO 
    266 !!gm new 
    267 !! 
    268 !!    ====>>>>  add below an wet-only calculation of u and v at t-point like in zdfsh2: 
    269 !!                   zmsku = ( 2. - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) 
    270 !!                   zmskv = ( 2. - vmask(ji,jj-1,mbkt(ji,jj)) * vmask(ji,jj,mbkt(ji,jj)) ) 
    271 !! 
    272 !! 
    273 !!    DO jj = 2, jpjm1 
    274 !!       DO ji = fs_2, fs_jpim1   ! vector opt. 
    275 !!          zmsku = ( 2. - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) 
    276 !!          zmskv = ( 2. - vmask(ji,jj-1,mbkt(ji,jj)) * vmask(ji,jj,mbkt(ji,jj)) ) 
    277 !!          !                             ! where 0.001875 = (rn_ebb0/rau0) * 0.5 = 3.75*0.5/1000. (CAUTION CdU<0) 
    278 !!          zebot = - 0.001875_wp * rCdU_bot(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mbkt(ji,jj))+ub(ji-1,jj,mbkt(ji,jj)) ) )**2 
    279 !!             &                                           + ( zmskv*( vb(ji,jj,mbkt(ji,jj))+vb(ji,jj-1,mbkt(ji,jj)) ) )**2  ) 
    280 !!          en(ji,jj,mbkt(ji,jj)+1) = MAX( zebot, rn_emin ) * ssmask(ji,jj) 
    281 !!       END DO 
    282 !!    END DO 
    283 !!    IF( ln_isfcav ) THEN       !top friction 
    284 !!       DO jj = 2, jpjm1 
    285 !!          DO ji = fs_2, fs_jpim1   ! vector opt. 
    286 !!             zmsku = ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) 
    287 !!             zmskv = ( 2. - vmask(ji,jj-1,mikt(ji,jj)) * vmask(ji,jj,mikt(ji,jj)) ) 
    288 !!             !                             ! where 0.001875 = (rn_ebb0/rau0) * 0.5 = 3.75*0.5/1000.  (CAUTION CdU<0) 
    289 !!             zebot = - 0.001875_wp * rCdU_top(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mikt(ji,jj))+ub(ji-1,jj,mikt(ji,jj)) ) )**2 
    290 !!                &                                           + ( zmskv*( vb(ji,jj,mikt(ji,jj))+vb(ji,jj-1,mikt(ji,jj)) ) )**2  ) 
    291 !!             en(ji,jj,mikt(ji,jj)+1) = MAX( zebot, rn_emin ) * (1._wp - tmask(ji,jj,1))   ! masked at ocean surface 
    292 !!          END DO 
    293 !!       END DO 
    294 !!    ENDIF 
    295 !! 
    296 !!bfr   - end commented area 
     248      !   en(bot)   = (ebb0/rau0)*0.5*sqrt(u_botfr^2+v_botfr^2) (min value rn_emin) 
     249      ! where ebb0 does not includes surface wave enhancement (i.e. ebb0=3.75) 
     250      ! Note that stress averaged is done using an wet-only calculation of u and v at t-point like in zdfsh2 
     251      ! 
     252      IF( ln_drg ) THEN       !== friction used as top/bottom boundary condition on TKE 
     253         ! 
     254         DO jj = 2, jpjm1           ! bottom friction 
     255            DO ji = fs_2, fs_jpim1     ! vector opt. 
     256               zmsku = ( 2. - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) 
     257               zmskv = ( 2. - vmask(ji,jj-1,mbkt(ji,jj)) * vmask(ji,jj,mbkt(ji,jj)) ) 
     258               !                       ! where 0.001875 = (rn_ebb0/rau0) * 0.5 = 3.75*0.5/1000. (CAUTION CdU<0) 
     259               zebot = - 0.001875_wp * rCdU_bot(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mbkt(ji,jj))+ub(ji-1,jj,mbkt(ji,jj)) ) )**2  & 
     260                  &                                           + ( zmskv*( vb(ji,jj,mbkt(ji,jj))+vb(ji,jj-1,mbkt(ji,jj)) ) )**2  ) 
     261               en(ji,jj,mbkt(ji,jj)+1) = MAX( zebot, rn_emin ) * ssmask(ji,jj) 
     262            END DO 
     263         END DO 
     264         IF( ln_isfcav ) THEN       ! top friction 
     265            DO jj = 2, jpjm1 
     266               DO ji = fs_2, fs_jpim1   ! vector opt. 
     267                  zmsku = ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) 
     268                  zmskv = ( 2. - vmask(ji,jj-1,mikt(ji,jj)) * vmask(ji,jj,mikt(ji,jj)) ) 
     269                  !                             ! where 0.001875 = (rn_ebb0/rau0) * 0.5 = 3.75*0.5/1000.  (CAUTION CdU<0) 
     270                  zetop = - 0.001875_wp * rCdU_top(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mikt(ji,jj))+ub(ji-1,jj,mikt(ji,jj)) ) )**2  & 
     271                     &                                           + ( zmskv*( vb(ji,jj,mikt(ji,jj))+vb(ji,jj-1,mikt(ji,jj)) ) )**2  ) 
     272                  en(ji,jj,mikt(ji,jj)) = MAX( zetop, rn_emin ) * (1._wp - tmask(ji,jj,1))   ! masked at ocean surface 
     273               END DO 
     274            END DO 
     275         ENDIF 
     276         ! 
     277      ENDIF 
    297278      ! 
    298279      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     
    426407      !                            !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    427408!!gm BUG : in the exp  remove the depth of ssh !!! 
     409!!gm       i.e. use gde3w in argument (pdepw) 
    428410       
    429411       
     
    678660      NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin  ,   & 
    679661         &                 rn_emin0, rn_bshear, nn_mxl , ln_mxl0  ,   & 
    680          &                 rn_mxl0 , nn_pdl   , ln_lc  , rn_lc    ,   & 
     662         &                 rn_mxl0 , nn_pdl   , ln_drg , ln_lc    , rn_lc    ,   & 
    681663         &                 nn_etau , nn_htau  , rn_efr    
    682664      !!---------------------------------------------------------------------- 
     
    703685         WRITE(numout,*) '      minimum value of tke                        rn_emin   = ', rn_emin 
    704686         WRITE(numout,*) '      surface minimum value of tke                rn_emin0  = ', rn_emin0 
     687         WRITE(numout,*) '      prandl number flag                          nn_pdl    = ', nn_pdl 
    705688         WRITE(numout,*) '      background shear (>0)                       rn_bshear = ', rn_bshear 
    706689         WRITE(numout,*) '      mixing length type                          nn_mxl    = ', nn_mxl 
    707          WRITE(numout,*) '      prandl number flag                          nn_pdl    = ', nn_pdl 
    708          WRITE(numout,*) '      surface mixing length = F(stress) or not    ln_mxl0   = ', ln_mxl0 
    709          WRITE(numout,*) '      surface  mixing length minimum value        rn_mxl0   = ', rn_mxl0 
    710          WRITE(numout,*) '      flag to take into acc.  Langmuir circ.      ln_lc     = ', ln_lc 
    711          WRITE(numout,*) '      coef to compute verticla velocity of LC     rn_lc     = ', rn_lc 
     690         WRITE(numout,*) '         surface mixing length = F(stress) or not    ln_mxl0   = ', ln_mxl0 
     691         WRITE(numout,*) '         surface  mixing length minimum value        rn_mxl0   = ', rn_mxl0 
     692         WRITE(numout,*) '      top/bottom friction forcing flag            ln_drg    = ', ln_drg 
     693         WRITE(numout,*) '      Langmuir cells parametrization              ln_lc     = ', ln_lc 
     694         WRITE(numout,*) '         coef to compute vertical velocity of LC     rn_lc  = ', rn_lc 
    712695         WRITE(numout,*) '      test param. to add tke induced by wind      nn_etau   = ', nn_etau 
    713          WRITE(numout,*) '      flag for computation of exp. tke profile    nn_htau   = ', nn_htau 
    714          WRITE(numout,*) '      fraction of en which pene. the thermocline  rn_efr    = ', rn_efr 
     696         WRITE(numout,*) '          type of tke penetration profile            nn_htau   = ', nn_htau 
     697         WRITE(numout,*) '          fraction of TKE that penetrates            rn_efr    = ', rn_efr 
    715698         WRITE(numout,*) 
    716          WRITE(numout,*) '      critical Richardson nb with your parameters  ri_cri = ', ri_cri 
     699         IF( ln_drg ) THEN 
     700            WRITE(numout,*) '   Namelist namdrg_top/_bot:   used values:' 
     701            WRITE(numout,*) '      top    ocean cavity roughness (m)          rn_z0(_top)= ', r_z0_top 
     702            WRITE(numout,*) '      Bottom seafloor     roughness (m)          rn_z0(_bot)= ', r_z0_bot 
     703         ENDIF 
     704         WRITE(numout,*) 
     705         WRITE(numout,*) 
     706         WRITE(numout,*) '   ==>> critical Richardson nb with your parameters  ri_cri = ', ri_cri 
    717707         WRITE(numout,*) 
    718708      ENDIF 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/step.F90

    r8093 r8143  
    129129 
    130130      !  VERTICAL PHYSICS 
    131                          CALL zdf_phy( kstp )         ! vertical physics update (bfr, avt, avs, avm + MLD) 
     131                         CALL zdf_phy( kstp )         ! vertical physics update (top/bot drag, avt, avs, avm + MLD) 
    132132 
    133133 
     
    211211      ENDIF 
    212212 
    213                          CALL dyn_bfr       ( kstp )  ! bottom friction 
     213      IF( .NOT.ln_drgimp)   CALL dyn_bfr       ( kstp )  ! bottom friction 
     214       
    214215                         CALL dyn_zdf       ( kstp )  ! vertical diffusion 
    215216 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/step_oce.F90

    r7990 r8143  
    77   !!             3.7  !  2014-01  (G. Madec) LDF simplication  
    88   !!---------------------------------------------------------------------- 
    9    USE oce              ! ocean dynamics and tracers variables 
    10    USE dom_oce          ! ocean space and time domain variables 
    11    USE zdf_oce          ! ocean vertical physics variables 
     9   USE oce             ! ocean dynamics and tracers variables 
     10   USE dom_oce         ! ocean space and time domain variables 
     11   USE zdf_oce         ! ocean vertical physics variables 
     12   USE zdfdrg  ,  ONLY : ln_drgimp   ! implicit top/bottom friction 
    1213 
    13    USE daymod           ! calendar                         (day     routine) 
     14   USE daymod          ! calendar                         (day     routine) 
    1415 
    15    USE sbc_oce          ! surface boundary condition: ocean 
    16    USE sbcmod           ! surface boundary condition       (sbc     routine) 
    17    USE sbcrnf           ! surface boundary condition: runoff variables 
    18    USE sbccpl           ! surface boundary condition: coupled formulation (call send at end of step) 
    19    USE sbcapr           ! surface boundary condition: atmospheric pressure 
    20    USE sbctide          ! Tide initialisation 
    21    USE sbcwave          ! Wave intialisation 
     16   USE sbc_oce         ! surface boundary condition: ocean 
     17   USE sbcmod          ! surface boundary condition       (sbc     routine) 
     18   USE sbcrnf          ! surface boundary condition: runoff variables 
     19   USE sbccpl          ! surface boundary condition: coupled formulation (call send at end of step) 
     20   USE sbcapr          ! surface boundary condition: atmospheric pressure 
     21   USE sbctide         ! Tide initialisation 
     22   USE sbcwave         ! Wave intialisation 
    2223 
    23    USE traqsr           ! solar radiation penetration      (tra_qsr routine) 
    24    USE trasbc           ! surface boundary condition       (tra_sbc routine) 
    25    USE trabbc           ! bottom boundary condition        (tra_bbc routine) 
    26    USE trabbl           ! bottom boundary layer            (tra_bbl routine) 
    27    USE tradmp           ! internal damping                 (tra_dmp routine) 
    28    USE traadv           ! advection scheme control     (tra_adv_ctl routine) 
    29    USE traldf           ! lateral mixing                   (tra_ldf routine) 
    30    USE trazdf           ! vertical mixing                  (tra_zdf routine) 
    31    USE tranxt           ! time-stepping                    (tra_nxt routine) 
    32    USE tranpc           ! non-penetrative convection       (tra_npc routine) 
     24   USE traqsr          ! solar radiation penetration      (tra_qsr routine) 
     25   USE trasbc          ! surface boundary condition       (tra_sbc routine) 
     26   USE trabbc          ! bottom boundary condition        (tra_bbc routine) 
     27   USE trabbl          ! bottom boundary layer            (tra_bbl routine) 
     28   USE tradmp          ! internal damping                 (tra_dmp routine) 
     29   USE traadv          ! advection scheme control     (tra_adv_ctl routine) 
     30   USE traldf          ! lateral mixing                   (tra_ldf routine) 
     31   USE trazdf          ! vertical mixing                  (tra_zdf routine) 
     32   USE tranxt          ! time-stepping                    (tra_nxt routine) 
     33   USE tranpc          ! non-penetrative convection       (tra_npc routine) 
    3334 
    34    USE eosbn2           ! equation of state                (eos_bn2 routine) 
     35   USE eosbn2          ! equation of state                (eos_bn2 routine) 
    3536 
    36    USE divhor           ! horizontal divergence            (div_hor routine) 
    37    USE dynadv           ! advection                        (dyn_adv routine) 
    38    USE dynbfr           ! Bottom friction terms            (dyn_bfr routine) 
    39    USE dynvor           ! vorticity term                   (dyn_vor routine) 
    40    USE dynhpg           ! hydrostatic pressure grad.       (dyn_hpg routine) 
    41    USE dynldf           ! lateral momentum diffusion       (dyn_ldf routine) 
    42    USE dynzdf           ! vertical diffusion               (dyn_zdf routine) 
    43    USE dynspg           ! surface pressure gradient        (dyn_spg routine) 
     37   USE divhor          ! horizontal divergence            (div_hor routine) 
     38   USE dynadv          ! advection                        (dyn_adv routine) 
     39   USE dynbfr          ! Bottom friction terms            (dyn_bfr routine) 
     40   USE dynvor          ! vorticity term                   (dyn_vor routine) 
     41   USE dynhpg          ! hydrostatic pressure grad.       (dyn_hpg routine) 
     42   USE dynldf          ! lateral momentum diffusion       (dyn_ldf routine) 
     43   USE dynzdf          ! vertical diffusion               (dyn_zdf routine) 
     44   USE dynspg          ! surface pressure gradient        (dyn_spg routine) 
    4445 
    45    USE dynnxt           ! time-stepping                    (dyn_nxt routine) 
     46   USE dynnxt          ! time-stepping                    (dyn_nxt routine) 
    4647 
    47    USE stopar           ! Stochastic parametrization       (sto_par routine) 
     48   USE stopar          ! Stochastic parametrization       (sto_par routine) 
    4849   USE stopts  
    4950 
    50    USE bdy_oce    , ONLY: ln_bdy 
    51    USE bdydta           ! open boundary condition data     (bdy_dta routine) 
    52    USE bdytra           ! bdy cond. for tracers            (bdy_tra routine) 
    53    USE bdydyn3d         ! bdy cond. for baroclinic vel.  (bdy_dyn3d routine) 
     51   USE bdy_oce  , ONLY : ln_bdy 
     52   USE bdydta          ! open boundary condition data     (bdy_dta routine) 
     53   USE bdytra          ! bdy cond. for tracers            (bdy_tra routine) 
     54   USE bdydyn3d        ! bdy cond. for baroclinic vel.  (bdy_dyn3d routine) 
    5455 
    55    USE sshwzv           ! vertical velocity and ssh        (ssh_nxt routine) 
     56   USE sshwzv          ! vertical velocity and ssh        (ssh_nxt routine) 
    5657   !                                                       (ssh_swp routine) 
    5758   !                                                       (wzv     routine) 
    58    USE domvvl           ! variable vertical scale factors  (dom_vvl_sf_nxt routine) 
     59   USE domvvl          ! variable vertical scale factors  (dom_vvl_sf_nxt routine) 
    5960   !                                                       (dom_vvl_sf_swp routine) 
    6061 
    61    USE ldfslp           ! iso-neutral slopes               (ldf_slp routine) 
    62    USE ldfdyn           ! lateral eddy viscosity coef.     (ldf_dyn routine) 
    63    USE ldftra           ! lateral eddy diffusive coef.     (ldf_tra routine) 
     62   USE ldfslp          ! iso-neutral slopes               (ldf_slp routine) 
     63   USE ldfdyn          ! lateral eddy viscosity coef.     (ldf_dyn routine) 
     64   USE ldftra          ! lateral eddy diffusive coef.     (ldf_tra routine) 
    6465 
    65    USE zdfphy         ! vertical physics manager      (zdf_phy_init routine) 
     66   USE zdfphy          ! vertical physics manager      (zdf_phy_init routine) 
    6667 
    6768   USE step_diu        ! Time stepping for diurnal sst 
     
    7071   USE sbc_oce         ! surface fluxes   
    7172    
    72    USE zpshde           ! partial step: hor. derivative     (zps_hde routine) 
     73   USE zpshde          ! partial step: hor. derivative     (zps_hde routine) 
    7374 
    74    USE diawri           ! Standard run outputs             (dia_wri routine) 
    75    USE diaptr           ! poleward transports              (dia_ptr routine) 
    76    USE diadct           ! sections transports              (dia_dct routine) 
    77    USE diaar5           ! AR5 diagnosics                   (dia_ar5 routine) 
    78    USE diahth           ! thermocline depth                (dia_hth routine) 
    79    USE diahsb           ! heat, salt and volume budgets    (dia_hsb routine) 
     75   USE diawri          ! Standard run outputs             (dia_wri routine) 
     76   USE diaptr          ! poleward transports              (dia_ptr routine) 
     77   USE diadct          ! sections transports              (dia_dct routine) 
     78   USE diaar5          ! AR5 diagnosics                   (dia_ar5 routine) 
     79   USE diahth          ! thermocline depth                (dia_hth routine) 
     80   USE diahsb          ! heat, salt and volume budgets    (dia_hsb routine) 
    8081   USE diaharm 
    8182   USE diacfl 
    82    USE flo_oce          ! floats variables 
    83    USE floats           ! floats computation               (flo_stp routine) 
     83   USE flo_oce         ! floats variables 
     84   USE floats          ! floats computation               (flo_stp routine) 
    8485 
    85    USE crsfld           ! Standard output on coarse grid   (crs_fld routine) 
     86   USE crsfld          ! Standard output on coarse grid   (crs_fld routine) 
    8687 
    87    USE asminc           ! assimilation increments      (tra_asm_inc routine) 
     88   USE asminc          ! assimilation increments      (tra_asm_inc routine) 
    8889   !                                                   (dyn_asm_inc routine) 
    8990   USE asmbkg 
    90    USE stpctl           ! time stepping control            (stp_ctl routine) 
    91    USE restart          ! ocean restart                    (rst_wri routine) 
    92    USE prtctl           ! Print control                    (prt_ctl routine) 
     91   USE stpctl          ! time stepping control            (stp_ctl routine) 
     92   USE restart         ! ocean restart                    (rst_wri routine) 
     93   USE prtctl          ! Print control                    (prt_ctl routine) 
    9394 
    94    USE diaobs           ! Observation operator 
     95   USE diaobs          ! Observation operator 
    9596 
    96    USE in_out_manager   ! I/O manager 
    97    USE iom              ! 
     97   USE in_out_manager  ! I/O manager 
     98   USE iom             ! 
    9899   USE lbclnk 
    99    USE timing           ! Timing 
     100   USE timing          ! Timing 
    100101 
    101102#if defined key_iomput 
    102    USE xios 
     103   USE xios            ! I/O server 
    103104#endif 
    104105#if defined key_agrif 
Note: See TracChangeset for help on using the changeset viewer.