Skip to content

Commit 8ca79b2

Browse files
knutarosknutaros
committed
individual coefficients for extra turbulence production term
1 parent 949f7c3 commit 8ca79b2

File tree

10 files changed

+119
-74
lines changed

10 files changed

+119
-74
lines changed

‎src/turbulence/cmue_a.F90‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -85,7 +85,7 @@ subroutine cmue_a(nlev)
8585
!
8686
! !USES:
8787
use turbulence, only: eps
88-
use turbulence, only: P,B,Pb,epsb
88+
use turbulence, only: P,B,Px,Pb,epsb
8989
use turbulence, only: an,as,at,r
9090
use turbulence, only: cmue1,cmue2,gam
9191
use turbulence, only: cm0
@@ -158,7 +158,7 @@ subroutine cmue_a(nlev)
158158

159159
do i=1,nlev-1
160160

161-
Pe = ( P(i) + B(i) )/eps(i)
161+
Pe = ( P(i) + Px(i) + B(i) )/eps(i)
162162
Pbeb = Pb(i)/epsb(i)
163163
r_i = 1./r(i)
164164

‎src/turbulence/dissipationeq.F90‎

Lines changed: 9 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -17,11 +17,15 @@ subroutine dissipationeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
1717
! =
1818
! {\cal D}_\epsilon
1919
! + \frac{\epsilon}{k} ( c_{\epsilon 1} P + c_{\epsilon 3} G
20+
! + c_{\epsilon x} P_x
21+
! + c_{\epsilon 4} P_s
2022
! - c_{\epsilon 2} \epsilon )
2123
! \comma
2224
! \end{equation}
2325
! where $\dot{\epsilon}$ denotes the material derivative of $\epsilon$.
24-
! The production terms $P$ and $G$ follow from \eq{PandG} and
26+
! The production terms $P$ and $G$ follow from \eq{PandG}.
27+
! $P_s$ is Stokes shear production defined in \eq{computePs}
28+
! and $P_x$ accounts for extra turbulence production.
2529
! ${\cal D}_\epsilon$ represents the sum of the viscous and turbulent
2630
! transport terms.
2731
!
@@ -62,9 +66,9 @@ subroutine dissipationeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
6266
! by setting {\tt length\_lim = .true.} in {\tt gotm.yaml}.
6367
!
6468
! !USES:
65-
use turbulence, only: P,B,PSTK,num
69+
use turbulence, only: P,B,Px,PSTK,num
6670
use turbulence, only: tke,tkeo,k_min,eps,eps_min,L
67-
use turbulence, only: ce1,ce2,ce3plus,ce3minus,ce4
71+
use turbulence, only: ce1,ce2,ce3plus,ce3minus,cex,ce4
6872
use turbulence, only: cm0,cde,galp,length_lim
6973
use turbulence, only: epsilon_bc, psi_ubc, psi_lbc, ubc_type, lbc_type
7074
use turbulence, only: sig_e,sig_e0,sig_peps
@@ -122,7 +126,7 @@ subroutine dissipationeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
122126
if (sig_peps) then ! With wave breaking
123127
sig_eff(nlev)=sig_e0
124128
do i=1,nlev-1
125-
peps=(P(i)+B(i))/eps(i)
129+
peps=(P(i)+Px(i)+B(i))/eps(i)
126130
if (peps .gt. 1.) peps=_ONE_
127131
sig_eff(i)=peps*sig_e+(_ONE_-peps)*sig_e0
128132
end do
@@ -147,7 +151,7 @@ subroutine dissipationeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
147151
end if
148152

149153
EpsOverTke = eps(i)/tkeo(i)
150-
prod = ce1*EpsOverTke*P(i) + ce4*EpsOverTke*PSTK(i)
154+
prod = EpsOverTke * ( ce1*P(i) + cex*Px(i) + ce4*PSTK(i) )
151155
buoyan = ce3*EpsOverTke*B(i)
152156
diss = ce2*EpsOverTke*eps(i)
153157

‎src/turbulence/genericeq.F90‎

