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

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 7953 for branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC – NEMO

Ignore:
Timestamp:
2017-04-23T09:30:41+02:00 (7 years ago)
Author:
gm
Message:

#1880 (HPC-09): add zdfphy (the ZDF manager) + remove all key_...

Location:
branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC
Files:
1 deleted
17 edited
1 moved

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ASM/asmbkg.F90

    r6140 r7953  
    2828   USE ldfslp             ! Lateral diffusion: slopes of neutral surfaces 
    2929   USE tradmp             ! Tracer damping 
    30 #if defined key_zdftke 
    3130   USE zdftke             ! TKE vertical physics 
    32 #endif 
    3331   USE eosbn2             ! Equation of state (eos_bn2 routine) 
    3432   USE zdfmxl             ! Mixed layer depth 
     
    9492            IF( nitbkg_r == nit000 - 1 ) THEN      ! Treat special case when nitbkg = 0 
    9593               zdate = REAL( ndastp ) 
    96 #if defined key_zdftke 
    97                ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
    98                IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...' 
    99                CALL tke_rst( nit000, 'READ' )               ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
    100  
    101 #endif 
     94               IF( ln_zdftke ) THEN                   ! read turbulent kinetic energy ( en ) 
     95                  IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...' 
     96                  CALL tke_rst( nit000, 'READ' ) 
     97               ENDIF 
    10298            ELSE 
    10399               zdate = REAL( ndastp ) 
     
    111107            CALL iom_rstput( kt, nitbkg_r, inum, 'sn'     , tsn(:,:,:,jp_sal) ) 
    112108            CALL iom_rstput( kt, nitbkg_r, inum, 'sshn'   , sshn              ) 
    113 #if defined key_zdftke 
    114             CALL iom_rstput( kt, nitbkg_r, inum, 'en'     , en                ) 
    115 #endif 
     109            IF( ln_zdftke )   CALL iom_rstput( kt, nitbkg_r, inum, 'en'     , en                ) 
    116110            ! 
    117111            CALL iom_close( inum ) 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/C1D/step_c1d.F90

    r7931 r7953  
    7474                         CALL bn2( tsb, rab_b, rn2b ) ! before Brunt-Vaisala frequency 
    7575                         CALL bn2( tsn, rab_n, rn2  ) ! now    Brunt-Vaisala frequency 
    76       !  VERTICAL PHYSICS    
    77                          CALL zdf_bfr( kstp )         ! bottom friction 
    78       !                                               ! Vertical eddy viscosity and diffusivity coefficients 
    79       IF( lk_zdfric  )   CALL zdf_ric( kstp )            ! Richardson number dependent Kz 
    80       IF( lk_zdftke  )   CALL zdf_tke( kstp )            ! TKE closure scheme for Kz 
    81       IF( lk_zdfgls  )   CALL zdf_gls( kstp )            ! GLS closure scheme for Kz 
    82       IF( lk_zdfcst  )   THEN                            ! Constant Kz (reset avt, avm[uv] to the background value) 
    83          avt (:,:,:) = rn_avt0 * tmask(:,:,:) 
    84          avmu(:,:,:) = rn_avm0 * umask(:,:,:) 
    85          avmv(:,:,:) = rn_avm0 * vmask(:,:,:) 
    86       ENDIF 
    87  
    88       IF( ln_rnf_mouth ) THEN                         ! increase diffusivity at rivers mouths 
    89          DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2.e0 * rn_avt_rnf * rnfmsk(:,:)   ;   END DO 
    90       ENDIF 
    91       IF( ln_zdfevd  )   CALL zdf_evd( kstp )         ! enhanced vertical eddy diffusivity 
    92       IF( lk_zdftmx  )   CALL zdf_tmx( kstp )         ! tidal vertical mixing 
    93       IF( ln_zdfddm  )   CALL zdf_ddm( kstp )         ! double diffusive mixing 
    94                          CALL zdf_mxl( kstp )         ! mixed layer depth 
    95  
    96                                                       ! write tke information in the restart file 
    97       IF( lrst_oce .AND. lk_zdftke )   CALL tke_rst( kstp, 'WRITE' ) 
    98                                                       ! write gls information in the restart file 
    99       IF( lrst_oce .AND. lk_zdfgls )   CALL gls_rst( kstp, 'WRITE' ) 
     76       
     77      !  VERTICAL PHYSICS 
     78                         CALL zdf_phy( kstp )         ! vertical physics update (bfr, avt, avs, avm + MLD) 
    10079 
    10180      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/CRS/crs.F90

    r6140 r7953  
    151151 
    152152      ! Vertical diffusion 
    153       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)  ::  avt_crs           !: vert. diffusivity coef. [m2/s] at w-point for temp   
    154 # if defined key_zdfddm 
    155153      REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)  ::  avs_crs           !: salinity vertical diffusivity coeff. [m2/s] at w-point 
    156 # endif 
    157154 
    158155      ! Mixing and Mixed Layer Depth 
     
    239236         &     fr_i_crs(jpi_crs,jpj_crs), sfx_crs(jpi_crs ,jpj_crs),  STAT=ierr(12)  ) 
    240237 
    241      ALLOCATE( tsn_crs(jpi_crs,jpj_crs,jpk,jpts), avt_crs(jpi_crs,jpj_crs,jpk),    & 
    242 # if defined key_zdfddm 
    243          &      avs_crs(jpi_crs,jpj_crs,jpk),    & 
    244 # endif 
    245          &      STAT=ierr(13) ) 
     238     ALLOCATE( tsn_crs(jpi_crs,jpj_crs,jpk,jpts), avs_crs(jpi_crs,jpj_crs,jpk), STAT=ierr(13) ) 
    246239 
    247240      ALLOCATE( nmln_crs(jpi_crs,jpj_crs) , hmld_crs(jpi_crs,jpj_crs) , & 
    248241         &      hmlp_crs(jpi_crs,jpj_crs) , hmlpt_crs(jpi_crs,jpj_crs) , STAT=ierr(14) ) 
    249242          
    250       ALLOCATE( nimppt_crs(jpnij) , nlcit_crs(jpnij) , nldit_crs(jpnij) , nleit_crs(jpnij), & 
    251        &  nimppt_full(jpnij) , nlcit_full(jpnij) , nldit_full(jpnij) , nleit_full(jpnij),   & 
    252                 njmppt_crs(jpnij) , nlcjt_crs(jpnij) , nldjt_crs(jpnij) , nlejt_crs(jpnij), & 
    253        &  njmppt_full(jpnij) , nlcjt_full(jpnij) , nldjt_full(jpnij) , nlejt_full(jpnij)  , STAT=ierr(15) ) 
    254  
    255           
     243      ALLOCATE( nimppt_crs (jpnij) , nlcit_crs (jpnij) , nldit_crs (jpnij) , nleit_crs (jpnij),   & 
     244         &      nimppt_full(jpnij) , nlcit_full(jpnij) , nldit_full(jpnij) , nleit_full(jpnij),   & 
     245                njmppt_crs (jpnij) , nlcjt_crs (jpnij) , nldjt_crs (jpnij) , nlejt_crs (jpnij),   & 
     246         &      njmppt_full(jpnij) , nlcjt_full(jpnij) , nldjt_full(jpnij) , nlejt_full(jpnij)  , STAT=ierr(15) ) 
     247    
    256248      crs_dom_alloc = MAXVAL(ierr) 
    257  
     249      ! 
    258250   END FUNCTION crs_dom_alloc 
    259251 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/CRS/crsfld.F90

    r6140 r7953  
    8484         vn_crs   (:,:,:  ) = 0._wp    ! v-velocity 
    8585         wn_crs   (:,:,:  ) = 0._wp    ! w 
    86          avt_crs  (:,:,:  ) = 0._wp    ! avt 
     86         avs_crs  (:,:,:  ) = 0._wp    ! avt 
    8787         hdivn_crs(:,:,:  ) = 0._wp    ! hdiv 
    8888         rke_crs  (:,:,:  ) = 0._wp    ! rke 
     
    200200      SELECT CASE ( nn_crs_kz ) 
    201201         CASE ( 0 ) 
    202             CALL crs_dom_ope( avt, 'VOL', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     202            CALL crs_dom_ope( avt, 'VOL', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
    203203         CASE ( 1 ) 
    204             CALL crs_dom_ope( avt, 'MAX', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     204            CALL crs_dom_ope( avt, 'MAX', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
    205205         CASE ( 2 ) 
    206             CALL crs_dom_ope( avt, 'MIN', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
     206            CALL crs_dom_ope( avt, 'MIN', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 ) 
    207207      END SELECT 
    208208      ! 
    209       CALL iom_put( "avt", avt_crs )   !  Kz 
     209      CALL iom_put( "avt", avs_crs )   !  Kz 
    210210       
    211211      !  sbc fields   
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90

    r7753 r7953  
    88   USE oce             ! ocean dynamics and tracers variables 
    99   USE dom_oce         ! ocean space and time domain 
     10   USE zdf_oce         ! ocean vertical physics     
     11   USE zdfgls, ONLY: mxln 
    1012   USE in_out_manager  ! I/O units 
    1113   USE iom             ! I/0 library 
    12    USE wrk_nemo        ! working arrays 
    13 #if defined key_zdftke  
    14    USE zdf_oce, ONLY: en 
    15 #endif 
    16    USE zdf_oce, ONLY: avt, avm 
    17 #if defined key_zdfgls 
    18    USE zdf_oce, ONLY: en 
    19    USE zdfgls, ONLY: mxln 
    20 #endif 
     14   USE wrk_nemo        ! work arrays 
    2115 
    2216   IMPLICIT NONE 
    2317   PRIVATE 
    2418 
    25    LOGICAL , PUBLIC ::   ln_dia25h     !:  25h mean output 
    2619   PUBLIC   dia_25h_init               ! routine called by nemogcm.F90 
    2720   PUBLIC   dia_25h                    ! routine called by diawri.F90 
    2821 
    29   !! * variables for calculating 25-hourly means 
    30    REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   tn_25h  , sn_25h 
    31    REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:)   ::   sshn_25h  
    32    REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   un_25h  , vn_25h  , wn_25h 
    33    REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   avt_25h , avm_25h 
    34 #if defined key_zdfgls || key_zdftke 
    35    REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   en_25h 
    36 #endif 
    37 #if defined key_zdfgls  
    38    REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   rmxln_25h 
    39 #endif 
    40    INTEGER, SAVE :: cnt_25h     ! Counter for 25 hour means 
    41  
    42  
     22   LOGICAL, PUBLIC ::   ln_dia25h      !:  25h mean output 
     23 
     24   ! variables for calculating 25-hourly means 
     25   INTEGER , SAVE ::   cnt_25h     ! Counter for 25 hour means 
     26   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   tn_25h  , sn_25h 
     27   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:)   ::   sshn_25h  
     28   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   un_25h  , vn_25h  , wn_25h 
     29   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   avt_25h , avm_25h 
     30   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   en_25h  , rmxln_25h 
    4331 
    4432   !!---------------------------------------------------------------------- 
     
    5644      !!         
    5745      !! ** Method : Read namelist 
    58       !!   History 
    59       !!   3.6  !  08-14  (E. O'Dea) Routine to initialize dia_25h 
    6046      !!--------------------------------------------------------------------------- 
    61       !! 
    6247      INTEGER ::   ios                 ! Local integer output status for namelist read 
    6348      INTEGER ::   ierror              ! Local integer for memory allocation 
     
    7964         WRITE(numout,*) 'dia_25h_init : Output 25 hour mean diagnostics' 
    8065         WRITE(numout,*) '~~~~~~~~~~~~' 
    81          WRITE(numout,*) 'Namelist nam_dia25h : set 25h outputs ' 
    82          WRITE(numout,*) 'Switch for 25h diagnostics (T) or not (F)  ln_dia25h  = ', ln_dia25h 
     66         WRITE(numout,*) '   Namelist nam_dia25h : set 25h outputs ' 
     67         WRITE(numout,*) '      Switch for 25h diagnostics (T) or not (F)  ln_dia25h  = ', ln_dia25h 
    8368      ENDIF 
    8469      IF( .NOT. ln_dia25h )   RETURN 
     
    8671      ! 1 - Allocate memory ! 
    8772      ! ------------------- ! 
    88       ALLOCATE( tn_25h(jpi,jpj,jpk), STAT=ierror ) 
     73      !                                ! ocean arrays 
     74      ALLOCATE( tn_25h (jpi,jpj,jpk), sn_25h (jpi,jpj,jpk), sshn_25h(jpi,jpj)  ,     & 
     75         &      un_25h (jpi,jpj,jpk), vn_25h (jpi,jpj,jpk), wn_25h(jpi,jpj,jpk),     & 
     76         &      avt_25h(jpi,jpj,jpk), avm_25h(jpi,jpj,jpk),                      STAT=ierror ) 
    8977      IF( ierror > 0 ) THEN 
    90          CALL ctl_stop( 'dia_25h: unable to allocate tn_25h' )   ;   RETURN 
    91       ENDIF 
    92       ALLOCATE( sn_25h(jpi,jpj,jpk), STAT=ierror ) 
    93       IF( ierror > 0 ) THEN 
    94          CALL ctl_stop( 'dia_25h: unable to allocate sn_25h' )   ;   RETURN 
    95       ENDIF 
    96       ALLOCATE( un_25h(jpi,jpj,jpk), STAT=ierror ) 
    97       IF( ierror > 0 ) THEN 
    98          CALL ctl_stop( 'dia_25h: unable to allocate un_25h' )   ;   RETURN 
    99       ENDIF 
    100       ALLOCATE( vn_25h(jpi,jpj,jpk), STAT=ierror ) 
    101       IF( ierror > 0 ) THEN 
    102          CALL ctl_stop( 'dia_25h: unable to allocate vn_25h' )   ;   RETURN 
    103       ENDIF 
    104       ALLOCATE( wn_25h(jpi,jpj,jpk), STAT=ierror ) 
    105       IF( ierror > 0 ) THEN 
    106          CALL ctl_stop( 'dia_25h: unable to allocate wn_25h' )   ;   RETURN 
    107       ENDIF 
    108       ALLOCATE( avt_25h(jpi,jpj,jpk), STAT=ierror ) 
    109       IF( ierror > 0 ) THEN 
    110          CALL ctl_stop( 'dia_25h: unable to allocate avt_25h' )   ;   RETURN 
    111       ENDIF 
    112       ALLOCATE( avm_25h(jpi,jpj,jpk), STAT=ierror ) 
    113       IF( ierror > 0 ) THEN 
    114          CALL ctl_stop( 'dia_25h: unable to allocate avm_25h' )   ;   RETURN 
    115       ENDIF 
    116 # if defined key_zdfgls || defined key_zdftke 
    117       ALLOCATE( en_25h(jpi,jpj,jpk), STAT=ierror ) 
    118       IF( ierror > 0 ) THEN 
    119          CALL ctl_stop( 'dia_25h: unable to allocate en_25h' )   ;   RETURN 
    120       ENDIF 
    121 #endif 
    122 # if defined key_zdfgls  
    123       ALLOCATE( rmxln_25h(jpi,jpj,jpk), STAT=ierror ) 
    124       IF( ierror > 0 ) THEN 
    125          CALL ctl_stop( 'dia_25h: unable to allocate rmxln_25h' )   ;   RETURN 
    126       ENDIF 
    127 #endif 
    128       ALLOCATE( sshn_25h(jpi,jpj), STAT=ierror ) 
    129       IF( ierror > 0 ) THEN 
    130          CALL ctl_stop( 'dia_25h: unable to allocate sshn_25h' )   ;   RETURN 
     78         CALL ctl_stop( 'dia_25h: unable to allocate ocean arrays' )   ;   RETURN 
     79      ENDIF 
     80      IF( ln_zdftke ) THEN             ! TKE physics 
     81         ALLOCATE( en_25h(jpi,jpj,jpk), STAT=ierror ) 
     82         IF( ierror > 0 ) THEN 
     83            CALL ctl_stop( 'dia_25h: unable to allocate en_25h' )   ;   RETURN 
     84         ENDIF 
     85      ENDIF 
     86      IF( ln_zdfgls ) THEN             ! GLS physics 
     87         ALLOCATE( en_25h(jpi,jpj,jpk), rmxln_25h(jpi,jpj,jpk), STAT=ierror ) 
     88         IF( ierror > 0 ) THEN 
     89            CALL ctl_stop( 'dia_25h: unable to allocate en_25h and rmxln_25h' )   ;   RETURN 
     90         ENDIF 
    13191      ENDIF 
    13292      ! ------------------------- ! 
     
    142102      avt_25h(:,:,:) = avt(:,:,:) 
    143103      avm_25h(:,:,:) = avm(:,:,:) 
    144 # if defined key_zdfgls || defined key_zdftke 
     104      IF( ln_zdftke ) THEN 
    145105         en_25h(:,:,:) = en(:,:,:) 
    146 #endif 
    147 # if defined key_zdfgls 
     106      ENDIF 
     107      IF( ln_zdfgls ) THEN 
     108         en_25h(:,:,:) = en(:,:,:) 
    148109         rmxln_25h(:,:,:) = mxln(:,:,:) 
    149 #endif 
     110      ENDIF 
    150111#if defined key_lim3 || defined key_lim2 
    151112         CALL ctl_stop('STOP', 'dia_25h not setup yet to do tidemean ice') 
    152113#endif  
    153  
    154       ! -------------------------- ! 
    155       ! 3 - Return to dia_wri      ! 
    156       ! -------------------------- ! 
    157  
    158  
     114      ! 
    159115   END SUBROUTINE dia_25h_init 
    160116 
     
    164120      !!                 ***  ROUTINE dia_25h  *** 
    165121      !!          
    166       !! 
    167       !!-------------------------------------------------------------------- 
    168       !!                    
    169122      !! ** Purpose :   Write diagnostics with M2/S2 tide removed 
    170123      !! 
    171       !! ** Method  :    
    172       !!      25hr mean outputs for shelf seas 
     124      !! ** Method  :   25hr mean outputs for shelf seas 
     125      !!---------------------------------------------------------------------- 
     126      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    173127      !! 
    174       !! History : 
    175       !!   ?.0  !  07-04  (A. Hines) New routine, developed from dia_wri_foam 
    176       !!   3.4  !  02-13  (J. Siddorn) Routine taken from old dia_wri_foam 
    177       !!   3.6  !  08-14  (E. O'Dea) adapted for VN3.6 
    178       !!---------------------------------------------------------------------- 
    179       !! * Modules used 
    180  
    181       IMPLICIT NONE 
    182  
    183       !! * Arguments 
    184       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    185  
    186  
    187       !! * Local declarations 
    188128      INTEGER ::   ji, jj, jk 
    189  
     129      INTEGER                          ::   iyear0, nimonth0,iday0            ! start year,imonth,day 
    190130      LOGICAL ::   ll_print = .FALSE.    ! =T print and flush numout 
    191       REAL(wp)                         ::   zsto, zout, zmax, zjulian, zmdi       ! temporary reals 
    192       INTEGER                          ::   i_steps                               ! no of timesteps per hour 
    193       REAL(wp), DIMENSION(jpi,jpj    ) ::   zw2d, un_dm, vn_dm                    ! temporary workspace 
    194       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zw3d                                  ! temporary workspace 
    195       REAL(wp), DIMENSION(jpi,jpj,3)   ::   zwtmb                                 ! temporary workspace 
    196       INTEGER                          ::   iyear0, nimonth0,iday0                ! start year,imonth,day 
    197  
     131      REAL(wp)                         ::   zsto, zout, zmax, zjulian, zmdi   ! local scalars 
     132      INTEGER                          ::   i_steps                           ! no of timesteps per hour 
     133      REAL(wp), DIMENSION(jpi,jpj    ) ::   zw2d, un_dm, vn_dm                ! workspace 
     134      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zw3d                              ! workspace 
     135      REAL(wp), DIMENSION(jpi,jpj,3)   ::   zwtmb                             ! workspace 
    198136      !!---------------------------------------------------------------------- 
    199137 
     
    207145      ENDIF 
    208146 
    209 #if defined key_lim3 || defined key_lim2 
    210       CALL ctl_stop('STOP', 'dia_wri_tide not setup yet to do tidemean ice') 
    211 #endif 
    212  
    213147      ! local variable for debugging 
    214148      ll_print = ll_print .AND. lwp 
    215149 
    216       ! Sum of 25 hourly instantaneous values to give a 25h mean from 24hours 
    217       ! every day 
    218       IF( MOD( kt, i_steps ) == 0  .and. kt .ne. nn_it000 ) THEN 
     150      ! Sum of 25 hourly instantaneous values to give a 25h mean from 24hours every day 
     151      IF( MOD( kt, i_steps ) == 0  .AND. kt /= nn_it000 ) THEN 
    219152 
    220153         IF (lwp) THEN 
     
    231164         avt_25h(:,:,:)       = avt_25h(:,:,:) + avt(:,:,:) 
    232165         avm_25h(:,:,:)       = avm_25h(:,:,:) + avm(:,:,:) 
    233 # if defined key_zdfgls || defined key_zdftke 
    234          en_25h(:,:,:)        = en_25h(:,:,:) + en(:,:,:) 
    235 #endif 
    236 # if defined key_zdfgls 
    237          rmxln_25h(:,:,:)      = rmxln_25h(:,:,:) + mxln(:,:,:) 
    238 #endif 
     166         IF( ln_zdftke ) THEN 
     167            en_25h(:,:,:)        = en_25h(:,:,:) + en(:,:,:) 
     168         ENDIF 
     169         IF( ln_zdfgls ) THEN 
     170            en_25h(:,:,:)        = en_25h(:,:,:) + en(:,:,:) 
     171            rmxln_25h(:,:,:)      = rmxln_25h(:,:,:) + mxln(:,:,:) 
     172         ENDIF 
    239173         cnt_25h = cnt_25h + 1 
    240  
     174         ! 
    241175         IF (lwp) THEN 
    242176            WRITE(numout,*) 'dia_tide : Summed the following number of hourly values so far',cnt_25h 
    243177         ENDIF 
    244  
     178         ! 
    245179      ENDIF ! MOD( kt, i_steps ) == 0 
    246180 
    247          ! Write data for 25 hour mean output streams 
    248       IF( cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps*24) == 0 .AND. kt .NE. nn_it000 ) THEN 
    249  
    250             IF(lwp) THEN 
    251                WRITE(numout,*) 'dia_wri_tide : Writing 25 hour mean tide diagnostics at timestep', kt 
    252                WRITE(numout,*) '~~~~~~~~~~~~ ' 
    253             ENDIF 
    254  
    255             tn_25h(:,:,:)        = tn_25h(:,:,:) / 25.0_wp 
    256             sn_25h(:,:,:)        = sn_25h(:,:,:) / 25.0_wp 
    257             sshn_25h(:,:)        = sshn_25h(:,:) / 25.0_wp 
    258             un_25h(:,:,:)        = un_25h(:,:,:) / 25.0_wp 
    259             vn_25h(:,:,:)        = vn_25h(:,:,:) / 25.0_wp 
    260             wn_25h(:,:,:)        = wn_25h(:,:,:) / 25.0_wp 
    261             avt_25h(:,:,:)       = avt_25h(:,:,:) / 25.0_wp 
    262             avm_25h(:,:,:)       = avm_25h(:,:,:) / 25.0_wp 
    263 # if defined key_zdfgls || defined key_zdftke 
     181      ! Write data for 25 hour mean output streams 
     182      IF( cnt_25h == 25 .AND.  MOD( kt, i_steps*24) == 0 .AND. kt /= nn_it000 ) THEN 
     183         ! 
     184         IF(lwp) THEN 
     185            WRITE(numout,*) 'dia_wri_tide : Writing 25 hour mean tide diagnostics at timestep', kt 
     186            WRITE(numout,*) '~~~~~~~~~~~~ ' 
     187         ENDIF 
     188         ! 
     189         tn_25h(:,:,:)        = tn_25h(:,:,:) / 25.0_wp 
     190         sn_25h(:,:,:)        = sn_25h(:,:,:) / 25.0_wp 
     191         sshn_25h(:,:)        = sshn_25h(:,:) / 25.0_wp 
     192         un_25h(:,:,:)        = un_25h(:,:,:) / 25.0_wp 
     193         vn_25h(:,:,:)        = vn_25h(:,:,:) / 25.0_wp 
     194         wn_25h(:,:,:)        = wn_25h(:,:,:) / 25.0_wp 
     195         avt_25h(:,:,:)       = avt_25h(:,:,:) / 25.0_wp 
     196         avm_25h(:,:,:)       = avm_25h(:,:,:) / 25.0_wp 
     197         IF( ln_zdftke ) THEN 
    264198            en_25h(:,:,:)        = en_25h(:,:,:) / 25.0_wp 
    265 #endif 
    266 # if defined key_zdfgls 
     199         ENDIF 
     200         IF( ln_zdfgls ) THEN 
     201            en_25h(:,:,:)        = en_25h(:,:,:) / 25.0_wp 
    267202            rmxln_25h(:,:,:)       = rmxln_25h(:,:,:) / 25.0_wp 
    268 #endif 
    269  
    270             IF (lwp)  WRITE(numout,*) 'dia_wri_tide : Mean calculated by dividing 25 hour sums and writing output' 
    271             zmdi=1.e+20 !missing data indicator for masking 
    272             ! write tracers (instantaneous) 
    273             zw3d(:,:,:) = tn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    274             CALL iom_put("temper25h", zw3d)   ! potential temperature 
    275             zw3d(:,:,:) = sn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    276             CALL iom_put( "salin25h", zw3d  )   ! salinity 
    277             zw2d(:,:) = sshn_25h(:,:)*tmask(:,:,1) + zmdi*(1.0-tmask(:,:,1)) 
    278             CALL iom_put( "ssh25h", zw2d )   ! sea surface  
    279  
    280  
    281             ! Write velocities (instantaneous) 
    282             zw3d(:,:,:) = un_25h(:,:,:)*umask(:,:,:) + zmdi*(1.0-umask(:,:,:)) 
    283             CALL iom_put("vozocrtx25h", zw3d)    ! i-current 
    284             zw3d(:,:,:) = vn_25h(:,:,:)*vmask(:,:,:) + zmdi*(1.0-vmask(:,:,:)) 
    285             CALL iom_put("vomecrty25h", zw3d  )   ! j-current 
    286  
    287             zw3d(:,:,:) = wn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    288             CALL iom_put("vomecrtz25h", zw3d )   ! k-current 
    289             zw3d(:,:,:) = avt_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    290             CALL iom_put("avt25h", zw3d )   ! diffusivity 
    291             zw3d(:,:,:) = avm_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    292             CALL iom_put("avm25h", zw3d)   ! viscosity 
    293 #if defined key_zdftke || defined key_zdfgls  
     203         ENDIF 
     204         ! 
     205         IF(lwp)  WRITE(numout,*) 'dia_wri_tide : Mean calculated by dividing 25 hour sums and writing output' 
     206         zmdi=1.e+20 !missing data indicator for masking 
     207         ! write tracers (instantaneous) 
     208         zw3d(:,:,:) = tn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
     209         CALL iom_put("temper25h", zw3d)   ! potential temperature 
     210         zw3d(:,:,:) = sn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
     211         CALL iom_put( "salin25h", zw3d  )   ! salinity 
     212         zw2d(:,:) = sshn_25h(:,:)*tmask(:,:,1) + zmdi*(1.0-tmask(:,:,1)) 
     213         CALL iom_put( "ssh25h", zw2d )   ! sea surface  
     214         ! Write velocities (instantaneous) 
     215         zw3d(:,:,:) = un_25h(:,:,:)*umask(:,:,:) + zmdi*(1.0-umask(:,:,:)) 
     216         CALL iom_put("vozocrtx25h", zw3d)    ! i-current 
     217         zw3d(:,:,:) = vn_25h(:,:,:)*vmask(:,:,:) + zmdi*(1.0-vmask(:,:,:)) 
     218         CALL iom_put("vomecrty25h", zw3d  )   ! j-current 
     219         zw3d(:,:,:) = wn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
     220         CALL iom_put("vomecrtz25h", zw3d )   ! k-current 
     221         ! Write vertical physics 
     222         zw3d(:,:,:) = avt_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
     223         CALL iom_put("avt25h", zw3d )   ! diffusivity 
     224         zw3d(:,:,:) = avm_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
     225         CALL iom_put("avm25h", zw3d)   ! viscosity 
     226         IF( ln_zdftke ) THEN 
    294227            zw3d(:,:,:) = en_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    295228            CALL iom_put("tke25h", zw3d)   ! tke 
    296 #endif 
    297 #if defined key_zdfgls  
     229         ENDIF 
     230         IF( ln_zdfgls ) THEN 
     231            zw3d(:,:,:) = en_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
     232            CALL iom_put("tke25h", zw3d)   ! tke 
    298233            zw3d(:,:,:) = rmxln_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 
    299234            CALL iom_put( "mxln25h",zw3d) 
    300 #endif 
    301  
    302             ! After the write reset the values to cnt=1 and sum values equal current value  
    303             tn_25h(:,:,:) = tsn(:,:,:,jp_tem) 
    304             sn_25h(:,:,:) = tsn(:,:,:,jp_sal) 
    305             sshn_25h(:,:) = sshn (:,:) 
    306             un_25h(:,:,:) = un(:,:,:) 
    307             vn_25h(:,:,:) = vn(:,:,:) 
    308             wn_25h(:,:,:) = wn(:,:,:) 
    309             avt_25h(:,:,:) = avt(:,:,:) 
    310             avm_25h(:,:,:) = avm(:,:,:) 
    311 # if defined key_zdfgls || defined key_zdftke 
     235         ENDIF 
     236         ! 
     237         ! After the write reset the values to cnt=1 and sum values equal current value  
     238         tn_25h(:,:,:) = tsn(:,:,:,jp_tem) 
     239         sn_25h(:,:,:) = tsn(:,:,:,jp_sal) 
     240         sshn_25h(:,:) = sshn (:,:) 
     241         un_25h(:,:,:) = un(:,:,:) 
     242         vn_25h(:,:,:) = vn(:,:,:) 
     243         wn_25h(:,:,:) = wn(:,:,:) 
     244         avt_25h(:,:,:) = avt(:,:,:) 
     245         avm_25h(:,:,:) = avm(:,:,:) 
     246         IF( ln_zdftke ) THEN 
    312247            en_25h(:,:,:) = en(:,:,:) 
    313 #endif 
    314 # if defined key_zdfgls 
     248         ENDIF 
     249         IF( ln_zdfgls ) THEN 
     250            en_25h(:,:,:) = en(:,:,:) 
    315251            rmxln_25h(:,:,:) = mxln(:,:,:) 
    316 #endif 
    317             cnt_25h = 1 
    318             IF (lwp)  WRITE(numout,*) 'dia_wri_tide : After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average',cnt_25h 
    319  
     252         ENDIF 
     253         cnt_25h = 1 
     254         IF (lwp)  WRITE(numout,*) 'dia_wri_tide : After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average',cnt_25h 
     255         ! 
    320256      ENDIF !  cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps * 24) == 0 .AND. kt .NE. nn_it000 
    321  
    322  
     257      ! 
    323258   END SUBROUTINE dia_25h  
    324259 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r7931 r7953  
    3939   USE zdfmxl          ! mixed layer 
    4040   USE dianam          ! build name of file (routine) 
    41    USE zdfddm          ! vertical  physics: double diffusion 
     41!   USE zdfddm          ! vertical  physics: double diffusion 
    4242   USE diahth          ! thermocline diagnostics 
    4343   USE wet_dry         ! wetting and drying 
     
    233233 
    234234      CALL iom_put( "avt" , avt                        )    ! T vert. eddy diff. coef. 
    235       CALL iom_put( "avs" , avs                        )    ! S vert. eddy diff. coef. (useful only with key_zdfddm) 
     235      CALL iom_put( "avs" , avs                        )    ! S vert. eddy diff. coef. 
    236236      CALL iom_put( "avm" , avmu                       )    ! T vert. eddy visc. coef. 
    237237 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf.F90

    r7753 r7953  
    9696      !!                 ***  ROUTINE dyn_zdf_init  *** 
    9797      !! 
    98       !! ** Purpose :   initializations of the vertical diffusion scheme 
     98      !! ** Purpose :   initialization of the vertical diffusion scheme 
    9999      !! 
    100100      !! ** Method  :   implicit (euler backward) scheme (default) 
     
    105105      !!---------------------------------------------------------------------- 
    106106      ! 
    107       ! Choice from ln_zdfexp read in namelist in zdfini 
     107      ! Choice from ln_zdfexp (namzdf namelist variable read in zdfphy module) 
    108108      IF( ln_zdfexp ) THEN   ;   nzdf = 0           ! use explicit scheme 
    109109      ELSE                   ;   nzdf = 1           ! use implicit scheme 
     
    111111      ! 
    112112      ! Force implicit schemes 
    113       IF( lk_zdftke .OR. lk_zdfgls   )   nzdf = 1   ! TKE or GLS physics 
     113      IF( ln_zdftke .OR. ln_zdfgls   )   nzdf = 1   ! TKE or GLS physics 
    114114      IF( ln_dynldf_iso              )   nzdf = 1   ! iso-neutral lateral physics 
    115115      IF( ln_dynldf_hor .AND. ln_sco )   nzdf = 1   ! horizontal lateral physics in s-coordinate 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/TRA/trazdf.F90

    r7931 r7953  
    124124      !!---------------------------------------------------------------------- 
    125125      ! 
    126       ! Choice from ln_zdfexp already read in namelist in zdfini module 
     126      ! Choice from ln_zdfexp (namzdf namelist variable read in zdfphy module) 
    127127      IF( ln_zdfexp ) THEN   ;   nzdf = 0           ! use explicit scheme 
    128128      ELSE                   ;   nzdf = 1           ! use implicit scheme 
     
    130130      ! 
    131131      ! Force implicit schemes 
    132       IF( lk_zdftke .OR. lk_zdfgls   )   nzdf = 1   ! TKE, or GLS physics 
     132      IF( ln_zdftke .OR. ln_zdfgls   )   nzdf = 1   ! TKE, or GLS physics 
    133133      IF( ln_traldf_iso              )   nzdf = 1   ! iso-neutral lateral physics 
    134134      IF( ln_traldf_hor .AND. ln_sco )   nzdf = 1   ! horizontal lateral physics in s-coordinate 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90

    r7931 r7953  
    1515 
    1616   PUBLIC  zdf_oce_alloc    ! Called in nemogcm.F90 
    17  
    18 #if defined key_zdfcst 
    19    LOGICAL, PARAMETER, PUBLIC ::   lk_zdfcst        = .TRUE.         !: constant vertical mixing flag 
    20 #else 
    21    LOGICAL, PARAMETER, PUBLIC ::   lk_zdfcst        = .FALSE.        !: constant vertical mixing flag 
    22 #endif 
    2317 
    2418   !                                 !!* namelist namzdf: vertical diffusion * 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r7646 r7953  
    88   !!             3.3  !  2010-10  (C. Bricaud)  Add in the reference 
    99   !!---------------------------------------------------------------------- 
    10 #if defined key_zdfgls 
    11    !!---------------------------------------------------------------------- 
    12    !!   'key_zdfgls'                 Generic Length Scale vertical physics 
     10 
    1311   !!---------------------------------------------------------------------- 
    1412   !!   zdf_gls       : update momentum and tracer Kz from a gls scheme 
     
    3937 
    4038   PUBLIC   zdf_gls        ! routine called in step module 
    41    PUBLIC   zdf_gls_init   ! routine called in opa module 
    42    PUBLIC   gls_rst        ! routine called in step module 
    43  
    44    LOGICAL , PUBLIC, PARAMETER ::   lk_zdfgls = .TRUE.   !: TKE vertical mixing flag 
     39   PUBLIC   zdf_gls_init   ! routine called in zdfphy module 
     40   PUBLIC   gls_rst        ! routine called in zdfphy module 
     41 
    4542   ! 
    4643   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   mxln    !: now mixing length 
     
    12111208   END SUBROUTINE gls_rst 
    12121209 
    1213 #else 
    1214    !!---------------------------------------------------------------------- 
    1215    !!   Dummy module :                                        NO TKE scheme 
    1216    !!---------------------------------------------------------------------- 
    1217    LOGICAL, PUBLIC, PARAMETER ::   lk_zdfgls = .FALSE.   !: TKE flag 
    1218 CONTAINS 
    1219    SUBROUTINE zdf_gls_init           ! Empty routine 
    1220       WRITE(*,*) 'zdf_gls_init: You should not have seen this print! error?' 
    1221    END SUBROUTINE zdf_gls_init 
    1222    SUBROUTINE zdf_gls( kt )          ! Empty routine 
    1223       WRITE(*,*) 'zdf_gls: You should not have seen this print! error?', kt 
    1224    END SUBROUTINE zdf_gls 
    1225    SUBROUTINE gls_rst( kt, cdrw )          ! Empty routine 
    1226       INTEGER         , INTENT(in) ::   kt         ! ocean time-step 
    1227       CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
    1228       WRITE(*,*) 'gls_rst: You should not have seen this print! error?', kt, cdrw 
    1229    END SUBROUTINE gls_rst 
    1230 #endif 
    1231  
    12321210   !!====================================================================== 
    12331211END MODULE zdfgls 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfphy.F90

    r7931 r7953  
    1 MODULE zdfini 
     1MODULE zdfphy 
    22   !!====================================================================== 
    3    !!                      ***  MODULE  zdfini  *** 
    4    !! Ocean physics :   read vertical mixing namelist and check consistancy 
     3   !!                      ***  MODULE  zdfphy  *** 
     4   !! Ocean physics :   manager of vertical mixing parametrizations 
    55   !!====================================================================== 
    6    !! History :  8.0  ! 1997-06  (G. Madec)  Original code from inimix 
    7    !!            1.0  ! 2002-08  (G. Madec)  F90 : free form 
    8    !!             -   ! 2005-06  (C. Ethe) KPP scheme 
    9    !!             -   ! 2009-07  (G. Madec) add avmb, avtb in restart for cen2 advection 
    10    !!            3.7  ! 2014-12  (G. Madec) remove KPP scheme 
    11    !!---------------------------------------------------------------------- 
    12  
    13    !!---------------------------------------------------------------------- 
    14    !!   zdf_init    : initialization, namelist read, and parameters control 
    15    !!---------------------------------------------------------------------- 
    16    USE par_oce         ! mesh and scale factors 
    17    USE zdf_oce         ! TKE vertical mixing           
    18    USE sbc_oce         ! surface module (only for nn_isf in the option compatibility test) 
    19    USE zdftke          ! TKE vertical mixing 
    20    USE zdfgls          ! GLS vertical mixing 
    21    USE zdfric          ! Richardson vertical mixing    
    22    USE zdfddm          ! double diffusion mixing       
    23    USE zdfevd          ! enhanced vertical diffusion   
    24    USE tranpc          ! convection: non penetrative adjustment 
    25    USE ldfslp          ! iso-neutral slopes 
     6   !! History :  4.0  !  2017-04  (G. Madec)  original code 
     7   !!---------------------------------------------------------------------- 
     8 
     9   !!---------------------------------------------------------------------- 
     10   !!   zdf_phy_init  : initialization of all vertical physics pakages 
     11   !!   zdf_phy       : upadate at each time-step the vertical mixing coeff.  
     12   !!---------------------------------------------------------------------- 
     13   USE par_oce        ! mesh and scale factors 
     14   USE zdf_oce        ! TKE vertical mixing           
     15   USE sbc_oce        ! surface module (only for nn_isf in the option compatibility test) 
     16   USE zdfbfr         ! bottom friction 
     17   USE zdftke         ! TKE vertical mixing 
     18   USE zdfgls         ! GLS vertical mixing 
     19   USE zdfric         ! Richardson vertical mixing    
     20   USE zdfddm         ! double diffusion mixing       
     21   USE zdfevd         ! enhanced vertical diffusion 
     22   USE zdftmx         ! internal tide-induced mixing 
     23   USE zdfqiao          !Qiao module wave induced mixing   (zdf_qiao routine) 
     24   USE zdfmxl           ! Mixed-layer depth                (zdf_mxl routine) 
     25   USE tranpc         ! convection: non penetrative adjustment 
     26   USE sbcrnf         ! surface boundary condition: runoff variables 
    2627   ! 
    27    USE in_out_manager  ! I/O manager 
    28    USE iom             ! IOM library 
    29    USE lib_mpp         ! distribued memory computing 
     28   USE in_out_manager ! I/O manager 
     29   USE iom            ! IOM library 
     30   USE lib_mpp        ! distribued memory computing 
    3031 
    3132   IMPLICIT NONE 
    3233   PRIVATE 
    3334 
    34    PUBLIC   zdf_init   ! routine called by opa.F90 
     35   PUBLIC   zdf_phy_init   ! routine called by nemogcm.F90 
     36   PUBLIC   zdf_phy        ! routine called by step.F90 
     37    
    3538    
    3639   !!---------------------------------------------------------------------- 
     
    4144CONTAINS 
    4245 
    43    SUBROUTINE zdf_init 
    44       !!---------------------------------------------------------------------- 
    45       !!                  ***  ROUTINE zdf_init  *** 
     46   SUBROUTINE zdf_phy_init 
     47      !!---------------------------------------------------------------------- 
     48      !!                  ***  ROUTINE zdf_phy_init  *** 
    4649      !!  
    4750      !! ** Purpose :   initializations of the vertical ocean physics 
     
    4952      !! ** Method  :   Read namelist namzdf, control logicals  
    5053      !!---------------------------------------------------------------------- 
    51       INTEGER ::   ioptio, ios       ! local integers 
     54      INTEGER ::   ioptio, ios   ! local integers 
    5255      !! 
    5356      NAMELIST/namzdf/ ln_zdfcst, ln_zdfric, ln_zdftke, ln_zdfgls,   &     ! type of closure scheme 
     
    7780      IF(lwp) THEN               !* Parameter print 
    7881         WRITE(numout,*) 
    79          WRITE(numout,*) 'zdf_init : vertical physics' 
     82         WRITE(numout,*) 'zdf_phy_init : vertical physics' 
    8083         WRITE(numout,*) '~~~~~~~~' 
    8184         WRITE(numout,*) '   Namelist namzdf : set vertical mixing mixing parameters' 
     
    106109      ENDIF 
    107110 
    108       IF(ln_zdfddm) THEN                    ! double diffusive mixing' 
    109          ALLOCATE( avs(jpi,jpj,jpk) ) 
    110          avs(:,:,:) = rn_avt0 * wmask(:,:,:)  
    111       ENDIF 
    112  
     111!!gm      IF(ln_zdfddm) THEN                    ! double diffusive mixing' 
     112!         avs(:,:,:) = rn_avt0 * wmask(:,:,:)  
     113!!gm      ENDIF 
    113114 
    114115      !                          !* Parameter & logical controls 
     
    122123      IF(lwp) WRITE(numout,*) '   vertical mixing option :' 
    123124      ioptio = 0 
    124       IF( lk_zdfcst ) THEN 
     125      IF( ln_zdfcst ) THEN 
    125126         IF(lwp) WRITE(numout,*) '      constant eddy diffusion coefficients' 
    126127         ioptio = ioptio+1 
    127128      ENDIF 
    128       IF( lk_zdfric ) THEN 
     129      IF( ln_zdfric ) THEN 
    129130         IF(lwp) WRITE(numout,*) '      Richardson dependent eddy coefficients' 
    130131         ioptio = ioptio+1 
    131132      ENDIF 
    132       IF( lk_zdftke ) THEN 
     133      IF( ln_zdftke ) THEN 
    133134         IF(lwp) WRITE(numout,*) '      TKE dependent eddy coefficients' 
    134135         ioptio = ioptio+1 
    135136      ENDIF 
    136       IF( lk_zdfgls ) THEN 
     137      IF( ln_zdfgls ) THEN 
    137138         IF(lwp) WRITE(numout,*) '      GLS dependent eddy coefficients' 
    138139         ioptio = ioptio+1 
     
    140141      IF( ioptio == 0 .OR. ioptio > 1 )   & 
    141142         &   CALL ctl_stop( ' one and only one vertical diffusion option has to be defined ' ) 
    142       IF( ( lk_zdfric .OR. lk_zdfgls ) .AND. ln_isfcav )   & 
     143      IF( ( ln_zdfric .OR. ln_zdfgls ) .AND. ln_isfcav )   & 
    143144         &   CALL ctl_stop( ' only zdfcst and zdftke were tested with ice shelves cavities ' ) 
    144145      ! 
     
    148149      ! 
    149150#if defined key_top 
    150       IF( ln_zdfnpc )   CALL ctl_stop( ' zdf_init: npc scheme is not working with key_top' ) 
     151      IF( ln_zdfnpc )   CALL ctl_stop( ' zdf_phy_init: npc scheme is not working with key_top' ) 
    151152#endif 
    152153      ! 
     
    160161         ioptio = ioptio+1 
    161162      ENDIF 
    162       IF( lk_zdftke ) THEN 
     163      IF( ln_zdftke ) THEN 
    163164         IF(lwp) WRITE(numout,*) '      use the 1.5 turbulent closure' 
    164165      ENDIF 
    165       IF( lk_zdfgls ) THEN 
     166      IF( ln_zdfgls ) THEN 
    166167         IF(lwp) WRITE(numout,*) '      use the GLS closure scheme' 
    167168      ENDIF 
    168169      IF ( ioptio > 1 )   CALL ctl_stop( ' chose between ln_zdfnpc and ln_zdfevd' ) 
    169       IF( ioptio == 0 .AND. .NOT.( lk_zdftke .OR. lk_zdfgls ) )           & 
     170      IF( ioptio == 0 .AND. .NOT.( ln_zdftke .OR. ln_zdfgls ) )           & 
    170171         CALL ctl_stop( ' except for TKE or GLS physics, a convection scheme is',   & 
    171172         &              ' required: ln_zdfevd or ln_zdfnpc logicals' ) 
     
    202203      ENDIF 
    203204      ! 
    204    END SUBROUTINE zdf_init 
     205 
     206!!gm moved into zdf_phy_init 
     207! 
     208                            CALL zdf_bfr_init      ! bottom friction 
     209 
     210      ioptio = 0                 !==  type of vertical turbulent closure  ==!    (set nzdfphy) 
     211      ! 
     212!      IF( ln_zdfcst ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_CST   ;   ENDIF 
     213!      IF( ln_zdfric ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_RIC   ;   CALL zdf_ric_init   ;   ENDIF 
     214!      IF( ln_zdftke ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_TKE   ;   CALL zdf_tke_init   ;   ENDIF 
     215!      IF( ln_zdfgls ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_GLS   ;   CALL zdf_gls_init   ;   ENDIF 
     216 
     217 
     218      ! 
     219      IF( ln_zdfric     )   CALL zdf_ric_init      ! Richardson number dependent Kz 
     220      IF( ln_zdftke     )   CALL zdf_tke_init      ! TKE closure scheme 
     221      IF( ln_zdfgls     )   CALL zdf_gls_init      ! GLS closure scheme 
     222      IF( ln_zdftmx     )   CALL zdf_tmx_init      ! tidal vertical mixing 
     223!!gm 
     224      ! 
     225   END SUBROUTINE zdf_phy_init 
     226 
     227 
     228   SUBROUTINE zdf_phy( kstp ) 
     229      !!---------------------------------------------------------------------- 
     230      !!                     ***  ROUTINE zdf_phy  *** 
     231      !! 
     232      !! ** Purpose :  Update ocean physics at each time-step 
     233      !! 
     234      !! ** Method  :  
     235      !! 
     236      !! ** Action  :   avm, avt vertical eddy viscosity and diffusivity at w-points 
     237      !!                nmld ??? mixed layer depth in level and meters   <<<<====verifier ! 
     238      !!                bottom stress.....                               <<<<====verifier ! 
     239      !!---------------------------------------------------------------------- 
     240      INTEGER, INTENT(in) ::   kstp   ! ocean time-step index 
     241      ! 
     242      INTEGER ::   ji, jj, jk    ! dummy loop indice 
     243      !!---------------------------------------------------------------------- 
     244      ! 
     245                         CALL zdf_bfr( kstp )         ! bottom friction (if quadratic) 
     246      !                                               ! Vertical eddy viscosity and diffusivity coefficients 
     247      IF( ln_zdfric  )   CALL zdf_ric ( kstp )             ! Richardson number dependent Kz 
     248      IF( ln_zdftke  )   CALL zdf_tke ( kstp )             ! TKE closure scheme for Kz 
     249      IF( ln_zdfgls  )   CALL zdf_gls ( kstp )             ! GLS closure scheme for Kz 
     250      IF( ln_zdfqiao )   CALL zdf_qiao( kstp )             ! Qiao vertical mixing  
     251      ! 
     252      IF( ln_zdfcst  ) THEN                                ! Constant Kz (reset avt, avm[uv] to the background value) 
     253         avt (:,:,:) = rn_avt0 * wmask (:,:,:) 
     254         avm (:,:,:) = rn_avm0 * wmask (:,:,:) 
     255         avmu(:,:,:) = rn_avm0 * wumask(:,:,:) 
     256         avmv(:,:,:) = rn_avm0 * wvmask(:,:,:) 
     257      ENDIF 
     258      ! 
     259      IF( ln_rnf_mouth ) THEN                         ! increase diffusivity at rivers mouths 
     260         DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2._wp * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO 
     261      ENDIF 
     262      ! 
     263      IF( ln_zdfevd  )   CALL zdf_evd( kstp )         ! enhanced vertical eddy diffusivity 
     264      ! 
     265      IF( ln_zdfddm  ) THEN                           ! double diffusive mixing 
     266                         CALL zdf_ddm( kstp ) 
     267      ELSE                                            ! avs=avt 
     268         DO jk = 2, jpkm1   ;   avs(:,:,jk) = avt(:,:,jk)   ;   END DO 
     269      ENDIF 
     270      ! 
     271      IF( ln_zdftmx  )   CALL zdf_tmx( kstp )         ! tidal vertical mixing 
     272 
     273                         CALL zdf_mxl( kstp )         ! mixed layer depth 
     274 
     275                                                      ! write TKE or GLS information in the restart file 
     276      IF( lrst_oce .AND. ln_zdftke )   CALL tke_rst( kstp, 'WRITE' ) 
     277      IF( lrst_oce .AND. ln_zdfgls )   CALL gls_rst( kstp, 'WRITE' ) 
     278      ! 
     279   END SUBROUTINE zdf_phy 
    205280 
    206281   !!====================================================================== 
    207 END MODULE zdfini 
     282END MODULE zdfphy 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfqiao.F90

    r7646 r7953  
    6767      !--------------------------------------------------------------------------------- 
    6868      ! 
     69!!gm Comment: I don't understand the use of min of 4 gdepw_n to define a quantity at w-point 
     70!!gm                       ==>> this is an error.... 
    6971      DO jk = 1, jpk 
    7072         DO jj = 1, jpjm1 
     
    101103      !------------------------------- 
    102104      ! 
     105!!gm  with double diffusion activated, avs is not updated... 
     106!!gm                    =====>>> BUG 
    103107      DO jk = 1, jpkm1 
    104108         DO jj = 1, jpj 
     
    112116      ! 
    113117   END SUBROUTINE zdf_qiao 
     118 
    114119 
    115120   INTEGER FUNCTION zdf_qiao_alloc() 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfric.F90

    r7646 r7953  
    1313   !!            3.3.1! 2011-09  (P. Oddo) Mixed layer depth parameterization 
    1414   !!---------------------------------------------------------------------- 
    15 #if defined key_zdfric 
    16    !!---------------------------------------------------------------------- 
    17    !!   'key_zdfric'                                             Kz = f(Ri) 
    18    !!---------------------------------------------------------------------- 
    19    !!   zdf_ric       : update momentum and tracer Kz from the Richardson 
    20    !!                  number computation 
     15 
     16   !!---------------------------------------------------------------------- 
     17   !!   zdf_ric       : update momentum and tracer Kz from the Richardson number 
    2118   !!   zdf_ric_init  : initialization, namelist read, & parameters control 
    2219   !!---------------------------------------------------------------------- 
     
    3835   PUBLIC   zdf_ric         ! called by step.F90 
    3936   PUBLIC   zdf_ric_init    ! called by opa.F90 
    40  
    41    LOGICAL, PUBLIC, PARAMETER ::   lk_zdfric = .TRUE.   !: Richardson vertical mixing flag 
    4237 
    4338   !                           !!* Namelist namzdf_ric : Richardson number dependent Kz * 
     
    108103      !!      namelist 
    109104      !!        N.B. the mask are required for implicit scheme, and surface 
    110       !!      and bottom value already set in zdfini.F90 
     105      !!      and bottom value already set in zdfphy.F90 
    111106      !! 
    112107      !! References : Pacanowski & Philander 1981, JPO, 1441-1451. 
     
    183178            zrhos          = rhop(ji,jj,1) + zflageos * ( 1. - tmask(ji,jj,1) ) 
    184179            zustar         = SQRT( taum(ji,jj) / ( zrhos +  rsmall ) ) 
    185             ekm_dep(ji,jj) = rn_ekmfc * zustar / ( ABS( ff(ji,jj) ) + rsmall ) 
     180            ekm_dep(ji,jj) = rn_ekmfc * zustar / ( ABS( ff_t(ji,jj) ) + rsmall ) 
    186181            ekm_dep(ji,jj) = MAX(ekm_dep(ji,jj),rn_mldmin) ! Minimun allowed 
    187182            ekm_dep(ji,jj) = MIN(ekm_dep(ji,jj),rn_mldmax) ! Maximum allowed 
     
    303298   END SUBROUTINE zdf_ric_init 
    304299 
    305 #else 
    306    !!---------------------------------------------------------------------- 
    307    !!   Dummy module :              NO Richardson dependent vertical mixing 
    308    !!---------------------------------------------------------------------- 
    309    LOGICAL, PUBLIC, PARAMETER ::   lk_zdfric = .FALSE.   !: Richardson mixing flag 
    310 CONTAINS 
    311    SUBROUTINE zdf_ric_init         ! Dummy routine 
    312    END SUBROUTINE zdf_ric_init 
    313    SUBROUTINE zdf_ric( kt )        ! Dummy routine 
    314       WRITE(*,*) 'zdf_ric: You should not have seen this print! error?', kt 
    315    END SUBROUTINE zdf_ric 
    316 #endif 
    317  
    318300   !!====================================================================== 
    319301END MODULE zdfric 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r7813 r7953  
    2727   !!            3.3  !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase 
    2828   !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability 
     29   !!            4.0  !  2017-04  (G. Madec)  Remove CPP keys 
    2930   !!---------------------------------------------------------------------- 
    30 #if defined key_zdftke 
    31    !!---------------------------------------------------------------------- 
    32    !!   'key_zdftke'                                   TKE vertical physics 
     31 
    3332   !!---------------------------------------------------------------------- 
    3433   !!   zdf_tke       : update momentum and tracer Kz from a tke scheme 
     
    6564   PUBLIC   tke_rst        ! routine called in step module 
    6665 
    67    LOGICAL , PUBLIC, PARAMETER ::   lk_zdftke = .TRUE.  !: TKE vertical mixing flag 
    68  
    6966   !                      !!** Namelist  namzdf_tke  ** 
    7067   LOGICAL  ::   ln_mxl0   ! mixing length scale surface value as function of wind stress or not 
     
    376373            DO ji = fs_2, fs_jpim1   ! vector opt. 
    377374               zcof   = zfact1 * tmask(ji,jj,jk) 
    378 # if defined key_zdftmx_new 
    379                ! key_zdftmx_new: New internal wave-driven param: set a minimum value for Kz on TKE (ensure numerical stability) 
    380                zzd_up = zcof * MAX( avm(ji,jj,jk+1) + avm(ji,jj,jk), 2.e-5_wp )   &  ! upper diagonal 
    381                   &          / (  e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk  )  ) 
    382                zzd_lw = zcof * MAX( avm(ji,jj,jk) + avm(ji,jj,jk-1), 2.e-5_wp )   &  ! lower diagonal 
    383                   &          / (  e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk  )  ) 
    384 # else 
    385                zzd_up = zcof * ( avm  (ji,jj,jk+1) + avm  (ji,jj,jk  ) )   &  ! upper diagonal 
    386                   &          / ( e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk  ) ) 
    387                zzd_lw = zcof * ( avm  (ji,jj,jk  ) + avm  (ji,jj,jk-1) )   &  ! lower diagonal 
    388                   &          / ( e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk  ) ) 
    389 # endif 
     375               !                                   ! A minimum of 2.e-5 m2/s is imposed on TKE vertical 
     376               !                                   ! eddy coefficient (ensure numerical stability) 
     377               zzd_up = zcof * MAX(   avm(ji,jj,jk+1) +   avm(ji,jj,jk  ) , 2.e-5_wp  )   &  ! upper diagonal 
     378                  &          /    ( e3t_n(ji,jj,jk  ) * e3w_n(ji,jj,jk  )  ) 
     379               zzd_lw = zcof * MAX(   avm(ji,jj,jk  ) +   avm(ji,jj,jk-1) , 2.e-5_wp  )   &  ! lower diagonal 
     380                  &          /    ( e3t_n(ji,jj,jk-1) * e3w_n(ji,jj,jk  )  ) 
     381               ! 
    390382               !                                   ! shear prod. at w-point weightened by mask 
    391383               zesh2  =  ( z3du(ji-1,jj,jk) + z3du(ji,jj,jk) ) / MAX( 1._wp , umask(ji-1,jj,jk) + umask(ji,jj,jk) )   & 
     
    741733      ! 
    742734      ri_cri   = 2._wp    / ( 2._wp + rn_ediss / rn_ediff )   ! resulting critical Richardson number 
    743 # if defined key_zdftmx_new 
    744       ! key_zdftmx_new: New internal wave-driven param: specified value of rn_emin & rmxl_min are used 
    745       rn_emin  = 1.e-10_wp 
    746       rmxl_min = 1.e-03_wp 
    747       IF(lwp) THEN                  ! Control print 
    748          WRITE(numout,*) 
    749          WRITE(numout,*) 'zdf_tke_init :  New tidal mixing case: force rn_emin = 1.e-10 and rmxl_min = 1.e-3 ' 
    750          WRITE(numout,*) '~~~~~~~~~~~~' 
    751       ENDIF 
    752 # else 
    753       rmxl_min = 1.e-6_wp / ( rn_ediff * SQRT( rn_emin ) )    ! resulting minimum length to recover molecular viscosity 
    754 # endif 
    755735      ! 
    756736      IF(lwp) THEN                    !* Control print 
     
    776756         WRITE(numout,*) 
    777757         WRITE(numout,*) '      critical Richardson nb with your parameters  ri_cri = ', ri_cri 
     758         WRITE(numout,*) 
     759      ENDIF 
     760      ! 
     761      IF( ln_zdftmx ) THEN          ! Internal wave driven mixing 
     762         !                          ! specific values of rn_emin & rmxl_min are used 
     763         rn_emin  = 1.e-10_wp 
     764         rmxl_min = 1.e-03_wp 
     765         IF(lwp) WRITE(numout,*) '      Internal wave-driven mixing case:   force   rn_emin = 1.e-10 and rmxl_min = 1.e-3 ' 
     766      ELSE 
     767         rmxl_min = 1.e-6_wp / ( rn_ediff * SQRT( rn_emin ) )    ! resulting minimum length to recover molecular viscosity 
     768         IF(lwp) WRITE(numout,*) '      minimum mixing length with your parameters rmxl_min = ', rmxl_min 
    778769      ENDIF 
    779770      ! 
     
    891882   END SUBROUTINE tke_rst 
    892883 
    893 #else 
    894    !!---------------------------------------------------------------------- 
    895    !!   Dummy module :                                        NO TKE scheme 
    896    !!---------------------------------------------------------------------- 
    897    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftke = .FALSE.   !: TKE flag 
    898 CONTAINS 
    899    SUBROUTINE zdf_tke_init           ! Dummy routine 
    900    END SUBROUTINE zdf_tke_init 
    901    SUBROUTINE zdf_tke( kt )          ! Dummy routine 
    902       WRITE(*,*) 'zdf_tke: You should not have seen this print! error?', kt 
    903    END SUBROUTINE zdf_tke 
    904    SUBROUTINE tke_rst( kt, cdrw ) 
    905      CHARACTER(len=*) ::   cdrw 
    906      WRITE(*,*) 'tke_rst: You should not have seen this print! error?', kt, cdwr 
    907    END SUBROUTINE tke_rst 
    908 #endif 
    909  
    910884   !!====================================================================== 
    911885END MODULE zdftke 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90

    r7931 r7953  
    22   !!======================================================================== 
    33   !!                       ***  MODULE  zdftmx  *** 
    4    !! Ocean physics: vertical tidal mixing coefficient 
     4   !! Ocean physics: Internal gravity wave-driven vertical mixing 
    55   !!======================================================================== 
    66   !! History :  1.0  !  2004-04  (L. Bessieres, G. Madec)  Original code 
    7    !!             -   !  2006-08  (A. Koch-Larrouy) Indonesian strait 
    8    !!            3.3  !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase 
     7   !!             -   !  2006-08  (A. Koch-Larrouy)  Indonesian strait 
     8   !!            3.3  !  2010-10  (C. Ethe, G. Madec)  reorganisation of initialisation phase 
     9   !!            3.6  !  2016-03  (C. de Lavergne)  New param: internal wave-driven mixing  
     10   !!            4.0  !  2017-04  (G. Madec)  Remove the old tidal mixing param. and key zdftmx(_new) 
    911   !!---------------------------------------------------------------------- 
    10 #if defined key_zdftmx 
    11    !!---------------------------------------------------------------------- 
    12    !!   'key_zdftmx'                                  Tidal vertical mixing 
    13    !!---------------------------------------------------------------------- 
    14    !!   zdf_tmx       : global     momentum & tracer Kz with tidal induced Kz 
    15    !!   tmx_itf       : Indonesian momentum & tracer Kz with tidal induced Kz  
    16    !!---------------------------------------------------------------------- 
    17    USE oce            ! ocean dynamics and tracers variables 
    18    USE dom_oce        ! ocean space and time domain variables 
    19    USE zdf_oce        ! ocean vertical physics variables 
    20    USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    21    USE eosbn2         ! ocean equation of state 
    22    USE phycst         ! physical constants 
    23    USE prtctl         ! Print control 
    24    USE in_out_manager ! I/O manager 
    25    USE iom            ! I/O Manager 
    26    USE lib_mpp        ! MPP library 
    27    USE wrk_nemo       ! work arrays 
    28    USE timing         ! Timing 
    29    USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
    30  
    31    IMPLICIT NONE 
    32    PRIVATE 
    33  
    34    PUBLIC   zdf_tmx         ! called in step module  
    35    PUBLIC   zdf_tmx_init    ! called in opa module  
    36    PUBLIC   zdf_tmx_alloc   ! called in nemogcm module 
    37  
    38    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftmx = .TRUE.    !: tidal mixing flag 
    39  
    40    !                       !!* Namelist  namzdf_tmx : tidal mixing * 
    41    REAL(wp) ::  rn_htmx     ! vertical decay scale for turbulence (meters) 
    42    REAL(wp) ::  rn_n2min    ! threshold of the Brunt-Vaisala frequency (s-1) 
    43    REAL(wp) ::  rn_tfe      ! tidal dissipation efficiency (St Laurent et al. 2002) 
    44    REAL(wp) ::  rn_me       ! mixing efficiency (Osborn 1980) 
    45    LOGICAL  ::  ln_tmx_itf  ! Indonesian Through Flow (ITF): Koch-Larrouy et al. (2007) parameterization 
    46    REAL(wp) ::  rn_tfe_itf  ! ITF tidal dissipation efficiency (St Laurent et al. 2002) 
    47  
    48    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   en_tmx     ! energy available for tidal mixing (W/m2) 
    49    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   mask_itf   ! mask to use over Indonesian area 
    50    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   az_tmx     ! coefficient used to evaluate the tidal induced Kz 
    51  
    52    !! * Substitutions 
    53 #  include "vectopt_loop_substitute.h90" 
    54    !!---------------------------------------------------------------------- 
    55    !! NEMO/OPA 3.7 , NEMO Consortium (2014) 
    56    !! $Id$ 
    57    !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    58    !!---------------------------------------------------------------------- 
    59 CONTAINS 
    60  
    61    INTEGER FUNCTION zdf_tmx_alloc() 
    62       !!---------------------------------------------------------------------- 
    63       !!                ***  FUNCTION zdf_tmx_alloc  *** 
    64       !!---------------------------------------------------------------------- 
    65       ALLOCATE(en_tmx(jpi,jpj), mask_itf(jpi,jpj), az_tmx(jpi,jpj,jpk), STAT=zdf_tmx_alloc ) 
    66       ! 
    67       IF( lk_mpp             )   CALL mpp_sum ( zdf_tmx_alloc ) 
    68       IF( zdf_tmx_alloc /= 0 )   CALL ctl_warn('zdf_tmx_alloc: failed to allocate arrays') 
    69    END FUNCTION zdf_tmx_alloc 
    70  
    71  
    72    SUBROUTINE zdf_tmx( kt ) 
    73       !!---------------------------------------------------------------------- 
    74       !!                  ***  ROUTINE zdf_tmx  *** 
    75       !!                    
    76       !! ** Purpose :   add to the vertical mixing coefficients the effect of 
    77       !!              tidal mixing (Simmons et al 2004). 
    78       !! 
    79       !! ** Method  : - tidal-induced vertical mixing is given by: 
    80       !!                  Kz_tides = az_tmx / max( rn_n2min, N^2 ) 
    81       !!              where az_tmx is a coefficient that specified the 3D space  
    82       !!              distribution of the faction of tidal energy taht is used 
    83       !!              for mixing. Its expression is set in zdf_tmx_init routine, 
    84       !!              following Simmons et al. 2004. 
    85       !!                NB: a specific bounding procedure is performed on av_tide 
    86       !!              so that the input tidal energy is actually almost used. The 
    87       !!              basic maximum value is 60 cm2/s, but values of 300 cm2/s  
    88       !!              can be reached in area where bottom stratification is too  
    89       !!              weak. 
    90       !! 
    91       !!              - update av_tide in the Indonesian Through Flow area 
    92       !!              following Koch-Larrouy et al. (2007) parameterisation 
    93       !!              (see tmx_itf routine). 
    94       !! 
    95       !!              - update the model vertical eddy viscosity and diffusivity:  
    96       !!                     avt  = avt  +    av_tides 
    97       !!                     avm  = avm  +    av_tides 
    98       !!                     avmu = avmu + mi(av_tides) 
    99       !!                     avmv = avmv + mj(av_tides) 
    100       !! 
    101       !! ** Action  :   avt, avm, avmu, avmv   increased by tidal mixing 
    102       !! 
    103       !! References : Simmons et al. 2004, Ocean Modelling, 6, 3-4, 245-263. 
    104       !!              Koch-Larrouy et al. 2007, GRL. 
    105       !!---------------------------------------------------------------------- 
    106       INTEGER, INTENT(in) ::   kt   ! ocean time-step  
    107       ! 
    108       INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    109       REAL(wp) ::   ztpc         ! scalar workspace 
    110       REAL(wp), POINTER, DIMENSION(:,:)   ::   zkz 
    111       REAL(wp), POINTER, DIMENSION(:,:,:) ::   zav_tide 
    112       !!---------------------------------------------------------------------- 
    113       ! 
    114       IF( nn_timing == 1 )  CALL timing_start('zdf_tmx') 
    115       ! 
    116       CALL wrk_alloc( jpi,jpj,       zkz ) 
    117       CALL wrk_alloc( jpi,jpj,jpk,   zav_tide ) 
    118       ! 
    119       !                          ! ----------------------- ! 
    120       !                          !  Standard tidal mixing  !  (compute zav_tide) 
    121       !                          ! ----------------------- ! 
    122       !                             !* First estimation (with n2 bound by rn_n2min) bounded by 60 cm2/s 
    123       zav_tide(:,:,:) = MIN(  60.e-4, az_tmx(:,:,:) / MAX( rn_n2min, rn2(:,:,:) )  ) 
    124  
    125       zkz(:,:) = 0.e0               !* Associated potential energy consummed over the whole water column 
    126       DO jk = 2, jpkm1 
    127          zkz(:,:) = zkz(:,:) + e3w_n(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zav_tide(:,:,jk) * wmask(:,:,jk) 
    128       END DO 
    129  
    130       DO jj = 1, jpj                !* Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx 
    131          DO ji = 1, jpi 
    132             IF( zkz(ji,jj) /= 0.e0 )   zkz(ji,jj) = en_tmx(ji,jj) / zkz(ji,jj) 
    133          END DO 
    134       END DO 
    135  
    136       DO jk = 2, jpkm1     !* Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zav_tide bound by 300 cm2/s 
    137          zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) * wmask(:,:,jk)  !kz max = 300 cm2/s 
    138       END DO 
    139  
    140       IF( kt == nit000 ) THEN       !* check at first time-step: diagnose the energy consumed by zav_tide 
    141          ztpc = 0._wp 
    142          DO jk= 1, jpk 
    143             DO jj= 1, jpj 
    144                DO ji= 1, jpi 
    145                   ztpc = ztpc + e3w_n(ji,jj,jk) * e1e2t(ji,jj)                  & 
    146                      &        * MAX( 0.e0, rn2(ji,jj,jk) ) * zav_tide(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj) 
    147                END DO 
    148             END DO 
    149          END DO 
    150          ztpc= rau0 / ( rn_tfe * rn_me ) * ztpc 
    151          IF( lk_mpp )   CALL mpp_sum( ztpc ) 
    152          IF(lwp) WRITE(numout,*)  
    153          IF(lwp) WRITE(numout,*) '          N Total power consumption by av_tide    : ztpc = ', ztpc * 1.e-12 ,'TW' 
    154       ENDIF 
    155         
    156       !                          ! ----------------------- ! 
    157       !                          !    ITF  tidal mixing    !  (update zav_tide) 
    158       !                          ! ----------------------- ! 
    159       IF( ln_tmx_itf )   CALL tmx_itf( kt, zav_tide ) 
    160  
    161       !                          ! ----------------------- ! 
    162       !                          !   Update  mixing coefs  !                           
    163       !                          ! ----------------------- ! 
    164       DO jk = 2, jpkm1              !* update momentum & tracer diffusivity with tidal mixing 
    165          avt(:,:,jk) = avs(:,:,jk) + zav_tide(:,:,jk) * wmask(:,:,jk) 
    166          avm(:,:,jk) = avm(:,:,jk) + zav_tide(:,:,jk) * wmask(:,:,jk) 
    167          DO jj = 2, jpjm1 
    168             DO ji = fs_2, fs_jpim1  ! vector opt. 
    169                avmu(ji,jj,jk) = avmu(ji,jj,jk) + 0.5 * ( zav_tide(ji,jj,jk) + zav_tide(ji+1,jj  ,jk) ) * wumask(ji,jj,jk) 
    170                avmv(ji,jj,jk) = avmv(ji,jj,jk) + 0.5 * ( zav_tide(ji,jj,jk) + zav_tide(ji  ,jj+1,jk) ) * wvmask(ji,jj,jk) 
    171             END DO 
    172          END DO 
    173       END DO 
    174       CALL lbc_lnk( avmu, 'U', 1. )   ;   CALL lbc_lnk( avmv, 'V', 1. )      ! lateral boundary condition 
    175  
    176       !                             !* output tidal mixing coefficient 
    177       CALL iom_put( "av_tide", zav_tide ) 
    178  
    179       IF(ln_ctl)   CALL prt_ctl(tab3d_1=zav_tide , clinfo1=' tmx - av_tide: ', tab3d_2=avt, clinfo2=' avt: ', ovlap=1, kdim=jpk) 
    180       ! 
    181       CALL wrk_dealloc( jpi,jpj,       zkz ) 
    182       CALL wrk_dealloc( jpi,jpj,jpk,   zav_tide ) 
    183       ! 
    184       IF( nn_timing == 1 )  CALL timing_stop('zdf_tmx') 
    185       ! 
    186    END SUBROUTINE zdf_tmx 
    187  
    188  
    189    SUBROUTINE tmx_itf( kt, pav ) 
    190       !!---------------------------------------------------------------------- 
    191       !!                  ***  ROUTINE tmx_itf  *** 
    192       !!                    
    193       !! ** Purpose :   modify the vertical eddy diffusivity coefficients  
    194       !!              (pav) in the Indonesian Through Flow area (ITF). 
    195       !! 
    196       !! ** Method  : - Following Koch-Larrouy et al. (2007), in the ITF defined 
    197       !!                by msk_itf (read in a file, see tmx_init), the tidal 
    198       !!                mixing coefficient is computed with : 
    199       !!                  * q=1 (i.e. all the tidal energy remains trapped in 
    200       !!                         the area and thus is used for mixing) 
    201       !!                  * the vertical distribution of the tifal energy is a 
    202       !!                    proportional to N above the thermocline (d(N^2)/dz > 0) 
    203       !!                    and to N^2 below the thermocline (d(N^2)/dz < 0) 
    204       !! 
    205       !! ** Action  :   av_tide   updated in the ITF area (msk_itf) 
    206       !! 
    207       !! References :  Koch-Larrouy et al. 2007, GRL  
    208       !!---------------------------------------------------------------------- 
    209       INTEGER , INTENT(in   )                         ::   kt   ! ocean time-step 
    210       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pav  ! Tidal mixing coef. 
    211       !!  
    212       INTEGER  ::   ji, jj, jk    ! dummy loop indices 
    213       REAL(wp) ::   zcoef, ztpc   ! temporary scalar 
    214       REAL(wp), DIMENSION(:,:)  , POINTER ::   zkz                        ! 2D workspace 
    215       REAL(wp), DIMENSION(:,:)  , POINTER ::   zsum1 , zsum2 , zsum       !  -      - 
    216       REAL(wp), DIMENSION(:,:,:), POINTER ::   zempba_3d_1, zempba_3d_2   ! 3D workspace 
    217       REAL(wp), DIMENSION(:,:,:), POINTER ::   zempba_3d  , zdn2dz        !  -      - 
    218       REAL(wp), DIMENSION(:,:,:), POINTER ::   zavt_itf                   !  -      - 
    219       !!---------------------------------------------------------------------- 
    220       ! 
    221       IF( nn_timing == 1 )  CALL timing_start('tmx_itf') 
    222       ! 
    223       CALL wrk_alloc( jpi,jpj, zkz, zsum1 , zsum2 , zsum ) 
    224       CALL wrk_alloc( jpi,jpj,jpk, zempba_3d_1, zempba_3d_2, zempba_3d, zdn2dz, zavt_itf ) 
    225  
    226       !                             ! compute the form function using N2 at each time step 
    227       zempba_3d_1(:,:,jpk) = 0.e0 
    228       zempba_3d_2(:,:,jpk) = 0.e0 
    229       DO jk = 1, jpkm1              
    230          zdn2dz     (:,:,jk) = rn2(:,:,jk) - rn2(:,:,jk+1)           ! Vertical profile of dN2/dz 
    231          zempba_3d_1(:,:,jk) = SQRT(  MAX( 0.e0, rn2(:,:,jk) )  )    !    -        -    of N 
    232          zempba_3d_2(:,:,jk) =        MAX( 0.e0, rn2(:,:,jk) )       !    -        -    of N^2 
    233       END DO 
    234       ! 
    235       zsum (:,:) = 0.e0 
    236       zsum1(:,:) = 0.e0 
    237       zsum2(:,:) = 0.e0 
    238       DO jk= 2, jpk 
    239          zsum1(:,:) = zsum1(:,:) + zempba_3d_1(:,:,jk) * e3w_n(:,:,jk) * wmask(:,:,jk) 
    240          zsum2(:,:) = zsum2(:,:) + zempba_3d_2(:,:,jk) * e3w_n(:,:,jk) * wmask(:,:,jk)                
    241       END DO 
    242       DO jj = 1, jpj 
    243          DO ji = 1, jpi 
    244             IF( zsum1(ji,jj) /= 0.e0 )   zsum1(ji,jj) = 1.e0 / zsum1(ji,jj) 
    245             IF( zsum2(ji,jj) /= 0.e0 )   zsum2(ji,jj) = 1.e0 / zsum2(ji,jj)                 
    246          END DO 
    247       END DO 
    248  
    249       DO jk= 1, jpk 
    250          DO jj = 1, jpj 
    251             DO ji = 1, jpi 
    252                zcoef = 0.5 - SIGN( 0.5, zdn2dz(ji,jj,jk) )       ! =0 if dN2/dz > 0, =1 otherwise  
    253                ztpc  = zempba_3d_1(ji,jj,jk) * zsum1(ji,jj) *        zcoef     & 
    254                   &  + zempba_3d_2(ji,jj,jk) * zsum2(ji,jj) * ( 1. - zcoef ) 
    255                ! 
    256                zempba_3d(ji,jj,jk) =               ztpc  
    257                zsum     (ji,jj)    = zsum(ji,jj) + ztpc * e3w_n(ji,jj,jk) 
    258             END DO 
    259          END DO 
    260        END DO 
    261        DO jj = 1, jpj 
    262           DO ji = 1, jpi 
    263              IF( zsum(ji,jj) > 0.e0 )   zsum(ji,jj) = 1.e0 / zsum(ji,jj)                 
    264           END DO 
    265        END DO 
    266  
    267       !                             ! first estimation bounded by 10 cm2/s (with n2 bounded by rn_n2min)  
    268       zcoef = rn_tfe_itf / ( rn_tfe * rau0 ) 
    269       DO jk = 1, jpk 
    270          zavt_itf(:,:,jk) = MIN(  10.e-4, zcoef * en_tmx(:,:) * zsum(:,:) * zempba_3d(:,:,jk)   & 
    271             &                                      / MAX( rn_n2min, rn2(:,:,jk) ) * tmask(:,:,jk)  ) 
    272       END DO            
    273  
    274       zkz(:,:) = 0.e0               ! Associated potential energy consummed over the whole water column 
    275       DO jk = 2, jpkm1 
    276          zkz(:,:) = zkz(:,:) + e3w_n(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zavt_itf(:,:,jk) * wmask(:,:,jk) 
    277       END DO 
    278  
    279       DO jj = 1, jpj                ! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx 
    280          DO ji = 1, jpi 
    281             IF( zkz(ji,jj) /= 0.e0 )   zkz(ji,jj) = en_tmx(ji,jj) * rn_tfe_itf / rn_tfe / zkz(ji,jj) 
    282          END DO 
    283       END DO 
    284  
    285       DO jk = 2, jpkm1              ! Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zavt_itf bound by 300 cm2/s 
    286          zavt_itf(:,:,jk) = zavt_itf(:,:,jk) * MIN( zkz(:,:), 120./10. ) * wmask(:,:,jk)   ! kz max = 120 cm2/s 
    287       END DO 
    288  
    289       IF( kt == nit000 ) THEN       ! diagnose the nergy consumed by zavt_itf 
    290          ztpc = 0.e0 
    291          DO jk= 1, jpk 
    292             DO jj= 1, jpj 
    293                DO ji= 1, jpi 
    294                   ztpc = ztpc + e1e2t(ji,jj) * e3w_n(ji,jj,jk) * MAX( 0.e0, rn2(ji,jj,jk) )   & 
    295                      &                       * zavt_itf(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj) 
    296                END DO 
    297             END DO 
    298          END DO 
    299          IF( lk_mpp )   CALL mpp_sum( ztpc ) 
    300          ztpc= rau0 * ztpc / ( rn_me * rn_tfe_itf ) 
    301          IF(lwp) WRITE(numout,*) '          N Total power consumption by zavt_itf: ztpc = ', ztpc * 1.e-12 ,'TW' 
    302       ENDIF 
    303  
    304       !                             ! Update pav with the ITF mixing coefficient 
    305       DO jk = 2, jpkm1 
    306          pav(:,:,jk) = pav     (:,:,jk) * ( 1.e0 - mask_itf(:,:) )   & 
    307             &        + zavt_itf(:,:,jk) *          mask_itf(:,:)  
    308       END DO 
    309       ! 
    310       CALL wrk_dealloc( jpi,jpj, zkz, zsum1 , zsum2 , zsum ) 
    311       CALL wrk_dealloc( jpi,jpj,jpk, zempba_3d_1, zempba_3d_2, zempba_3d, zdn2dz, zavt_itf ) 
    312       ! 
    313       IF( nn_timing == 1 )  CALL timing_stop('tmx_itf') 
    314       ! 
    315    END SUBROUTINE tmx_itf 
    316  
    317  
    318    SUBROUTINE zdf_tmx_init 
    319       !!---------------------------------------------------------------------- 
    320       !!                  ***  ROUTINE zdf_tmx_init  *** 
    321       !!                      
    322       !! ** Purpose :   Initialization of the vertical tidal mixing, Reading 
    323       !!              of M2 and K1 tidal energy in nc files 
    324       !! 
    325       !! ** Method  : - Read the namtmx namelist and check the parameters 
    326       !! 
    327       !!              - Read the input data in NetCDF files : 
    328       !!              M2 and K1 tidal energy. The total tidal energy, en_tmx,  
    329       !!              is the sum of M2, K1 and S2 energy where S2 is assumed  
    330       !!              to be: S2=(1/2)^2 * M2 
    331       !!              mask_itf, a mask array that determine where substituing  
    332       !!              the standard Simmons et al. (2005) formulation with the 
    333       !!              one of Koch_Larrouy et al. (2007). 
    334       !! 
    335       !!              - Compute az_tmx, a 3D coefficient that allows to compute 
    336       !!             the standard tidal-induced vertical mixing as follows: 
    337       !!                  Kz_tides = az_tmx / max( rn_n2min, N^2 ) 
    338       !!             with az_tmx a bottom intensified coefficient is given by: 
    339       !!                 az_tmx(z) = en_tmx / ( rau0 * rn_htmx ) * EXP( -(H-z)/rn_htmx ) 
    340       !!                                                  / ( 1. - EXP( - H   /rn_htmx ) )  
    341       !!             where rn_htmx the characteristic length scale of the bottom  
    342       !!             intensification, en_tmx the tidal energy, and H the ocean depth 
    343       !! 
    344       !! ** input   :   - Namlist namtmx 
    345       !!                - NetCDF file : M2_ORCA2.nc, K1_ORCA2.nc, and mask_itf.nc 
    346       !! 
    347       !! ** Action  : - Increase by 1 the nstop flag is setting problem encounter 
    348       !!              - defined az_tmx used to compute tidal-induced mixing 
    349       !! 
    350       !! References : Simmons et al. 2004, Ocean Modelling, 6, 3-4, 245-263. 
    351       !!              Koch-Larrouy et al. 2007, GRL. 
    352       !!---------------------------------------------------------------------- 
    353       INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    354       INTEGER  ::   inum         ! local integer 
    355       INTEGER  ::   ios 
    356       REAL(wp) ::   ztpc, ze_z   ! local scalars 
    357       REAL(wp), DIMENSION(:,:)  , POINTER ::  zem2, zek1     ! read M2 and K1 tidal energy 
    358       REAL(wp), DIMENSION(:,:)  , POINTER ::  zkz            ! total M2, K1 and S2 tidal energy 
    359       REAL(wp), DIMENSION(:,:)  , POINTER ::  zfact          ! used for vertical structure function 
    360       REAL(wp), DIMENSION(:,:)  , POINTER ::  zhdep          ! Ocean depth  
    361       REAL(wp), DIMENSION(:,:,:), POINTER ::  zpc, zav_tide  ! power consumption 
    362       !! 
    363       NAMELIST/namzdf_tmx/ rn_htmx, rn_n2min, rn_tfe, rn_me, ln_tmx_itf, rn_tfe_itf 
    364       !!---------------------------------------------------------------------- 
    365       ! 
    366       IF( nn_timing == 1 )  CALL timing_start('zdf_tmx_init') 
    367       ! 
    368       CALL wrk_alloc( jpi,jpj,       zem2, zek1, zkz, zfact, zhdep ) 
    369       CALL wrk_alloc( jpi,jpj,jpk,   zpc, zav_tide ) 
    370       ! 
    371       REWIND( numnam_ref )             ! Namelist namzdf_tmx in reference namelist : Tidal Mixing 
    372       READ  ( numnam_ref, namzdf_tmx, IOSTAT = ios, ERR = 901) 
    373 901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_tmx in reference namelist', lwp ) 
    374       ! 
    375       REWIND( numnam_cfg )             ! Namelist namzdf_tmx in configuration namelist : Tidal Mixing 
    376       READ  ( numnam_cfg, namzdf_tmx, IOSTAT = ios, ERR = 902 ) 
    377 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_tmx in configuration namelist', lwp ) 
    378       IF(lwm) WRITE ( numond, namzdf_tmx ) 
    379       ! 
    380       IF(lwp) THEN                     ! Control print 
    381          WRITE(numout,*) 
    382          WRITE(numout,*) 'zdf_tmx_init : tidal mixing' 
    383          WRITE(numout,*) '~~~~~~~~~~~~' 
    384          WRITE(numout,*) '   Namelist namzdf_tmx : set tidal mixing parameters' 
    385          WRITE(numout,*) '      Vertical decay scale for turbulence   = ', rn_htmx  
    386          WRITE(numout,*) '      Brunt-Vaisala frequency threshold     = ', rn_n2min 
    387          WRITE(numout,*) '      Tidal dissipation efficiency          = ', rn_tfe 
    388          WRITE(numout,*) '      Mixing efficiency                     = ', rn_me 
    389          WRITE(numout,*) '      ITF specific parameterisation         = ', ln_tmx_itf 
    390          WRITE(numout,*) '      ITF tidal dissipation efficiency      = ', rn_tfe_itf 
    391       ENDIF 
    392       !                                ! allocate tmx arrays 
    393       IF( zdf_tmx_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'zdf_tmx_init : unable to allocate tmx arrays' ) 
    394  
    395       IF( ln_tmx_itf ) THEN            ! read the Indonesian Through Flow mask 
    396          CALL iom_open('mask_itf',inum) 
    397          CALL iom_get (inum, jpdom_data, 'tmaskitf',mask_itf,1) !  
    398          CALL iom_close(inum) 
    399       ENDIF 
    400       !                                ! read M2 tidal energy flux : W/m2  ( zem2 < 0 ) 
    401       CALL iom_open('M2rowdrg',inum) 
    402       CALL iom_get (inum, jpdom_data, 'field',zem2,1) !  
    403       CALL iom_close(inum) 
    404       !                                ! read K1 tidal energy flux : W/m2  ( zek1 < 0 ) 
    405       CALL iom_open('K1rowdrg',inum) 
    406       CALL iom_get (inum, jpdom_data, 'field',zek1,1) !  
    407       CALL iom_close(inum) 
    408       !                                ! Total tidal energy ( M2, S2 and K1  with S2=(1/2)^2 * M2 ) 
    409       !                                ! only the energy available for mixing is taken into account, 
    410       !                                ! (mixing efficiency tidal dissipation efficiency) 
    411       en_tmx(:,:) = - rn_tfe * rn_me * ( zem2(:,:) * 1.25 + zek1(:,:) ) * ssmask(:,:) 
    412  
    413 !============ 
    414 !TG: Bug for VVL? Should this section be moved out of _init and be updated at every timestep? 
    415 !!gm : you are right, but tidal mixing acts in deep ocean (H>500m) where e3 is O(100m) 
    416 !!     the error is thus ~1% which I feel comfortable with, compared to uncertainties in tidal energy dissipation. 
    417       !                                ! Vertical structure (az_tmx) 
    418       DO jj = 1, jpj                         ! part independent of the level 
    419          DO ji = 1, jpi 
    420             zhdep(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1)       ! depth of the ocean 
    421             zfact(ji,jj) = rau0 * rn_htmx * ( 1. - EXP( -zhdep(ji,jj) / rn_htmx ) ) 
    422             IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = en_tmx(ji,jj) / zfact(ji,jj) 
    423          END DO 
    424       END DO 
    425       DO jk= 1, jpk                          ! complete with the level-dependent part 
    426          DO jj = 1, jpj 
    427             DO ji = 1, jpi 
    428                az_tmx(ji,jj,jk) = zfact(ji,jj) * EXP( -( zhdep(ji,jj)-gdepw_0(ji,jj,jk) ) / rn_htmx ) * tmask(ji,jj,jk) 
    429             END DO 
    430          END DO 
    431       END DO 
    432 !=========== 
    433       ! 
    434       IF( nprint == 1 .AND. lwp ) THEN 
    435          ! Control print 
    436          ! Total power consumption due to vertical mixing 
    437          ! zpc = rau0 * 1/rn_me * rn2 * zav_tide 
    438          zav_tide(:,:,:) = 0.e0 
    439          DO jk = 2, jpkm1 
    440             zav_tide(:,:,jk) = az_tmx(:,:,jk) / MAX( rn_n2min, rn2(:,:,jk) ) 
    441          END DO 
    442          ! 
    443          ztpc = 0._wp 
    444          zpc(:,:,:) = MAX(rn_n2min,rn2(:,:,:)) * zav_tide(:,:,:) 
    445          DO jk= 2, jpkm1 
    446             DO jj = 1, jpj 
    447                DO ji = 1, jpi 
    448                   ztpc = ztpc + e3w_n(ji,jj,jk) * e1e2t(ji,jj) * zpc(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 
    449                END DO 
    450             END DO 
    451          END DO 
    452          IF( lk_mpp )   CALL mpp_sum( ztpc ) 
    453          ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc 
    454          ! 
    455          WRITE(numout,*)  
    456          WRITE(numout,*) '          Total power consumption of the tidally driven part of Kz : ztpc = ', ztpc * 1.e-12 ,'TW' 
    457          ! 
    458          ! control print 2 
    459          zav_tide(:,:,:) = MIN( zav_tide(:,:,:), 60.e-4 )    
    460          zkz(:,:) = 0._wp 
    461          DO jk = 2, jpkm1 
    462                zkz(:,:) = zkz(:,:) + e3w_n(:,:,jk) * MAX(0.e0, rn2(:,:,jk)) * rau0 * zav_tide(:,:,jk) * wmask(:,:,jk) 
    463          END DO 
    464          ! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz 
    465          DO jj = 1, jpj 
    466             DO ji = 1, jpi 
    467                IF( zkz(ji,jj) /= 0.e0 )   THEN 
    468                    zkz(ji,jj) = en_tmx(ji,jj) / zkz(ji,jj) 
    469                ENDIF 
    470             END DO 
    471          END DO 
    472          ztpc = 1.e50 
    473          DO jj = 1, jpj 
    474             DO ji = 1, jpi 
    475                IF( zkz(ji,jj) /= 0.e0 )   THEN 
    476                    ztpc = Min( zkz(ji,jj), ztpc) 
    477                ENDIF 
    478             END DO 
    479          END DO 
    480          WRITE(numout,*) '          Min de zkz ', ztpc, ' Max = ', maxval(zkz(:,:) ) 
    481          ! 
    482          DO jk = 2, jpkm1 
    483             zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) * wmask(:,:,jk)  !kz max = 300 cm2/s 
    484          END DO 
    485          ztpc = 0._wp 
    486          zpc(:,:,:) = Max(0.e0,rn2(:,:,:)) * zav_tide(:,:,:) 
    487          DO jk= 1, jpk 
    488             DO jj = 1, jpj 
    489                DO ji = 1, jpi 
    490                   ztpc = ztpc + e3w_n(ji,jj,jk) * e1e2t(ji,jj) * zpc(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 
    491                END DO 
    492             END DO 
    493          END DO 
    494          IF( lk_mpp )   CALL mpp_sum( ztpc ) 
    495          ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc 
    496          WRITE(numout,*) '          2 Total power consumption of the tidally driven part of Kz : ztpc = ', ztpc * 1.e-12 ,'TW' 
    497 !!gm bug mpp  in these diagnostics 
    498          DO jk = 1, jpk 
    499             ze_z =                  SUM( e1e2t(:,:) * zav_tide(:,:,jk) * tmask_i(:,:) )   & 
    500                &     / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask   (:,:,jk) * tmask_i(:,:) ) ) 
    501             ztpc = 1.e50 
    502             DO jj = 1, jpj 
    503                DO ji = 1, jpi 
    504                   IF( zav_tide(ji,jj,jk) /= 0.e0 )   ztpc = MIN( ztpc, zav_tide(ji,jj,jk) ) 
    505                END DO 
    506             END DO 
    507             WRITE(numout,*) '            N2 min - jk= ', jk,'   ', ze_z * 1.e4,' cm2/s min= ',ztpc*1.e4,   & 
    508                &       'max= ', MAXVAL(zav_tide(:,:,jk) )*1.e4, ' cm2/s' 
    509          END DO 
    510  
    511          WRITE(numout,*) '          e_tide : ', SUM( e1e2t*en_tmx ) / ( rn_tfe * rn_me ) * 1.e-12, 'TW' 
    512          WRITE(numout,*)  
    513          WRITE(numout,*) '          Initial profile of tidal vertical mixing' 
    514          DO jk = 1, jpk 
    515             DO jj = 1,jpj 
    516                DO ji = 1,jpi 
    517                   zkz(ji,jj) = az_tmx(ji,jj,jk) /MAX( rn_n2min, rn2(ji,jj,jk) ) 
    518                END DO 
    519             END DO 
    520             ze_z =                  SUM( e1e2t(:,:) * zkz  (:,:)    * tmask_i(:,:) )   & 
    521                &     / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) ) 
    522             WRITE(numout,*) '                jk= ', jk,'   ', ze_z * 1.e4,' cm2/s' 
    523          END DO 
    524          DO jk = 1, jpk 
    525             zkz(:,:) = az_tmx(:,:,jk) /rn_n2min 
    526             ze_z =                  SUM( e1e2t(:,:) * zkz  (:,:)    * tmask_i(:,:) )   & 
    527                &     / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) ) 
    528             WRITE(numout,*)  
    529             WRITE(numout,*) '          N2 min - jk= ', jk,'   ', ze_z * 1.e4,' cm2/s min= ',MINVAL(zkz)*1.e4,   & 
    530                &       'max= ', MAXVAL(zkz)*1.e4, ' cm2/s' 
    531          END DO 
    532 !!gm  end bug mpp 
    533          ! 
    534       ENDIF 
    535       ! 
    536       CALL wrk_dealloc( jpi,jpj, zem2, zek1, zkz, zfact, zhdep ) 
    537       CALL wrk_dealloc( jpi,jpj,jpk, zpc, zav_tide ) 
    538       ! 
    539       IF( nn_timing == 1 )  CALL timing_stop('zdf_tmx_init') 
    540       ! 
    541    END SUBROUTINE zdf_tmx_init 
    542  
    543 #elif defined key_zdftmx_new 
    544    !!---------------------------------------------------------------------- 
    545    !!   'key_zdftmx_new'               Internal wave-driven vertical mixing 
     12 
    54613   !!---------------------------------------------------------------------- 
    54714   !!   zdf_tmx       : global     momentum & tracer Kz with wave induced Kz 
     
    56936   PUBLIC   zdf_tmx_init    ! called in nemogcm module  
    57037   PUBLIC   zdf_tmx_alloc   ! called in nemogcm module 
    571  
    572    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftmx = .TRUE.    !: wave-driven mixing flag 
    57338 
    57439   !                       !!* Namelist  namzdf_tmx : internal wave-driven mixing * 
     
    1027492   END SUBROUTINE zdf_tmx_init 
    1028493 
    1029 #else 
    1030    !!---------------------------------------------------------------------- 
    1031    !!   Default option          Dummy module                NO Tidal MiXing 
    1032    !!---------------------------------------------------------------------- 
    1033    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftmx = .FALSE.   !: tidal mixing flag 
    1034 CONTAINS 
    1035    SUBROUTINE zdf_tmx_init           ! Dummy routine 
    1036       WRITE(*,*) 'zdf_tmx: You should not have seen this print! error?' 
    1037    END SUBROUTINE zdf_tmx_init 
    1038    SUBROUTINE zdf_tmx( kt )          ! Dummy routine 
    1039       WRITE(*,*) 'zdf_tmx: You should not have seen this print! error?', kt 
    1040    END SUBROUTINE zdf_tmx 
    1041 #endif 
    1042  
    1043494   !!====================================================================== 
    1044495END MODULE zdftmx 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r7931 r7953  
    5555   USE ldfdyn         ! lateral viscosity setting      (ldfdyn_init routine) 
    5656   USE ldftra         ! lateral diffusivity setting    (ldftra_init routine) 
    57    USE zdfini         ! vertical physics setting          (zdf_init routine) 
     57!!gm   USE zdfphy         ! vertical physics manager      (zdf_phy_init routine) 
    5858   USE trdini         ! dyn/tra trends initialization     (trd_init routine) 
    5959   USE asminc         ! assimilation increments      
     
    429429      IF( ln_ctl        )   CALL prt_ctl_init   ! Print control 
    430430       
    431       CALL diurnal_sst_bulk_init            ! diurnal sst 
     431      CALL diurnal_sst_bulk_init             ! diurnal sst 
    432432      IF ( ln_diurnal ) CALL diurnal_sst_coolskin_init   ! cool skin    
    433433       
     
    455455                            CALL     sbc_init   ! surface boundary conditions (including sea-ice) 
    456456                            CALL     bdy_init   ! Open boundaries initialisation 
     457 
    457458      !                                      ! Ocean physics 
    458       !                                         ! Vertical physics 
    459                             CALL     zdf_init      ! namelist read 
    460                             CALL zdf_bfr_init      ! bottom friction 
    461       IF( lk_zdfric     )   CALL zdf_ric_init      ! Richardson number dependent Kz 
    462       IF( lk_zdftke     )   CALL zdf_tke_init      ! TKE closure scheme 
    463       IF( lk_zdfgls     )   CALL zdf_gls_init      ! GLS closure scheme 
    464       IF( lk_zdftmx     )   CALL zdf_tmx_init      ! tidal vertical mixing 
    465 !!gm      IF( ln_zdfddm     )   CALL zdf_ddm_init      ! double diffusive mixing 
    466           
     459                            CALL zdf_phy_init   ! Vertical physics 
     460                                      
    467461      !                                         ! Lateral physics 
    468462                            CALL ldf_tra_init      ! Lateral ocean tracer physics 
     
    470464                            CALL ldf_dyn_init      ! Lateral ocean momentum physics 
    471465 
    472       !                                         ! Active tracers 
     466      !                                      ! Active tracers 
    473467                            CALL tra_qsr_init      ! penetrative solar radiation qsr 
    474468                            CALL tra_bbc_init      ! bottom heat flux 
     
    479473                            CALL tra_zdf_init      ! vertical mixing and after tracer fields 
    480474 
    481       !                                         ! Dynamics 
     475      !                                      ! Dynamics 
    482476      IF( lk_c1d        )   CALL dyn_dmp_init      ! internal momentum damping 
    483477                            CALL dyn_adv_init      ! advection (vector or flux form) 
     
    511505      IF( ln_diaobs     )   CALL dia_obs( nit000 - 1 )   ! Observation operator for restart 
    512506 
    513       !                                         ! Assimilation increments 
     507      !                                      ! Assimilation increments 
    514508      IF( lk_asminc     )   CALL asm_inc_init   ! Initialize assimilation increments 
    515509      IF(lwp) WRITE(numout,*) 'Euler time step switch is ', neuler 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/step.F90

    r7931 r7953  
    7474      !!              -8- Outputs and diagnostics 
    7575      !!---------------------------------------------------------------------- 
    76       INTEGER ::   ji,jj,jk ! dummy loop indice 
    77       INTEGER ::   indic    ! error indicator if < 0 
    78       INTEGER ::   kcall    ! optional integer argument (dom_vvl_sf_nxt) 
     76      INTEGER ::   ji, jj, jk  ! dummy loop indice 
     77      INTEGER ::   indic        ! error indicator if < 0 
     78      INTEGER ::   kcall        ! optional integer argument (dom_vvl_sf_nxt) 
    7979      !! --------------------------------------------------------------------- 
    8080#if defined key_agrif 
     
    125125                         CALL bn2    ( tsn, rab_n, rn2  ) ! now    Brunt-Vaisala frequency 
    126126 
    127       ! 
    128127      !  VERTICAL PHYSICS 
    129                          CALL zdf_bfr( kstp )         ! bottom friction (if quadratic) 
    130       !                                               ! Vertical eddy viscosity and diffusivity coefficients 
    131       IF( lk_zdfric  )   CALL zdf_ric ( kstp )             ! Richardson number dependent Kz 
    132       IF( lk_zdftke  )   CALL zdf_tke ( kstp )             ! TKE closure scheme for Kz 
    133       IF( lk_zdfgls  )   CALL zdf_gls ( kstp )             ! GLS closure scheme for Kz 
    134       IF( ln_zdfqiao )   CALL zdf_qiao( kstp )             ! Qiao vertical mixing  
    135       ! 
    136       IF( lk_zdfcst  ) THEN                                ! Constant Kz (reset avt, avm[uv] to the background value) 
    137          avt (:,:,:) = rn_avt0 * wmask (:,:,:) 
    138          avm (:,:,:) = rn_avm0 * wmask (:,:,:) 
    139          avmu(:,:,:) = rn_avm0 * wumask(:,:,:) 
    140          avmv(:,:,:) = rn_avm0 * wvmask(:,:,:) 
    141       ENDIF 
    142       ! 
    143       IF( ln_rnf_mouth ) THEN                         ! increase diffusivity at rivers mouths 
    144          DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2._wp * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO 
    145       ENDIF 
    146       ! 
    147       IF( ln_zdfevd  )   CALL zdf_evd( kstp )         ! enhanced vertical eddy diffusivity 
    148       ! 
    149       IF( ln_zdfddm  ) THEN                           ! double diffusive mixing 
    150                          CALL zdf_ddm( kstp ) 
    151       ELSE                                            ! avs=avt 
    152          DO jk = 2, jpkm1   ;   avs(:,:,jk) = avt(:,:,jk)   ;   END DO 
    153       ENDIF 
    154       ! 
    155       IF( lk_zdftmx  )   CALL zdf_tmx( kstp )         ! tidal vertical mixing 
    156  
    157                          CALL zdf_mxl( kstp )         ! mixed layer depth 
    158  
    159                                                       ! write TKE or GLS information in the restart file 
    160       IF( lrst_oce .AND. lk_zdftke )   CALL tke_rst( kstp, 'WRITE' ) 
    161       IF( lrst_oce .AND. lk_zdfgls )   CALL gls_rst( kstp, 'WRITE' ) 
    162       ! 
     128                         CALL zdf_phy( kstp )         ! vertical physics update (bfr, avt, avs, avm + MLD) 
     129 
    163130      !  LATERAL  PHYSICS 
    164131      ! 
  • branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/step_oce.F90

    r7646 r7953  
    6363   USE ldftra           ! lateral eddy diffusive coef.     (ldf_tra routine) 
    6464 
     65   USE zdfphy         ! vertical physics manager      (zdf_phy_init routine) 
     66!!gm to be suppressed 
    6567   USE zdftmx           ! tide-induced vertical mixing     (zdf_tmx routine) 
    6668   USE zdfbfr           ! bottom friction                  (zdf_bfr routine) 
     
    7274   USE zdfmxl           ! Mixed-layer depth                (zdf_mxl routine) 
    7375   USE zdfqiao          !Qiao module wave induced mixing   (zdf_qiao routine) 
     76!!gm end 
    7477 
    7578   USE step_diu        ! Time stepping for diurnal sst 
Note: See TracChangeset for help on using the changeset viewer.