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 @ 247

Last change on this file since 247 was 247, checked in by opalod, 19 years ago

CL : Add CVS Header and CeCILL licence information
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.1 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, PARAMETER ::   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 , LOCEAN-IPSL (2005)
41   !! $Header$
42   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
43   !!----------------------------------------------------------------------
44
45CONTAINS
46
47   SUBROUTINE zdf_ddm( kt )
48      !!----------------------------------------------------------------------
49      !!                  ***  ROUTINE zdf_ddm  ***
50      !!                   
51      !! ** Purpose :   Add to the vertical eddy diffusivity coefficient the
52      !!      effect of salt fingering and diffusive convection.
53      !!
54      !! ** Method  :   Diapycnal mixing is increased in case of double
55      !!      diffusive mixing (i.e. salt fingering and diffusive layering)
56      !!      following Merryfield et al. (1999). The rate of double diffusive
57      !!      mixing depend on the buoyancy ratio: Rrau=alpha/beta dk[T]/dk[S]
58      !!      which is computed in rn2.F
59      !!         * salt fingering (Schmitt 1981):
60      !!      for Rrau > 1 and rn2 > 0 : zavfs = avts / ( 1 + (Rrau/hsbfr)^6 )
61      !!      for Rrau > 1 and rn2 > 0 : zavfs = O
62      !!      otherwise                : zavft = 0.7 zavs / Rrau
63      !!         * diffusive layering (Federov 1988):
64      !!      for 0< Rrau < 1 and rn2 > 0 : zavdt = 1.3635e-6 
65      !!                                 * exp( 4.6 exp(-0.54 (1/Rrau-1) ) )
66      !!      otherwise                   : zavdt = 0
67      !!      for .5 < Rrau < 1 and rn2 > 0 : zavds = zavdt (1.885 Rrau -0.85)
68      !!      for  0 < Rrau <.5 and rn2 > 0 : zavds = zavdt 0.15 Rrau     
69      !!      otherwise                     : zavds = 0
70      !!         * update the eddy diffusivity:
71      !!      avt = avt + zavft + zavdt
72      !!      avs = avs + zavfs + zavds
73      !!      avmu, avmv are required to remain at least above avt and avs.
74      !!     
75      !! ** Action  :   avt, avs : update vertical eddy diffusivity coef.
76      !!                           for temperature and salinity
77      !!
78      !! References :
79      !!      Merryfield et al., JPO, 29, 1124-1142, 1999.
80      !! History :
81      !!        !  00-08  (G. Madec)  double diffusive mixing
82      !!   8.5  !  02-06  (G. Madec)  F90: Free form and module
83      !!----------------------------------------------------------------------
84      !! * Arguments
85      INTEGER, INTENT( in ) ::   kt         ! ocean time-step indexocean time step
86
87      !! * Local declarations
88      INTEGER ::   ji, jj , jk              ! dummy loop indices
89      REAL(wp), DIMENSION(jpi,jpj) ::   &
90         zmsks, zmskf,                    & ! temporary workspace
91         zmskd1, zmskd2, zmskd3             !    "           "
92      REAL(wp) ::   &
93         zinr, zrr,                       & ! temporary scalars
94         zavft, zavfs,                    & !    "         "
95         zavdt, zavds                       !    "         "
96      !!----------------------------------------------------------------------
97
98
99      IF ( kt == nit000 )   CALL zdf_ddm_init          ! Initialization (first time-step only)
100
101
102      ! Compute avs
103      ! -----------
104      !                                                ! ===============
105      DO jk = 2, jpkm1                                 ! Horizontal slab
106         !                                             ! ===============
107         ! Define the mask
108         ! ---------------
109         ! only retains positive value of rrau
110         rrau(:,:,jk) = MAX( 1.e-20, rrau(:,:,jk) )
111
112         ! indicators:
113         DO jj = 1, jpj
114            DO ji = 1, jpi
115               ! stability indicator: msks=1 if rn2>0; 0 elsewhere
116               IF( rn2(ji,jj,jk) + 1.e-12  <= 0. ) THEN
117                  zmsks(ji,jj) = 0.e0
118               ELSE
119                  zmsks(ji,jj) = 1.e0
120               ENDIF
121               ! salt fingering indicator: msksf=1 if rrau>1; 0 elsewhere           
122               IF( rrau(ji,jj,jk) <= 1. ) THEN
123                  zmskf(ji,jj) = 0.e0
124               ELSE
125                  zmskf(ji,jj) = 1.e0
126               ENDIF
127               ! diffusive layering indicators:
128               !   mskdl1=1 if 0<rrau<1; 0 elsewhere
129               IF( rrau(ji,jj,jk) >= 1. ) THEN
130                  zmskd1(ji,jj) = 0.e0
131               ELSE
132                  zmskd1(ji,jj) = 1.e0
133               ENDIF
134               !   mskdl2=1 if 0<rrau<0.5; 0 elsewhere
135               IF( rrau(ji,jj,jk) >= 0.5 ) THEN
136                  zmskd2(ji,jj) = 0.e0
137               ELSE
138                  zmskd2(ji,jj) = 1.e0
139               ENDIF
140               !   mskdl3=1 if 0.5<rrau<1; 0 elsewhere
141               IF( rrau(ji,jj,jk) <= 0.5 .OR. rrau(ji,jj,jk) >= 1. ) THEN
142                  zmskd3(ji,jj) = 0.e0
143               ELSE
144                  zmskd3(ji,jj) = 1.e0
145               ENDIF
146            END DO
147         END DO
148         ! mask zmsk in order to have avt and avs masked
149         zmsks(:,:) = zmsks(:,:) * tmask(:,:,jk)
150
151
152         ! Update avt and avs
153         ! ------------------
154         ! Constant eddy coefficient: reset to the background value
155!CDIR NOVERRCHK
156         DO jj = 1, jpj
157!CDIR NOVERRCHK
158            DO ji = 1, jpi
159               zinr = 1./rrau(ji,jj,jk)
160               ! salt fingering
161               zrr = rrau(ji,jj,jk)/hsbfr
162               zrr = zrr * zrr
163               zavfs = avts / ( 1 + zrr*zrr*zrr ) * zmsks(ji,jj) *zmskf(ji,jj)
164               zavft = 0.7 * zavfs / rrau(ji,jj,jk)
165               ! diffusive layering
166               zavdt = 1.3635e-6 * EXP(4.6*EXP(-0.54*(zinr-1.) ) )   &
167                                 * zmsks(ji,jj) * zmskd1(ji,jj)
168               zavds = zavdt * zmsks(ji,jj)   &
169                     * ( (1.85 * rrau(ji,jj,jk) - 0.85 ) * zmskd3(ji,jj)   &
170                         + zavdt * 0.15 * rrau(ji,jj,jk) * zmskd2(ji,jj)  )
171               ! add to the eddy viscosity coef. previously computed
172               avs (ji,jj,jk) = avt(ji,jj,jk) + zavfs + zavds
173               avt (ji,jj,jk) = avt(ji,jj,jk) + zavft + zavdt
174            END DO
175         END DO
176
177
178         ! Increase avmu, avmv if necessary
179         ! --------------------------------
180         DO jj = 1, jpjm1
181            DO ji = 1, fs_jpim1   ! vector opt.
182               avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk),    &
183                                     avt(ji,jj,jk), avt(ji+1,jj,jk),   &
184                                     avs(ji,jj,jk), avs(ji+1,jj,jk) )   &
185                              * umask(ji,jj,jk)
186               avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk),    &
187                                     avt(ji,jj,jk), avt(ji,jj+1,jk),   &
188                                     avs(ji,jj,jk), avs(ji,jj+1,jk) )   &
189                              * vmask(ji,jj,jk)
190            END DO
191         END DO
192         !                                                ! ===============
193      END DO                                              !   End of slab
194      !                                                   ! ===============
195     
196      ! Lateral boundary conditions on ( avt, avs, avmu, avmv )   (unchanged sign)
197      ! -------------------------------========================
198      CALL lbc_lnk( avt , 'W', 1. )
199      CALL lbc_lnk( avs , 'W', 1. )
200      CALL lbc_lnk( avmu, 'U', 1. ) 
201      CALL lbc_lnk( avmv, 'V', 1. )
202
203      IF(l_ctl) THEN
204         WRITE(numout,*) ' ddm  t : ', SUM( avt (1:nictl+1,1:njctl+1,:) ), ' s : ', SUM( avs (1:nictl+1,1:njctl+1,:) )
205         WRITE(numout,*) '      u : ', SUM( avmu(1:nictl+1,1:njctl+1,:) ), ' v : ', SUM( avmv(1:nictl+1,1:njctl+1,:) )
206      ENDIF
207     
208   END SUBROUTINE zdf_ddm
209   
210   
211   SUBROUTINE zdf_ddm_init
212      !!----------------------------------------------------------------------
213      !!                  ***  ROUTINE zdf_ddm_init  ***
214      !!
215      !! ** Purpose :   Initialization of double diffusion mixing scheme
216      !!
217      !! ** Method  :   Read the nammbf namelist and check the parameter values
218      !!      called by zdf_ddm at the first timestep (nit000)
219      !!
220      !! History :
221      !!   8.5  !  02-08  (G. Madec)  Original code
222      !!----------------------------------------------------------------------
223      NAMELIST/namddm/ avts, hsbfr
224      !!----------------------------------------------------------------------
225
226      ! Read Namelist namddm : double diffusion mixing scheme
227      ! --------------------
228      REWIND ( numnam )
229      READ   ( numnam, namddm )
230
231
232      ! Parameter control and print
233      ! ---------------------------
234      IF(lwp) THEN
235         WRITE(numout,*)
236         WRITE(numout,*) 'zdf_ddm : double diffusive mixing'
237         WRITE(numout,*) '~~~~~~~'
238         WRITE(numout,*) '          Namelist namddm : set dd mixing parameter'
239         WRITE(numout,*) '             maximum avs for dd mixing      avts   = ', avts
240         WRITE(numout,*) '             heat/salt buoyancy flux ratio  hsbfr  = ', hsbfr
241         WRITE(numout,*)
242      ENDIF
243
244   END SUBROUTINE zdf_ddm_init
245
246#else
247   !!----------------------------------------------------------------------
248   !!   Default option :          Dummy module          No double diffusion
249   !!----------------------------------------------------------------------
250   LOGICAL, PUBLIC, PARAMETER ::   lk_zdfddm = .FALSE.   !: double diffusion flag
251CONTAINS
252   SUBROUTINE zdf_ddm( kt )           ! Dummy routine
253      WRITE(*,*) 'zdf_ddm: You should not have seen this print! error?', kt
254   END SUBROUTINE zdf_ddm
255#endif
256
257   !!======================================================================
258END MODULE zdfddm
Note: See TracBrowser for help on using the repository browser.