Lines changed: 7 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -43,12 +43,16 @@ subroutine genericeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
4343
! \label{generic}
4444
! \dot{\psi} = {\cal D}_\psi
4545
! + \frac{\psi}{k} ( c_{\psi_1} P + c_{\psi_3} G
46+
! + c_{\psi x} P_x
47+
! + c_{\psi 4} P_s
4648
! - c_{\psi 2} \epsilon )
4749
! \comma
4850
! \end{equation}
4951
! where $\dot{\psi}$ denotes the material derivative of $\psi$,
5052
! see \cite{UmlaufBurchard2003}.
5153
! The production terms $P$ and $G$ follow from \eq{PandG}.
54+
! $P_s$ is Stokes shear production defined in \eq{computePs}
55+
! and $P_x$ accounts for extra turbulence production.
5256
! ${\cal D}_\psi$ represents the sum of the viscous and turbulent
5357
! transport terms. The rate of dissipation can computed by solving
5458
! \eq{psi_l} for $l$ and inserting the result into \eq{epsilon}.
@@ -118,9 +122,9 @@ subroutine genericeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
118122
! implemented in GOTM and are described in \sect{sec:generate}.
119123
!
120124
! !USES:
121-
use turbulence, only: P,B,PSTK,num
125+
use turbulence, only: P,B,Px,PSTK,num
122126
use turbulence, only: tke,tkeo,k_min,eps,eps_min,L
123-
use turbulence, only: cpsi1,cpsi2,cpsi3plus,cpsi3minus,cpsi4,sig_psi
127+
use turbulence, only: cpsi1,cpsi2,cpsi3plus,cpsi3minus,cpsix,cpsi4,sig_psi
124128
use turbulence, only: gen_m,gen_n,gen_p
125129
use turbulence, only: cm0,cde,galp,length_lim
126130
use turbulence, only: psi_bc, psi_ubc, psi_lbc, ubc_type, lbc_type
@@ -199,7 +203,7 @@ subroutine genericeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
199203

200204
! compute production terms in psi-equation
201205
PsiOverTke = psi(i)/tkeo(i)
202-
prod = cpsi1*PsiOverTke*P(i) + cpsi4*PsiOverTke*PSTK(i)
206+
prod = PsiOverTke * ( cpsi1*P(i) + cpsix*Px(i) + cpsi4*PSTK(i) )
203207
buoyan = cpsi3*PsiOverTke*B(i)
204208
diss = cpsi2*PsiOverTke*eps(i)
205209

‎src/turbulence/lengthscaleeq.F90‎

Lines changed: 8 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -14,11 +14,14 @@ subroutine lengthscaleeq(nlev,dt,depth,u_taus,u_taub,z0s,z0b,h,NN,SS)
1414
! \begin{equation}
1515
! \label{MY}
1616
! \dot{\overline{q^2 l}}
17-
! = {\cal D}_l + l ( E_1 P + E_3 G - E_2 F \epsilon )
17+
! = {\cal D}_l + l ( E_1 P + E_3 G + E_x P_x + E_6 P_s - E_2 F \epsilon )
1818
! \comma
1919
! \end{equation}
2020
! where $\dot{\overline{q^2 l}}$ denotes the material derivative of $q^2 l$.
21-
! The production terms $P$ and $G$ follow from \eq{PandG}, and $\epsilon$
21+
! The production terms $P$ and $G$ follow from \eq{PandG}.
22+
! $P_s$ is Stokes shear production defined in \eq{computePs}
23+
! and $P_x$ accounts for extra turbulence production.
24+
! $\epsilon$
2225
! can be computed either directly from \eq{epsilonMY}, or from \eq{epsilon}
2326
! with the help \eq{B1}.
2427
!
@@ -67,9 +70,9 @@ subroutine lengthscaleeq(nlev,dt,depth,u_taus,u_taub,z0s,z0b,h,NN,SS)
6770
! by setting {\tt length\_lim = .true.} in {\tt gotm.yaml}.
6871
!
6972
! !USES:
70-
use turbulence, only: P,B, PSTK
73+
use turbulence, only: P,B,Px,PSTK
7174
use turbulence, only: tke,tkeo,k_min,eps,eps_min,L
72-
use turbulence, only: kappa,e1,e2,e3,e6,b1
75+
use turbulence, only: kappa,e1,e2,e3,ex,e6,b1
7376
use turbulence, only: MY_length,cm0,cde,galp,length_lim
7477
use turbulence, only: q2l_bc, psi_ubc, psi_lbc, ubc_type, lbc_type
7578
use turbulence, only: sl_var
@@ -159,7 +162,7 @@ subroutine lengthscaleeq(nlev,dt,depth,u_taus,u_taub,z0s,z0b,h,NN,SS)
159162
avh(i) = sl_var(i) * sqrt(2.*tkeo(i))*L(i)
160163

161164
! compute production terms in q^2 l - equation
162-
prod = e1*L(i)*P(i) + e6*L(i)*PSTK(i)
165+
prod = L(i) * ( e1*P(i) + ex*Px(i) + e6*PSTK(i) )
163166
buoyan = e3*L(i)*B(i)
164167
diss = q3(i)/b1*(1.+e2*(L(i)/Lz(i))*(L(i)/Lz(i)))
165168

‎src/turbulence/omegaeq.F90‎

