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 5581 for branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfddm.F90 – NEMO

Ignore:
Timestamp:
2015-07-10T13:28:53+02:00 (9 years ago)
Author:
timgraham
Message:

Merged head of trunk into branch

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfddm.F90

    r4624 r5581  
    66   !! History :  OPA  ! 2000-08  (G. Madec)  double diffusive mixing 
    77   !!   NEMO     1.0  ! 2002-06  (G. Madec)  F90: Free form and module 
    8    !!            3.3  !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase 
     8   !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase 
     9   !!            3.6  ! 2013-04  (G. Madec, F. Roquet) zrau compute locally using interpolation of alpha & beta 
    910   !!---------------------------------------------------------------------- 
    1011#if defined key_zdfddm   ||   defined key_esopa 
     
    1819   USE dom_oce         ! ocean space and time domain variables  
    1920   USE zdf_oce         ! ocean vertical physics variables 
     21   USE eosbn2         ! equation of state 
     22   ! 
    2023   USE in_out_manager  ! I/O manager 
    2124   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
     
    3437   LOGICAL , PUBLIC, PARAMETER ::   lk_zdfddm = .TRUE.  !: double diffusive mixing flag 
    3538 
    36    REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   avs    !: salinity vertical diffusivity coeff. at w-point 
    37    REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   rrau   !: heat/salt buoyancy flux ratio 
    38  
    39    !                      !!* Namelist namzdf_ddm : double diffusive mixing * 
    40    REAL(wp) ::   rn_avts   ! maximum value of avs for salt fingering 
    41    REAL(wp) ::   rn_hsbfr  ! heat/salt buoyancy flux ratio 
     39   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   avs   !: salinity vertical diffusivity coeff. at w-point 
     40 
     41   !                       !!* Namelist namzdf_ddm : double diffusive mixing * 
     42   REAL(wp) ::   rn_avts    ! maximum value of avs for salt fingering 
     43   REAL(wp) ::   rn_hsbfr   ! heat/salt buoyancy flux ratio 
    4244 
    4345   !! * Substitutions 
     46#  include "domzgr_substitute.h90" 
    4447#  include "vectopt_loop_substitute.h90" 
    4548   !!---------------------------------------------------------------------- 
    46    !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
     49   !! NEMO/OPA 3.7 , NEMO Consortium (2014) 
    4750   !! $Id$ 
    4851   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    5457      !!                ***  ROUTINE zdf_ddm_alloc  *** 
    5558      !!---------------------------------------------------------------------- 
    56       ALLOCATE( avs(jpi,jpj,jpk), rrau(jpi,jpj,jpk), STAT= zdf_ddm_alloc ) 
    57       ! 
     59      ALLOCATE( avs(jpi,jpj,jpk) , STAT= zdf_ddm_alloc ) 
    5860      IF( lk_mpp             )   CALL mpp_sum ( zdf_ddm_alloc ) 
    5961      IF( zdf_ddm_alloc /= 0 )   CALL ctl_warn('zdf_ddm_alloc: failed to allocate arrays') 
     
    7173      !!      diffusive mixing (i.e. salt fingering and diffusive layering) 
    7274      !!      following Merryfield et al. (1999). The rate of double diffusive  
    73       !!      mixing depend on the buoyancy ratio: Rrau=alpha/beta dk[T]/dk[S] 
    74       !!      which is computed in rn2.F 
     75      !!      mixing depend on the buoyancy ratio (R=alpha/beta dk[T]/dk[S]): 
    7576      !!         * salt fingering (Schmitt 1981): 
    76       !!      for Rrau > 1 and rn2 > 0 : zavfs = rn_avts / ( 1 + (Rrau/rn_hsbfr)^6 ) 
    77       !!      for Rrau > 1 and rn2 > 0 : zavfs = O 
    78       !!      otherwise                : zavft = 0.7 zavs / Rrau 
     77      !!      for R > 1 and rn2 > 0 : zavfs = rn_avts / ( 1 + (R/rn_hsbfr)^6 ) 
     78      !!      for R > 1 and rn2 > 0 : zavfs = O 
     79      !!      otherwise                : zavft = 0.7 zavs / R 
    7980      !!         * diffusive layering (Federov 1988): 
    80       !!      for 0< Rrau < 1 and rn2 > 0 : zavdt = 1.3635e-6   
    81       !!                                 * exp( 4.6 exp(-0.54 (1/Rrau-1) ) ) 
     81      !!      for 0< R < 1 and N^2 > 0 : zavdt = 1.3635e-6 * exp( 4.6 exp(-0.54 (1/R-1) ) ) 
    8282      !!      otherwise                   : zavdt = 0  
    83       !!      for .5 < Rrau < 1 and rn2 > 0 : zavds = zavdt (1.885 Rrau -0.85) 
    84       !!      for  0 < Rrau <.5 and rn2 > 0 : zavds = zavdt 0.15 Rrau       
     83      !!      for .5 < R < 1 and N^2 > 0 : zavds = zavdt (1.885 R -0.85) 
     84      !!      for  0 < R <.5 and N^2 > 0 : zavds = zavdt 0.15 R       
    8585      !!      otherwise                     : zavds = 0  
    8686      !!         * update the eddy diffusivity: 
     
    9696      ! 
    9797      INTEGER  ::   ji, jj , jk     ! dummy loop indices 
    98       REAL(wp) ::   zinr, zrr       ! temporary scalars 
    99       REAL(wp) ::   zavft, zavfs    !    -         - 
    100       REAL(wp) ::   zavdt, zavds    !    -         - 
    101       REAL(wp), POINTER, DIMENSION(:,:) ::   zmsks, zmskf, zmskd1, zmskd2, zmskd3 
     98      REAL(wp) ::   zaw, zbw, zrw   ! local scalars 
     99      REAL(wp) ::   zdt, zds 
     100      REAL(wp) ::   zinr, zrr       !   -      - 
     101      REAL(wp) ::   zavft, zavfs    !   -      - 
     102      REAL(wp) ::   zavdt, zavds    !   -      - 
     103      REAL(wp), POINTER, DIMENSION(:,:) ::   zrau, zmsks, zmskf, zmskd1, zmskd2, zmskd3 
    102104      !!---------------------------------------------------------------------- 
    103105      ! 
    104106      IF( nn_timing == 1 )  CALL timing_start('zdf_ddm') 
    105107      ! 
    106       CALL wrk_alloc( jpi,jpj, zmsks, zmskf, zmskd1, zmskd2, zmskd3 ) 
    107  
     108      CALL wrk_alloc( jpi,jpj, zrau, zmsks, zmskf, zmskd1, zmskd2, zmskd3 ) 
     109      ! 
    108110      !                                                ! =============== 
    109111      DO jk = 2, jpkm1                                 ! Horizontal slab 
     
    111113         ! Define the mask  
    112114         ! --------------- 
    113          rrau(:,:,jk) = MAX( 1.e-20, rrau(:,:,jk) )         ! only retains positive value of rrau 
     115         DO jj = 1, jpj                                ! R=zrau = (alpha / beta) (dk[t] / dk[s]) 
     116            DO ji = 1, jpi 
     117               zrw =   ( fsdepw(ji,jj,jk  ) - fsdept(ji,jj,jk) )   & 
     118                  &  / ( fsdept(ji,jj,jk-1) - fsdept(ji,jj,jk) )  
     119               ! 
     120               zaw = (  rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem) * zrw  )  & 
     121                   &    * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) 
     122               zbw = (  rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal) * zrw  )  & 
     123                   &    * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) 
     124               ! 
     125               zdt = zaw * ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) 
     126               zds = zbw * ( tsn(ji,jj,jk-1,jp_sal) - tsn(ji,jj,jk,jp_sal) )  
     127               IF( ABS( zds) <= 1.e-20_wp )   zds = 1.e-20_wp 
     128               zrau(ji,jj) = MAX(  1.e-20, zdt / zds  )    ! only retains positive value of zrau 
     129            END DO 
     130         END DO 
    114131 
    115132         DO jj = 1, jpj                                     ! indicators: 
     
    119136               ELSE                                       ;   zmsks(ji,jj) = 1._wp 
    120137               ENDIF 
    121                ! salt fingering indicator: msksf=1 if rrau>1; 0 elsewhere             
    122                IF( rrau(ji,jj,jk) <= 1.          ) THEN   ;   zmskf(ji,jj) = 0._wp 
     138               ! salt fingering indicator: msksf=1 if R>1; 0 elsewhere             
     139               IF( zrau(ji,jj) <= 1.             ) THEN   ;   zmskf(ji,jj) = 0._wp 
    123140               ELSE                                       ;   zmskf(ji,jj) = 1._wp 
    124141               ENDIF 
    125142               ! diffusive layering indicators:  
    126                !     ! mskdl1=1 if 0<rrau<1; 0 elsewhere 
    127                IF( rrau(ji,jj,jk) >= 1.          ) THEN   ;   zmskd1(ji,jj) = 0._wp 
     143               !     ! mskdl1=1 if 0< R <1; 0 elsewhere 
     144               IF( zrau(ji,jj) >= 1.             ) THEN   ;   zmskd1(ji,jj) = 0._wp 
    128145               ELSE                                       ;   zmskd1(ji,jj) = 1._wp 
    129146               ENDIF 
    130                !     ! mskdl2=1 if 0<rrau<0.5; 0 elsewhere 
    131                IF( rrau(ji,jj,jk) >= 0.5         ) THEN   ;   zmskd2(ji,jj) = 0._wp 
     147               !     ! mskdl2=1 if 0< R <0.5; 0 elsewhere 
     148               IF( zrau(ji,jj) >= 0.5            ) THEN   ;   zmskd2(ji,jj) = 0._wp 
    132149               ELSE                                       ;   zmskd2(ji,jj) = 1._wp 
    133150               ENDIF 
    134                !   mskdl3=1 if 0.5<rrau<1; 0 elsewhere 
    135                IF( rrau(ji,jj,jk) <= 0.5 .OR. rrau(ji,jj,jk) >= 1. ) THEN   ;   zmskd3(ji,jj) = 0._wp 
    136                ELSE                                                         ;   zmskd3(ji,jj) = 1._wp 
     151               !   mskdl3=1 if 0.5< R <1; 0 elsewhere 
     152               IF( zrau(ji,jj) <= 0.5 .OR. zrau(ji,jj) >= 1. ) THEN   ;   zmskd3(ji,jj) = 0._wp 
     153               ELSE                                                   ;   zmskd3(ji,jj) = 1._wp 
    137154               ENDIF 
    138155            END DO 
    139156         END DO 
    140157         ! mask zmsk in order to have avt and avs masked 
    141          zmsks(:,:) = zmsks(:,:) * tmask(:,:,jk) 
     158         zmsks(:,:) = zmsks(:,:) * wmask(:,:,jk) 
    142159 
    143160 
     
    149166!CDIR NOVERRCHK 
    150167            DO ji = 1, jpi 
    151                zinr = 1./rrau(ji,jj,jk) 
     168               zinr = 1._wp / zrau(ji,jj) 
    152169               ! salt fingering 
    153                zrr = rrau(ji,jj,jk)/rn_hsbfr 
     170               zrr = zrau(ji,jj) / rn_hsbfr 
    154171               zrr = zrr * zrr 
    155172               zavfs = rn_avts / ( 1 + zrr*zrr*zrr ) * zmsks(ji,jj) * zmskf(ji,jj) 
     
    157174               ! diffusive layering 
    158175               zavdt = 1.3635e-6 * EXP(  4.6 * EXP( -0.54*(zinr-1.) )  ) * zmsks(ji,jj) * zmskd1(ji,jj) 
    159                zavds = zavdt * zmsks(ji,jj) * (  (1.85 * rrau(ji,jj,jk) - 0.85 ) * zmskd3(ji,jj)   & 
    160                   &                            +  0.15 * rrau(ji,jj,jk)          * zmskd2(ji,jj)  ) 
     176               zavds = zavdt * zmsks(ji,jj) * (  ( 1.85 * zrau(ji,jj) - 0.85 ) * zmskd3(ji,jj)   & 
     177                  &                             +  0.15 * zrau(ji,jj)          * zmskd2(ji,jj)  ) 
    161178               ! add to the eddy viscosity coef. previously computed 
    162179               avs (ji,jj,jk) = avt(ji,jj,jk) + zavfs + zavds 
     
    174191               avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk),    & 
    175192                  &                  avt(ji,jj,jk), avt(ji+1,jj,jk),   & 
    176                   &                  avs(ji,jj,jk), avs(ji+1,jj,jk) )  * umask(ji,jj,jk) 
     193                  &                  avs(ji,jj,jk), avs(ji+1,jj,jk) )  * wumask(ji,jj,jk) 
    177194               avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk),    & 
    178195                  &                  avt(ji,jj,jk), avt(ji,jj+1,jk),   & 
    179                   &                  avs(ji,jj,jk), avs(ji,jj+1,jk) )  * vmask(ji,jj,jk) 
     196                  &                  avs(ji,jj,jk), avs(ji,jj+1,jk) )  * wvmask(ji,jj,jk) 
    180197            END DO 
    181198         END DO 
     
    196213      ENDIF 
    197214      ! 
    198       CALL wrk_dealloc( jpi,jpj, zmsks, zmskf, zmskd1, zmskd2, zmskd3 ) 
     215      CALL wrk_dealloc( jpi,jpj, zrau, zmsks, zmskf, zmskd1, zmskd2, zmskd3 ) 
    199216      ! 
    200217      IF( nn_timing == 1 )  CALL timing_stop('zdf_ddm') 
     
    212229      !!              called by zdf_ddm at the first timestep (nit000) 
    213230      !!---------------------------------------------------------------------- 
     231      INTEGER ::   ios   ! local integer 
     232      !! 
    214233      NAMELIST/namzdf_ddm/ rn_avts, rn_hsbfr 
    215       INTEGER  ::   ios                 ! Local integer output status for namelist read 
    216234      !!---------------------------------------------------------------------- 
    217235      ! 
     
    237255      IF( zdf_ddm_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'zdf_ddm_init : unable to allocate arrays' ) 
    238256      !                               ! initialization to masked Kz 
    239       avs(:,:,:) = rn_avt0 * tmask(:,:,:)  
     257      avs(:,:,:) = rn_avt0 * wmask(:,:,:)  
    240258      ! 
    241259   END SUBROUTINE zdf_ddm_init 
Note: See TracChangeset for help on using the changeset viewer.