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.
zdfddm.F90 in trunk/NEMO/OPA_SRC/ZDF – NEMO

source: trunk/NEMO/OPA_SRC/ZDF/zdfddm.F90 @ 3

Last change on this file since 3 was 3, checked in by opalod, 20 years ago

Initial revision
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 10.7 KB
Line 
1MODULE zdfddm
2   !!======================================================================
3   !!                       ***  MODULE  zdfddm  ***
4   !! Ocean physics : double diffusion mixing parameterization
5   !!======================================================================
6#if defined key_zdfddm   ||   defined key_esopa
7   !!----------------------------------------------------------------------
8   !!   'key_zdfddm' :                                     double diffusion
9   !!----------------------------------------------------------------------
10   !!   zdf_ddm       : compute the Ks for salinity
11   !!   zdf_ddm_init  : read namelist and control the parameters
12   !!----------------------------------------------------------------------
13   !! * Modules used
14   USE oce             ! ocean dynamics and tracers variables
15   USE dom_oce         ! ocean space and time domain variables
16   USE zdf_oce         ! ocean vertical physics variables
17   USE in_out_manager  ! I/O manager
18   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
19
20   IMPLICIT NONE
21   PRIVATE
22
23   !! * Routine accessibility
24   PUBLIC zdf_ddm     ! called by step.F90
25
26   !! * Shared module variables
27   LOGICAL, PUBLIC ::   lk_zdfddm = .TRUE.    !: double diffusive mixing flag
28   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &
29      avs ,               &  !: salinity vertical diffusivity coeff. at w-point
30      rrau                   !: heat/salt buoyancy flux ratio
31
32   !! * Module variables
33   REAL(wp) ::            & !!! * double diffusive mixing namelist *
34      avts  = 1.e-4_wp ,  &  ! maximum value of avs for salt fingering
35      hsbfr = 1.6_wp         ! heat/salt buoyancy flux ratio
36
37   !! * Substitutions
38#  include "vectopt_loop_substitute.h90"
39   !!----------------------------------------------------------------------
40   !!   OPA 9.0 , LODYC-IPSL  (2003)
41   !!----------------------------------------------------------------------
42
43CONTAINS
44
45   SUBROUTINE zdf_ddm( kt )
46      !!----------------------------------------------------------------------
47      !!                  ***  ROUTINE zdf_ddm  ***
48      !!                   
49      !! ** Purpose :   Add to the vertical eddy diffusivity coefficient the
50      !!      effect of salt fingering and diffusive convection.
51      !!
52      !! ** Method  :   Diapycnal mixing is increased in case of double
53      !!      diffusive mixing (i.e. salt fingering and diffusive layering)
54      !!      following Merryfield et al. (1999). The rate of double diffusive
55      !!      mixing depend on the buoyancy ratio: Rrau=alpha/beta dk[T]/dk[S]
56      !!      which is computed in rn2.F
57      !!         * salt fingering (Schmitt 1981):
58      !!      for Rrau > 1 and rn2 > 0 : zavfs = avts / ( 1 + (Rrau/hsbfr)^6 )
59      !!      for Rrau > 1 and rn2 > 0 : zavfs = O
60      !!      otherwise                : zavft = 0.7 zavs / Rrau
61      !!         * diffusive layering (Federov 1988):
62      !!      for 0< Rrau < 1 and rn2 > 0 : zavdt = 1.3635e-6 
63      !!                                 * exp( 4.6 exp(-0.54 (1/Rrau-1) ) )
64      !!      otherwise                   : zavdt = 0
65      !!      for .5 < Rrau < 1 and rn2 > 0 : zavds = zavdt (1.885 Rrau -0.85)
66      !!      for  0 < Rrau <.5 and rn2 > 0 : zavds = zavdt 0.15 Rrau     
67      !!      otherwise                     : zavds = 0
68      !!         * update the eddy diffusivity:
69      !!      avt = avt + zavft + zavdt
70      !!      avs = avs + zavfs + zavds
71      !!      avmu, avmv are required to remain at least above avt and avs.
72      !!     
73      !! ** Action  :   avt, avs : update vertical eddy diffusivity coef.
74      !!                           for temperature and salinity
75      !!
76      !! References :
77      !!      Merryfield et al., JPO, 29, 1124-1142, 1999.
78      !! History :
79      !!        !  00-08  (G. Madec)  double diffusive mixing
80      !!   8.5  !  02-06  (G. Madec)  F90: Free form and module
81      !!----------------------------------------------------------------------
82      !! * Arguments
83      INTEGER, INTENT( in ) ::   kt         ! ocean time-step indexocean time step
84
85      !! * Local declarations
86      INTEGER ::   ji, jj , jk              ! dummy loop indices
87      REAL(wp), DIMENSION(jpi,jpj) ::   &
88         zmsks, zmskf,                    & ! temporary workspace
89         zmskd1, zmskd2, zmskd3             !    "           "
90      REAL(wp) ::   &
91         zinr, zrr,                       & ! temporary scalars
92         zavft, zavfs,                    & !    "         "
93         zavdt, zavds                       !    "         "
94      !!----------------------------------------------------------------------
95
96
97      IF ( kt == nit000 )   CALL zdf_ddm_init          ! Initialization (first time-step only)
98
99
100      ! Compute avs
101      ! -----------
102      !                                                ! ===============
103      DO jk = 2, jpkm1                                 ! Horizontal slab
104         !                                             ! ===============
105         ! Define the mask
106         ! ---------------
107         ! only retains positive value of rrau
108         rrau(:,:,jk) = MAX( 1.e-20, rrau(:,:,jk) )
109
110         ! indicators:
111         DO jj = 1, jpj
112            DO ji = 1, jpi
113               ! stability indicator: msks=1 if rn2>0; 0 elsewhere
114               IF( rn2(ji,jj,jk) + 1.e-12  <= 0. ) THEN
115                  zmsks(ji,jj) = 0.e0
116               ELSE
117                  zmsks(ji,jj) = 1.e0
118               ENDIF
119               ! salt fingering indicator: msksf=1 if rrau>1; 0 elsewhere           
120               IF( rrau(ji,jj,jk) <= 1. ) THEN
121                  zmskf(ji,jj) = 0.e0
122               ELSE
123                  zmskf(ji,jj) = 1.e0
124               ENDIF
125               ! diffusive layering indicators:
126               !   mskdl1=1 if 0<rrau<1; 0 elsewhere
127               IF( rrau(ji,jj,jk) >= 1. ) THEN
128                  zmskd1(ji,jj) = 0.e0
129               ELSE
130                  zmskd1(ji,jj) = 1.e0
131               ENDIF
132               !   mskdl2=1 if 0<rrau<0.5; 0 elsewhere
133               IF( rrau(ji,jj,jk) >= 0.5 ) THEN
134                  zmskd2(ji,jj) = 0.e0
135               ELSE
136                  zmskd2(ji,jj) = 1.e0
137               ENDIF
138               !   mskdl3=1 if 0.5<rrau<1; 0 elsewhere
139               IF( rrau(ji,jj,jk) <= 0.5 .OR. rrau(ji,jj,jk) >= 1. ) THEN
140                  zmskd3(ji,jj) = 0.e0
141               ELSE
142                  zmskd3(ji,jj) = 1.e0
143               ENDIF
144            END DO
145         END DO
146         ! mask zmsk in order to have avt and avs masked
147         zmsks(:,:) = zmsks(:,:) * tmask(:,:,jk)
148
149
150         ! Update avt and avs
151         ! ------------------
152         ! Constant eddy coefficient: reset to the background value
153!CDIR NOVERRCHK
154         DO jj = 1, jpj
155!CDIR NOVERRCHK
156            DO ji = 1, jpi
157               zinr = 1./rrau(ji,jj,jk)
158               ! salt fingering
159               zrr = rrau(ji,jj,jk)/hsbfr
160               zrr = zrr * zrr
161               zavfs = avts / ( 1 + zrr*zrr*zrr ) * zmsks(ji,jj) *zmskf(ji,jj)
162               zavft = 0.7 * zavfs / rrau(ji,jj,jk)
163               ! diffusive layering
164               zavdt = 1.3635e-6 * EXP(4.6*EXP(-0.54*(zinr-1.) ) )   &
165                                 * zmsks(ji,jj) * zmskd1(ji,jj)
166               zavds = zavdt * zmsks(ji,jj)   &
167                     * ( (1.85 * rrau(ji,jj,jk) - 0.85 ) * zmskd3(ji,jj)   &
168                         + zavdt * 0.15 * rrau(ji,jj,jk) * zmskd2(ji,jj)  )
169               ! add to the eddy viscosity coef. previously computed
170               avs (ji,jj,jk) = avt(ji,jj,jk) + zavfs + zavds
171               avt (ji,jj,jk) = avt(ji,jj,jk) + zavft + zavdt
172            END DO
173         END DO
174
175
176         ! Increase avmu, avmv if necessary
177         ! --------------------------------
178         DO jj = 1, jpjm1
179            DO ji = 1, fs_jpim1   ! vector opt.
180               avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk),    &
181                                     avt(ji,jj,jk), avt(ji+1,jj,jk),   &
182                                     avs(ji,jj,jk), avs(ji+1,jj,jk) )   &
183                              * umask(ji,jj,jk)
184               avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk),    &
185                                     avt(ji,jj,jk), avt(ji,jj+1,jk),   &
186                                     avs(ji,jj,jk), avs(ji,jj+1,jk) )   &
187                              * vmask(ji,jj,jk)
188            END DO
189         END DO
190         !                                                ! ===============
191      END DO                                              !   End of slab
192      !                                                   ! ===============
193     
194      ! Lateral boundary conditions on ( avt, avs, avmu, avmv )   (unchanged sign)
195      ! -------------------------------========================
196      CALL lbc_lnk( avt , 'W', 1. )
197      CALL lbc_lnk( avs , 'W', 1. )
198      CALL lbc_lnk( avmu, 'U', 1. ) 
199      CALL lbc_lnk( avmv, 'V', 1. )
200     
201   END SUBROUTINE zdf_ddm
202   
203   
204   SUBROUTINE zdf_ddm_init
205      !!----------------------------------------------------------------------
206      !!                  ***  ROUTINE zdf_ddm_init  ***
207      !!
208      !! ** Purpose :   Initialization of double diffusion mixing scheme
209      !!
210      !! ** Method  :   Read the nammbf namelist and check the parameter values
211      !!      called by zdf_ddm at the first timestep (nit000)
212      !!
213      !! History :
214      !!   8.5  !  02-08  (G. Madec)  Original code
215      !!----------------------------------------------------------------------
216      NAMELIST/namddm/ avts, hsbfr
217      !!----------------------------------------------------------------------
218
219      ! Read Namelist namddm : double diffusion mixing scheme
220      ! --------------------
221      REWIND ( numnam )
222      READ   ( numnam, namddm )
223
224
225      ! Parameter control and print
226      ! ---------------------------
227      IF(lwp) THEN
228         WRITE(numout,*)
229         WRITE(numout,*) 'zdf_ddm : double diffusive mixing'
230         WRITE(numout,*) '~~~~~~~'
231         WRITE(numout,*) '          Namelist namddm : set dd mixing parameter'
232         WRITE(numout,*)
233         WRITE(numout,*) '             maximum avs for dd mixing      avts   = ', avts
234         WRITE(numout,*) '             heat/salt buoyancy flux ratio  hsbfr  = ', hsbfr
235         WRITE(numout,*)
236      ENDIF
237
238   END SUBROUTINE zdf_ddm_init
239
240#else
241   !!----------------------------------------------------------------------
242   !!   Default option :          Dummy module          No double diffusion
243   !!----------------------------------------------------------------------
244   LOGICAL, PUBLIC ::   lk_zdfddm = .FALSE.   !: double diffusion flag
245CONTAINS
246   SUBROUTINE zdf_ddm( kt )           ! Dummy routine
247      WRITE(*,*) kt                          ! avoid compil warning
248   END SUBROUTINE zdf_ddm
249#endif
250
251   !!======================================================================
252END MODULE zdfddm
Note: See TracBrowser for help on using the repository browser.