Lines changed: 8 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -17,11 +17,15 @@ subroutine omegaeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
1717
! =
1818
! {\cal D}_\omega
1919
! + \frac{\omega}{k} ( c_{\omega 1} P + c_{\omega 3} G
20+
! + c_{\omega x} P_x
21+
! + c_{\omega 4} P_s
2022
! - c_{\omega 2} \varepsilon )
2123
! \comma
2224
! \end{equation}
2325
! where $\dot{\omega}$ denotes the material derivative of $\omega$.
24-
! The production terms $P$ and $G$ follow from \eq{PandG} and
26+
! The production terms $P$ and $G$ follow from \eq{PandG}.
27+
! $P_s$ is Stokes shear production defined in \eq{computePs}
28+
! and $P_x$ accounts for extra turbulence production.
2529
! ${\cal D}_\omega$ represents the sum of the viscous and turbulent
2630
! transport terms.
2731
!
@@ -57,9 +61,9 @@ subroutine omegaeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
5761
! by setting {\tt length\_lim = .true.} in {\tt gotm.yaml}.
5862
!
5963
! !USES:
60-
use turbulence, only: P,B,PSTK,num
64+
use turbulence, only: P,B,Px,PSTK,num
6165
use turbulence, only: tke,tkeo,k_min,eps,eps_min,L
62-
use turbulence, only: cw1,cw2,cw3plus,cw3minus,cw4
66+
use turbulence, only: cw1,cw2,cw3plus,cw3minus,cwx,cw4
6367
use turbulence, only: cm0,cde,galp,length_lim
6468
use turbulence, only: omega_bc, psi_ubc, psi_lbc, ubc_type, lbc_type
6569
use turbulence, only: sig_w
@@ -130,7 +134,7 @@ subroutine omegaeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
130134
end if
131135

132136
OmgOverTke = omega(i)/tkeo(i)
133-
prod = cw1*OmgOverTke*P(i) + cw4*OmgOverTke*PSTK(i)
137+
prod = OmgOverTke * ( cw1*P(i) + cwx*Px(i) + cw4*PSTK(i) )
134138
buoyan = cw3*OmgOverTke*B(i)
135139
diss = cw2*OmgOverTke*eps(i)
136140

‎src/turbulence/production.F90‎

Lines changed: 26 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -9,20 +9,32 @@ subroutine production(nlev,NN,SS,xP, SSCSTK, SSSTK)
99
!
1010
! !DESCRIPTION:
1111
! This subroutine calculates the production terms of turbulent kinetic
12-
! energy as defined in \eq{PandG} and the production of buoayancy
12+
! energy as defined in \eq{PandG} and the production of buoyancy
1313
! variance as defined in \eq{Pbvertical}.
14-
! The shear-production is computed according to
14+
! The Eulerian shear-production is computed according to
1515
! \begin{equation}
1616
! \label{computeP}
17-
! P = \nu_t (M^2 + \alpha_w N^2) + X_P
17+
! P = \nu_t (M^2 + \alpha_w N^2) + \nu^S_t S_c^2
1818
! \comma
1919
! \end{equation}
2020
! with the turbulent diffusivity of momentum, $\nu_t$, defined in
2121
! \eq{nu}. The shear-frequency, $M$, is discretised as described
2222
! in \sect{sec:shear}.
2323
! The term multiplied by $\alpha_w$ traces back to
2424
! a parameterisation of breaking internal waves suggested by
25-
! \cite{Mellor89}. $X_P$ is an extra production term, connected for
25+
! \cite{Mellor89}.
26+
! The turbulent momentum fluxes due to Stokes velocities induce the
27+
! Stokes-Eulerian cross-shear term
28+
! $S_c^2 = \frac{\partial u}{\partial z}\frac{\partial u_s}{\partial z} + \frac{\partial v}{\partial z}\frac{\partial v_s}{\partial z}$
29+
! with corresponding diffusivity $\nu^S_t$, and the additional
30+
! Stokes shear-production
31+
! \begin{equation}
32+
! \label{computePs}
33+
! P_s = \nu_t S_c^2 + \nu^S_t S_s^2
34+
! \end{equation}
35+
! with squared Stokes shear
36+
! $S_s^2 = \frac{\partial u_s}{\partial z}^2 + \frac{\partial v_s}{\partial z}^2$.
37+
! $X_P$ is an extra production term, connected for
2638
! example with turbulence production caused by sea-grass, see
2739
! \eq{sgProduction} in \sect{sec:seagrass}. {\tt xP} is an {\tt optional}
2840
! argument in the FORTRAN code.
@@ -51,7 +63,7 @@ subroutine production(nlev,NN,SS,xP, SSCSTK, SSSTK)
5163
! kinetic and potential energy in \eq{tkeA} and \eq{kbeq}, respectively.
5264
!
5365
! !USES:
54-
use turbulence, only: P,B,Pb, PSTK
66+
use turbulence, only: P,B,Pb,Px,PSTK
5567
use turbulence, only: num,nuh, nucl
5668
use turbulence, only: alpha,iw_model
5769
IMPLICIT NONE
@@ -92,19 +104,17 @@ subroutine production(nlev,NN,SS,xP, SSCSTK, SSSTK)
92104
alpha_eff=alpha
93105
end if
94106

