Newer
Older
MODULE zdfddm
!!======================================================================
!! *** MODULE zdfddm ***
!! Ocean physics : double diffusion mixing parameterization
!!======================================================================
!! History : OPA ! 2000-08 (G. Madec) double diffusive mixing
!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
!! 3.6 ! 2013-04 (G. Madec, F. Roquet) zrau compute locally using interpolation of alpha & beta
!! 4.0 ! 2017-04 (G. Madec) remove CPP ddm key & avm at t-point only
!!----------------------------------------------------------------------
!!----------------------------------------------------------------------
!! zdf_ddm : compute the Kz for salinity
!!----------------------------------------------------------------------
USE oce ! ocean dynamics and tracers variables
USE dom_oce ! ocean space and time domain variables
USE zdf_oce ! ocean vertical physics variables
USE eosbn2 ! equation of state
!
USE in_out_manager ! I/O manager
USE lbclnk ! ocean lateral boundary conditions (or mpp link)
USE prtctl ! Print control
USE lib_mpp ! MPP library
IMPLICIT NONE
PRIVATE
PUBLIC zdf_ddm ! called by step.F90
!! * Substitutions
# include "do_loop_substitute.h90"
sparonuz
committed
# include "single_precision_substitute.h90"
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: zdfddm.F90 14853 2021-05-12 13:07:30Z smasson $
!! Software governed by the CeCILL license (see ./LICENSE)
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE zdf_ddm( kt, Kmm, p_avm, p_avt, p_avs )
!!----------------------------------------------------------------------
!! *** ROUTINE zdf_ddm ***
!!
!! ** Purpose : Add to the vertical eddy diffusivity coefficient the
!! effect of salt fingering and diffusive convection.
!!
!! ** Method : Diapycnal mixing is increased in case of double
!! diffusive mixing (i.e. salt fingering and diffusive layering)
!! following Merryfield et al. (1999). The rate of double diffusive
!! mixing depend on the buoyancy ratio (R=alpha/beta dk[T]/dk[S]):
!! * salt fingering (Schmitt 1981):
!! for R > 1 and rn2 > 0 : zavfs = rn_avts / ( 1 + (R/rn_hsbfr)^6 )
!! for R > 1 and rn2 > 0 : zavfs = O
!! otherwise : zavft = 0.7 zavs / R
!! * diffusive layering (Federov 1988):
!! for 0< R < 1 and N^2 > 0 : zavdt = 1.3635e-6 * exp( 4.6 exp(-0.54 (1/R-1) ) )
!! otherwise : zavdt = 0
!! for .5 < R < 1 and N^2 > 0 : zavds = zavdt (1.885 R -0.85)
!! for 0 < R <.5 and N^2 > 0 : zavds = zavdt 0.15 R
!! otherwise : zavds = 0
!! * update the eddy diffusivity:
!! avt = avt + zavft + zavdt
!! avs = avs + zavfs + zavds
!! avm is required to remain at least above avt and avs.
!!
!! ** Action : avt, avs : updated vertical eddy diffusivity coef. for T & S
!!
!! References : Merryfield et al., JPO, 29, 1124-1142, 1999.
!!----------------------------------------------------------------------
INTEGER, INTENT(in ) :: kt ! ocean time-step index
INTEGER, INTENT(in ) :: Kmm ! ocean time level index
REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: p_avm ! Kz on momentum (w-points)
REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: p_avt ! Kz on temperature (w-points)
REAL(wp), DIMENSION(:,:,:), INTENT( out) :: p_avs ! Kz on salinity (w-points)
!
INTEGER :: ji, jj , jk ! dummy loop indices
REAL(wp) :: zaw, zbw, zrw ! local scalars
REAL(wp) :: zdt, zds
REAL(wp) :: zinr ! - -
REAL(dp) :: zrr ! - -
REAL(wp) :: zavft ! - -
REAL(dp) :: zavfs ! - -
REAL(wp) :: zavdt, zavds ! - -
REAL(wp), DIMENSION(A2D(nn_hls)) :: zrau, zmsks, zmskf, zmskd1, zmskd2, zmskd3
!!----------------------------------------------------------------------
!
! ! ===============
DO jk = 2, jpkm1 ! Horizontal slab
! ! ===============
! Define the mask
! ---------------
!!gm WORK to be done: change the code from vector optimisation to scalar one.
!!gm ==>>> test in the loop instead of use of mask arrays
!!gm and many acces in memory
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) !== R=zrau = (alpha / beta) (dk[t] / dk[s]) ==!
zrw = ( gdepw(ji,jj,jk ,Kmm) - gdept(ji,jj,jk,Kmm) ) &
!!gm please, use e3w at Kmm below
& / ( gdept(ji,jj,jk-1,Kmm) - gdept(ji,jj,jk,Kmm) )
!
zaw = ( rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem) * zrw ) &
& * tmask(ji,jj,jk) * tmask(ji,jj,jk-1)
zbw = ( rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal) * zrw ) &
& * tmask(ji,jj,jk) * tmask(ji,jj,jk-1)
!
zdt = zaw * ( ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) )
zds = zbw * ( ts(ji,jj,jk-1,jp_sal,Kmm) - ts(ji,jj,jk,jp_sal,Kmm) )
IF( ABS( zds) <= 1.e-20_wp ) zds = 1.e-20_wp
zrau(ji,jj) = MAX( 1.e-20, zdt / zds ) ! only retains positive value of zrau
END_2D
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) !== indicators ==!
! stability indicator: msks=1 if rn2>0; 0 elsewhere
IF( rn2(ji,jj,jk) + 1.e-12 <= 0. ) THEN ; zmsks(ji,jj) = 0._wp
ELSE ; zmsks(ji,jj) = 1._wp * wmask(ji,jj,jk) ! mask so avt and avs masked
ENDIF
! salt fingering indicator: msksf=1 if R>1; 0 elsewhere
IF( zrau(ji,jj) <= 1. ) THEN ; zmskf(ji,jj) = 0._wp
ELSE ; zmskf(ji,jj) = 1._wp
ENDIF
! diffusive layering indicators:
! ! mskdl1=1 if 0< R <1; 0 elsewhere
IF( zrau(ji,jj) >= 1. ) THEN ; zmskd1(ji,jj) = 0._wp
ELSE ; zmskd1(ji,jj) = 1._wp
ENDIF
! ! mskdl2=1 if 0< R <0.5; 0 elsewhere
IF( zrau(ji,jj) >= 0.5 ) THEN ; zmskd2(ji,jj) = 0._wp
ELSE ; zmskd2(ji,jj) = 1._wp
ENDIF
! mskdl3=1 if 0.5< R <1; 0 elsewhere
IF( zrau(ji,jj) <= 0.5 .OR. zrau(ji,jj) >= 1. ) THEN ; zmskd3(ji,jj) = 0._wp
ELSE ; zmskd3(ji,jj) = 1._wp
ENDIF
END_2D
! Update avt and avs
! ------------------
! Constant eddy coefficient: reset to the background value
DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zinr = 1._wp / zrau(ji,jj)
! salt fingering
zrr = zrau(ji,jj) / rn_hsbfr
zrr = zrr * zrr
zavfs = rn_avts / ( 1. + zrr*zrr*zrr ) * zmsks(ji,jj) * zmskf(ji,jj)
zavft = 0.7 * zavfs * zinr
! diffusive layering
zavdt = 1.3635e-6 * EXP( 4.6 * EXP( -0.54*(zinr-1.) ) ) * zmsks(ji,jj) * zmskd1(ji,jj)
zavds = zavdt * zmsks(ji,jj) * ( ( 1.85 * zrau(ji,jj) - 0.85 ) * zmskd3(ji,jj) &
& + 0.15 * zrau(ji,jj) * zmskd2(ji,jj) )
! add to the eddy viscosity coef. previously computed
p_avs(ji,jj,jk) = p_avt(ji,jj,jk) + zavfs + zavds
p_avt(ji,jj,jk) = p_avt(ji,jj,jk) + zavft + zavdt
p_avm(ji,jj,jk) = p_avm(ji,jj,jk) + MAX( zavft + zavdt, zavfs + zavds )
END_2D
! ! ===============
END DO ! End of slab
! ! ===============
!
IF(sn_cfctl%l_prtctl) THEN
sparonuz
committed
CALL prt_ctl(tab3d_1=CASTDP(avt), clinfo1=' ddm - t: ', tab3d_2=CASTDP(avs) , clinfo2=' s: ')
ENDIF
!
END SUBROUTINE zdf_ddm
!!======================================================================
END MODULE zdfddm