107+
do i=0,nlev
108+
P(i) = num(i)*( SS(i)+alpha_eff*NN(i) )
109+
B(i) = -nuh(i)*NN(i)
110+
Pb(i) = - B(i)*NN(i)
111+
end do
112+
95113
if ( PRESENT(xP) ) then
96114
do i=0,nlev
97-
P(i) = num(i)*( SS(i)+alpha_eff*NN(i) ) + xP(i)
98-
B(i) = -nuh(i)*NN(i)
99-
Pb(i) = - B(i)*NN(i)
100-
enddo
101-
else
102-
do i=0,nlev
103-
P(i) = num(i)*( SS(i)+alpha_eff*NN(i) )
104-
B(i) = -nuh(i)*NN(i)
105-
Pb(i) = - B(i)*NN(i)
106-
enddo
107-
endif
115+
Px(i) = xP(i)
116+
end do
117+
end if
108118

109119
! P is -<u'w'>du/dz for q2l production with e1,
110120
! PSTK is -<u'w'>dus/dz for q2l prodcution with e6,

‎src/turbulence/q2over2eq.F90‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -13,7 +13,7 @@ subroutine q2over2eq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
1313
! \label{tkeB}
1414
! \dot{\overline{q^2/2}}
1515
! =
16-
! {\cal D}_q + P + G - \epsilon
16+
! {\cal D}_q + P + G + P_x + P_s - \epsilon
1717
! \comma
1818
! \end{equation}
1919
! where $\dot{\overline{q^2/2}}$ denotes the material derivative of $q^2/2$.
@@ -44,7 +44,7 @@ subroutine q2over2eq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
4444
! \end{equation}
4545
!
4646
! !USES:
47-
use turbulence, only: P,B, PSTK
47+
use turbulence, only: P,B,Px,PSTK
4848
use turbulence, only: tke,tkeo,k_min,eps,L
4949
use turbulence, only: q2over2_bc, k_ubc, k_lbc, ubc_type, lbc_type
5050
use turbulence, only: sq_var
@@ -103,7 +103,7 @@ subroutine q2over2eq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
103103
avh(i) = sq_var(i) * sqrt( 2.*tke(i) )*L(i)
104104

105105
! compute production terms in q^2/2-equation
106-
prod = P(i) + PSTK(i)
106+
prod = P(i) + Px(i) + PSTK(i)
107107
buoyan = B(i)
108108
diss = eps(i)
109109

‎src/turbulence/tkeeq.F90‎

Lines changed: 7 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -16,12 +16,15 @@ subroutine tkeeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
1616
! \label{tkeA}
1717
! \dot{k}
1818
! =
19-
! {\cal D}_k + P + G - \epsilon
19+
! {\cal D}_k + P + G + P_x + P_s - \epsilon
2020
! \comma
2121
! \end{equation}
2222
! where $\dot{k}$ denotes the material derivative of $k$. $P$ and $G$ are
2323
! the production of $k$ by mean shear and buoyancy, respectively, and
24-
! $\epsilon$ the rate of dissipation. ${\cal D}_k$ represents the sum of
24+
! $\epsilon$ the rate of dissipation.
25+
! $P_s$ is Stokes shear production defined in \eq{computePs}
26+
! and $P_x$ accounts for extra turbulence production.
27+
! ${\cal D}_k$ represents the sum of
2528
! the viscous and turbulent transport terms.
2629
! For horizontally homogeneous flows, the transport term ${\cal D}_k$
2730
! appearing in \eq{tkeA} is presently expressed by a simple
@@ -58,7 +61,7 @@ subroutine tkeeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
5861
! where $c_\mu^0$ is a constant of the model.
5962
!
6063
! !USES:
61-
use turbulence, only: P,B,PSTK,num
64+
use turbulence, only: P,B,Px,PSTK,num
6265
use turbulence, only: tke,tkeo,k_min,eps
6366
use turbulence, only: k_bc, k_ubc, k_lbc, ubc_type, lbc_type
6467
use turbulence, only: sig_k
@@ -118,7 +121,7 @@ subroutine tkeeq(nlev,dt,u_taus,u_taub,z0s,z0b,h,NN,SS)
118121
avh(i) = num(i)/sig_k
119122

120123
! compute production terms in k-equation
121-
prod = P(i) + PSTK(i)
124+
prod = P(i) + Px(i) + PSTK(i)
122125
buoyan = B(i)
123126
diss = eps(i)
124127

0 commit comments

Comments
 (0)