Context Navigation

zdfmfc.F90 @ 15540

Last change on this file since 15540 was 15540, checked in by sparonuz, 3 years ago
Mixed precision version, tested up to 30 years on ORCA2.
File size: 22.9 KB

Line
1	MODULE zdfmfc
2	!!======================================================================
3	!! * MODULE zdfmfc *
4	!! Ocean physics: Mass-Flux scheme parameterization of Convection:
5	!! Non-local transport for the convective ocean boundary
6	!! layer. Subgrid-scale large eddies are represented by a
7	!! mass-flux contribution (ln_zdfmfc = .TRUE.)
8	!!======================================================================
9	!! History : NEMO !
10	!! 3.6 ! 2016-06 (H. Giordani, R. Bourdallé-Badie) Original code
11	!! 4.2 ! 2020-12 (H. Giordani, R. Bourdallé-Badie) adapt to NEM04.2
12	!!----------------------------------------------------------------------
13	!!----------------------------------------------------------------------
14	!! tra_mfc : Compute the Mass Flux and trends of T/S
15	!! diag_mfc : Modify diagonal of trazdf Matrix
16	!! rhs_mfc : Modify RHS of trazdf Matrix
17	!! zdf_mfc_init : initialization, namelist read, and parameters control
18	!!----------------------------------------------------------------------
19	!
20	USE oce ! ocean dynamics and active tracers
21	USE dom_oce ! ocean space and time domain
22	USE domvvl ! ocean space and time domain : variable volume layer
23	USE domzgr
24	USE zdf_oce ! ocean vertical physics
25	USE sbc_oce ! surface boundary condition: ocean
26	USE phycst ! physical constants
27	USE eosbn2 ! equation of state (eos routine)
28	USE zdfmxl ! mixed layer
29	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
30	USE lib_mpp ! MPP manager
31	USE prtctl ! Print control
32	USE in_out_manager ! I/O manager
33	USE iom ! I/O manager library
34	USE timing ! Timing
35	USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
36
37	IMPLICIT NONE
38	PRIVATE
39
40	PUBLIC tra_mfc ! routine called in step module
41	PUBLIC diag_mfc ! routine called in trazdf module
42	PUBLIC rhs_mfc ! routine called in trazdf module
43	PUBLIC zdf_mfc_init ! routine called in nemo module
44	!
45	REAL(dp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: edmfa, edmfb, edmfc !: diagonal term of the matrix.
46	REAL(dp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: edmftra !: y term for matrix inversion
47	REAL(dp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: edmfm !: y term for matrix inversion
48	!
49	!! Namelist namzdf_edmf
50	REAL(dp) :: rn_cemf ! entrain of T/S
51	REAL(dp) :: rn_cwmf ! detrain of T/S
52	REAL(dp) :: rn_cent ! entrain of the convective mass flux
53	REAL(dp) :: rn_cdet ! detrain of the convective mass flux
54	REAL(dp) :: rn_cap ! Factor of computation for convective area (negative => area constant)
55	REAL(dp) :: App_max ! Maximum of the convective area
56	LOGICAL, PUBLIC, SAVE :: ln_edmfuv !: EDMF flag for velocity !
57	!
58	!! * Substitutions
59	# include "do_loop_substitute.h90"
60	# include "domzgr_substitute.h90"
61	# include "single_precision_substitute.h90"
62	!!----------------------------------------------------------------------
63	!! NEMO/OCE 4.2 , NEMO Consortium (2018)
64	!! $Id: zdfmfc.F90 13783 2020-20-02 15:30:22Z rbourdal $
65	!! Software governed by the CeCILL license (see ./LICENSE)
66	!!----------------------------------------------------------------------
67	CONTAINS
68
69	INTEGER FUNCTION zdf_mfc_alloc()
70	!!----------------------------------------------------------------------
71	!! * FUNCTION zdf_edmf_alloc *
72	!!----------------------------------------------------------------------
73	ALLOCATE( edmfa(jpi,jpj,jpk) , edmfb(jpi,jpj,jpk) , edmfc(jpi,jpj,jpk) &
74	& , edmftra(jpi,jpj,jpk,2), edmfm(jpi,jpj,jpk) , STAT= zdf_mfc_alloc )
75	!
76	IF( lk_mpp ) CALL mpp_sum ( 'zdfmfc', zdf_mfc_alloc )
77	IF( zdf_mfc_alloc /= 0 ) CALL ctl_warn('zdf_mfc_alloc: failed to allocate arrays')
78	END FUNCTION zdf_mfc_alloc
79
80
81	SUBROUTINE tra_mfc( kt, Kmm, pts, Krhs )
82	!!----------------------------------------------------------------------
83	!! * ROUTINE zdf_mfc *
84	!!
85	!! ** Purpose : Compute a mass flux, depending on surface flux, over
86	!! the instable part of the water column.
87	!!
88	!! ** Method : Compute surface instability and mix tracers until stable level
89	!!
90	!!
91	!! ** Action : Compute convection plume and (ta,sa)-trends for trazdf (EDMF scheme)
92	!!
93	!! References :
94	!! Giordani, Bourdallé-Badie and Madec JAMES 2020
95	!!----------------------------------------------------------------------
96	!!----------------------------------------------------------------------
97	INTEGER , INTENT(in) :: Kmm, Krhs ! time level indices
98	REAL(dp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation
99	REAL(dp), DIMENSION(A2D(nn_hls),jpk,2) :: ztsp ! T/S of the plume
100	REAL(dp), DIMENSION(A2D(nn_hls),jpk,2) :: ztse ! T/S at W point
101	REAL(dp), DIMENSION(A2D(nn_hls),jpk) :: zrwp !
102	REAL(dp), DIMENSION(A2D(nn_hls),jpk) :: zrwp2 !
103	REAL(dp), DIMENSION(A2D(nn_hls),jpk) :: zapp !
104	REAL(dp), DIMENSION(A2D(nn_hls),jpk) :: zedmf !
105	REAL(dp), DIMENSION(A2D(nn_hls),jpk) :: zepsT, zepsW !
106	!
107	REAL(dp), DIMENSION(A2D(nn_hls)) :: zustar, zustar2 !
108	REAL(dp), DIMENSION(A2D(nn_hls)) :: zuws, zvws, zsws, zfnet !
109	REAL(dp), DIMENSION(A2D(nn_hls)) :: zfbuo, zrautbm1
110	REAL(dp), DIMENSION(A2D(nn_hls)) :: zrautb, zraupl
111	REAL(dp), DIMENSION(A2D(nn_hls)) :: zwpsurf !
112	REAL(dp), DIMENSION(A2D(nn_hls)) :: zop0 , zsp0 !
113	REAL(dp), DIMENSION(A2D(nn_hls)) :: zrwp_0, zrwp2_0 !
114	REAL(dp), DIMENSION(A2D(nn_hls)) :: zapp0 !
115	REAL(dp), DIMENSION(A2D(nn_hls)) :: zphp, zph, zphpm1, zphm1, zNHydro
116	REAL(dp), DIMENSION(A2D(nn_hls)) :: zhcmo !
117	!
118	REAL(dp), DIMENSION(A2D(nn_hls),jpk) :: zn2 ! N^2
119	REAL(dp), DIMENSION(A2D(nn_hls),2 ) :: zab, zabm1, zabp ! alpha and beta
120
121	REAL(dp), PARAMETER :: zepsilon = 1.e-30 ! local small value
122
123	REAL(dp) :: zrho, zrhop
124	REAL(dp) :: zcnh, znum, zden, zcoef1, zcoef2
125	REAL(dp) :: zca, zcb, zcd, zrw, zxl, zcdet, zctre
126	REAL(dp) :: zaw, zbw, zxw
127	REAL(dp) :: alpha
128	!
129	INTEGER, INTENT(in ) :: kt ! ocean time-step index !
130	!
131	INTEGER :: ji, jj, jk ! dummy loop arguments
132	!
133	!------------------------------------------------------------------
134	! Initialisation of coefficients
135	!------------------------------------------------------------------
136	zca = 1._wp
137	zcb = 1._wp
138	zcd = 1._wp
139
140	DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
141	!------------------------------------------------------------------
142	! Surface boundary condition
143	!------------------------------------------------------------------
144	! surface Stress
145	!--------------------
146	zuws(ji,jj) = utau(ji,jj) * r1_rho0
147	zvws(ji,jj) = vtau(ji,jj) * r1_rho0
148	zustar2(ji,jj) = SQRT(zuws(ji,jj)zuws(ji,jj)+zvws(ji,jj)zvws(ji,jj))
149	zustar(ji,jj) = SQRT(zustar2(ji,jj))
150
151	! Heat Flux
152	!--------------------
153	zfnet(ji,jj) = qns(ji,jj) + qsr(ji,jj)
154	zfnet(ji,jj) = zfnet(ji,jj) / (rho0 * rcp)
155
156	! Water Flux
157	!---------------------
158	zsws(ji,jj) = emp(ji,jj)
159
160	!-------------------------------------------
161	! Initialisation of prognostic variables
162	!-------------------------------------------
163	zrwp (ji,jj,:) = 0._wp ; zrwp2(ji,jj,:) = 0._wp ; zedmf(ji,jj,:) = 0._wp
164	zph (ji,jj) = 0._wp ; zphm1(ji,jj) = 0._wp ; zphpm1(ji,jj) = 0._wp
165	ztsp(ji,jj,:,:)= 0._wp
166
167	! Tracers inside plume (ztsp) and environment (ztse)
168	ztsp(ji,jj,1,jp_tem) = pts(ji,jj,1,jp_tem,Kmm) * tmask(ji,jj,1)
169	ztsp(ji,jj,1,jp_sal) = pts(ji,jj,1,jp_sal,Kmm) * tmask(ji,jj,1)
170	ztse(ji,jj,1,jp_tem) = pts(ji,jj,1,jp_tem,Kmm) * tmask(ji,jj,1)
171	ztse(ji,jj,1,jp_sal) = pts(ji,jj,1,jp_sal,Kmm) * tmask(ji,jj,1)
172	END_2D
173
174	CALL eos( ztse(:,:,1,:) , zrautb(:,:) )
175	CALL eos( ztsp(:,:,1,:) , zraupl(:,:) )
176
177	!-------------------------------------------
178	! Boundary Condition of Mass Flux (plume velo.; convective area, entrain/detrain)
179	!-------------------------------------------
180	zhcmo(:,:) = e3t(A1Di(nn_hls),A1Dj(nn_hls),1,Kmm)
181	zfbuo(:,:) = 0._wp
182	WHERE ( ABS(zrautb(:,:)) > 1.e-20 ) zfbuo(:,:) = &
183	& grav * ( 2.e-4_wp *zfnet(:,:) &
184	& - 7.6E-4_wp*pts(A2D(nn_hls),1,jp_sal,Kmm) &
185	& * zsws(:,:)/zrautb(:,:)) * zhcmo(:,:)
186
187	zedmf(:,:,1) = -0.065_wp(ABS(zfbuo(:,:)))(1._wp/3._wp)SIGN(1._wp,zfbuo(:,:))
188	zedmf(:,:,1) = MAX(0., zedmf(:,:,1))
189
190	zwpsurf(:,:) = 2._wp/3._wpzustar(:,:) + 2._wp/3._wpABS(zfbuo(:,:))**(1._wp/3._wp)
191	zwpsurf(:,:) = MAX(1.e-5_wp,zwpsurf(:,:))
192	zwpsurf(:,:) = MIN(1.,zwpsurf(:,:))
193
194	zapp(:,:,:) = App_max
195	WHERE(zwpsurf .NE. 0.) zapp(:,:,1) = MIN(MAX(0.,zedmf(:,:,1)/zwpsurf(:,:)), App_max)
196
197	zedmf(:,:,1) = 0._wp
198	zrwp (:,:,1) = 0._wp
199	zrwp2(:,:,1) = 0._wp
200	zepsT(:,:,:) = 0.001_wp
201	zepsW(:,:,:) = 0.001_wp
202
203
204	!--------------------------------------------------------------
205	! Compute plume properties
206	! In the same loop on vert. levels computation of:
207	! - Vertical velocity: zWp
208	! - Convective Area: zAp
209	! - Tracers properties inside the plume (if necessary): ztp
210	!---------------------------------------------------------------
211
212	DO jk= 2, jpk
213
214	! Compute the buoyancy acceleration on T-points at jk-1
215	zrautbm1(:,:) = zrautb(:,:)
216	CALL eos( pts (:,:,jk ,:,Kmm) , zrautb(:,:) )
217	CALL eos( ztsp(:,:,jk-1,: ) , zraupl(:,:) )
218
219	DO_2D( 0, 0, 0, 0 )
220	zphm1(ji,jj) = zphm1(ji,jj) + grav * zrautbm1(ji,jj) * e3t(ji,jj,jk-1, Kmm)
221	zphpm1(ji,jj) = zphpm1(ji,jj) + grav * zraupl(ji,jj) * e3t(ji,jj,jk-1, Kmm)
222	zph(ji,jj) = zphm1(ji,jj) + grav * zrautb(ji,jj) * e3t(ji,jj,jk , Kmm)
223	zph(ji,jj) = MAX( zph(ji,jj), zepsilon)
224	END_2D
225
226	WHERE(zrautbm1 .NE. 0.) zfbuo(:,:) = grav * (zraupl(:,:) - zrautbm1(:,:)) / zrautbm1(:,:)
227
228	DO_2D( 0, 0, 0, 0 )
229
230	! Compute Environment of Plume. Interpolation T/S (before time step) on W-points
231	zrw = (gdept(ji,jj,jk,Kmm) - gdepw(ji,jj,jk,Kmm)) &
232	& / (gdept(ji,jj,jk,Kmm) - gdept(ji,jj,jk-1,Kmm))
233	ztse(ji,jj,jk,:) = (pts(ji,jj,jk,:,Kmm) * zrw + pts(ji,jj,jk-1,:,Kmm)(1._wp - zrw) )tmask(ji,jj,jk)
234
235	!---------------------------------------------------------------
236	! Compute the vertical velocity on W-points
237	!---------------------------------------------------------------
238
239	! Non-hydrostatic pressure terms in the wp2 equation
240	zcnh = 0.2_wp
241	znum = 0.5_wp + zcnh - &
242	(zcnhgravzraupl(ji,jj)/zph(ji,jj)+zcb*zepsW(ji,jj,jk-1)) &
243	e3t(ji,jj,jk-1,Kmm)0.5_wp
244	zden = 0.5_wp + zcnh + &
245	(zcnhgravzraupl(ji,jj)/zph(ji,jj)+zcb*zepsW(ji,jj,jk-1)) &
246	e3t(ji,jj,jk-1,Kmm)0.5_wp
247
248	zcoef1 = zca*e3t(ji,jj,jk-1,Kmm) / zden
249	zcoef2 = znum/zden
250
251	! compute wp2
252	zrwp2(ji,jj,jk) = zcoef1*zfbuo(ji,jj) &
253	+ zcoef2*zrwp2(ji,jj,jk-1)
254	zrwp2(ji,jj,jk) = MAX ( zrwp2(ji,jj,jk)*wmask(ji,jj,jk) , 0.)
255	zrwp (ji,jj,jk) = SQRT( zrwp2(ji,jj,jk) )
256
257	!----------------------------------------------------------------------------------
258	! Compute convective area on W-point
259	! Compute vertical profil of the convective area with mass conservation hypothesis
260	! If rn_cap negative => constant value on the water column.
261	!----------------------------------------------------------------------------------
262	IF( rn_cap .GT. 0. ) THEN
263
264	zxw = MAX(zrwp(ji,jj,jk-1), zrwp(ji,jj,jk) )
265	IF( zxw > 0. ) THEN
266
267	zxl = (zrwp(ji,jj,jk-1)-zrwp(ji,jj,jk))/(e3t(ji,jj,jk-1,Kmm)*zxw)
268	IF (zxl .LT. 0._wp) THEN
269	zctre = -1.rn_capzxl
270	zcdet = 0._wp
271	ELSE
272	zctre = 0._wp
273	zcdet = rn_cap*zxl
274	END IF
275	zapp(ji,jj,jk) = zapp(ji,jj,jk-1)* &
276	& (1._wp + (zxl + zctre - zcdet )*e3t(ji,jj,jk-1,Kmm))
277	ELSE
278	zapp(ji,jj,jk) = App_max
279	END IF
280	zapp(ji,jj,jk) = MIN( MAX(zapp(ji,jj,jk),0.), App_max)
281	ELSE
282	zapp(ji,jj,jk) = -1. * rn_cap
283	END IF
284
285	! Compute Mass Flux on W-point
286	zedmf(ji,jj,jk) = -zapp(ji,jj,jk) * zrwp(ji,jj,jk)* wmask(ji,jj,jk)
287
288	! Compute Entrainment coefficient
289	IF(rn_cemf .GT. 0.) THEN
290	zxw = 0.5_wp*(zrwp(ji,jj,jk-1)+ zrwp(ji,jj,jk) )
291	zepsT(ji,jj,jk) = 0.01_wp
292	IF( zxw > 0. ) THEN
293	zepsT(ji,jj,jk) = zepsT(ji,jj,jk) + &
294	& ABS( zrwp(ji,jj,jk-1)-zrwp(ji,jj,jk) ) &
295	& / ( e3t(ji,jj,jk-1,Kmm) * zxw )
296	zepsT(ji,jj,jk) = zepsT(ji,jj,jk) * rn_cemf * wmask(ji,jj,jk)
297	ENDIF
298	ELSE
299	zepsT(ji,jj,jk) = -rn_cemf
300	ENDIF
301
302	! Compute the detrend coef for velocity (on W-point and not T-points, bug ???)
303	IF(rn_cwmf .GT. 0.) THEN
304	zepsW(ji,jj,jk) = rn_cwmf * zepsT(ji,jj,jk)
305	ELSE
306	zepsW(ji,jj,jk) = -rn_cwmf
307	ENDIF
308
309	!---------------------------------------------------------------
310	! Compute the plume properties on T-points
311	!---------------------------------------------------------------
312	IF(zrwp (ji,jj,jk) .LT. 1.e-12_wp .AND. zrwp (ji,jj,jk-1) .LT. 1.e-12_wp) THEN
313	ztsp(ji,jj,jk-1,jp_tem) = pts(ji,jj,jk-1,jp_tem,Kmm)
314	ztsp(ji,jj,jk-1,jp_sal) = pts(ji,jj,jk-1,jp_sal,Kmm)
315	ENDIF
316
317	zcoef1 = (1._wp-zepsT(ji,jj,jk)(1._wp-zrw)e3w(ji,jj,jk,Kmm)*wmask(ji,jj,jk ) ) &
318	& / (1._wp+zepsT(ji,jj,jk)zrwe3w(ji,jj,jk,Kmm)*wmask(ji,jj,jk) )
319	!
320	zcoef2 = zepsT(ji,jj,jk)e3w(ji,jj,jk,Kmm)wmask(ji,jj,jk) &
321	& / (1._wp+zepsT(ji,jj,jk)zrwe3w(ji,jj,jk,Kmm)*wmask(ji,jj,jk))
322	!
323	ztsp(ji,jj,jk,jp_tem) = (zcoef1 * ztsp(ji,jj,jk-1,jp_tem) + &
324	& zcoef2 * ztse(ji,jj,jk ,jp_tem) )*tmask(ji,jj,jk)
325	ztsp(ji,jj,jk,jp_sal) = (zcoef1 * ztsp(ji,jj,jk-1,jp_sal) + &
326	& zcoef2 * ztse(ji,jj,jk ,jp_sal) )*tmask(ji,jj,jk)
327
328	END_2D
329	END DO ! end of loop on jpk
330
331	! Compute Mass Flux on T-point
332	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
333	edmfm(ji,jj,jk) = (zedmf(ji,jj,jk+1) + zedmf(ji,jj,jk) )*0.5_wp
334	END_3D
335	DO_2D( 0, 0, 0, 0 )
336	edmfm(ji,jj,jpk) = zedmf(ji,jj,jpk)
337	END_2D
338
339	! Save variable (on T point)
340	CALL iom_put( "mf_Tp" , ztsp(:,:,:,jp_tem) ) ! Save plume temperature
341	CALL iom_put( "mf_Sp" , ztsp(:,:,:,jp_sal) ) ! Save plume salinity
342	CALL iom_put( "mf_mf" , edmfm(:,:,:) ) ! Save Mass Flux
343	! Save variable (on W point)
344	CALL iom_put( "mf_wp" , zrwp (:,:,:) ) ! Save convective velocity in the plume
345	CALL iom_put( "mf_app", zapp (:,:,:) ) ! Save convective area
346
347	!=================================================================================
348	! Computation of a tridiagonal matrix and right hand side terms of the linear system
349	!=================================================================================
350	DO_3D( 0, 0, 0, 0, 1, jpk )
351	edmfa(ji,jj,jk) = 0._wp
352	edmfb(ji,jj,jk) = 0._wp
353	edmfc(ji,jj,jk) = 0._wp
354	edmftra(ji,jj,jk,:) = 0._wp
355	END_3D
356
357	!---------------------------------------------------------------
358	! Diagonal terms
359	!---------------------------------------------------------------
360	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
361	edmfa(ji,jj,jk) = 0._wp
362	edmfb(ji,jj,jk) = -edmfm(ji,jj,jk ) / e3w(ji,jj,jk+1,Kmm)
363	edmfc(ji,jj,jk) = edmfm(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm)
364	END_3D
365	DO_2D( 0, 0, 0, 0 )
366	edmfa(ji,jj,jpk) = -edmfm(ji,jj,jpk-1) / e3w(ji,jj,jpk,Kmm)
367	edmfb(ji,jj,jpk) = edmfm(ji,jj,jpk ) / e3w(ji,jj,jpk,Kmm)
368	edmfc(ji,jj,jpk) = 0._wp
369	END_2D
370
371	!---------------------------------------------------------------
372	! right hand side term for Temperature
373	!---------------------------------------------------------------
374	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
375	edmftra(ji,jj,jk,1) = - edmfm(ji,jj,jk ) * ztsp(ji,jj,jk ,jp_tem) / e3w(ji,jj,jk+1,Kmm) &
376	& + edmfm(ji,jj,jk+1) * ztsp(ji,jj,jk+1,jp_tem) / e3w(ji,jj,jk+1,Kmm)
377	END_3D
378	DO_2D( 0, 0, 0, 0 )
379	edmftra(ji,jj,jpk,1) = - edmfm(ji,jj,jpk-1) * ztsp(ji,jj,jpk-1,jp_tem) / e3w(ji,jj,jpk,Kmm) &
380	& + edmfm(ji,jj,jpk ) * ztsp(ji,jj,jpk ,jp_tem) / e3w(ji,jj,jpk,Kmm)
381	END_2D
382
383	!---------------------------------------------------------------
384	! Right hand side term for Salinity
385	!---------------------------------------------------------------
386	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
387	edmftra(ji,jj,jk,2) = - edmfm(ji,jj,jk ) * ztsp(ji,jj,jk ,jp_sal) / e3w(ji,jj,jk+1,Kmm) &
388	& + edmfm(ji,jj,jk+1) * ztsp(ji,jj,jk+1,jp_sal) / e3w(ji,jj,jk+1,Kmm)
389	END_3D
390	DO_2D( 0, 0, 0, 0 )
391	edmftra(ji,jj,jpk,2) = - edmfm(ji,jj,jpk-1) * ztsp(ji,jj,jpk-1,jp_sal) / e3w(ji,jj,jpk,Kmm) &
392	& + edmfm(ji,jj,jpk ) * ztsp(ji,jj,jpk ,jp_sal) / e3w(ji,jj,jpk,Kmm)
393	END_2D
394	!
395	END SUBROUTINE tra_mfc
396
397
398	SUBROUTINE diag_mfc( zdiagi, zdiagd, zdiags, p2dt, Kaa )
399
400	REAL(dp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: zdiagi, zdiagd, zdiags ! inout: tridaig. terms
401	REAL(dp) , INTENT(in ) :: p2dt ! tracer time-step
402	INTEGER , INTENT(in ) :: Kaa ! ocean time level indices
403
404	INTEGER :: ji, jj, jk ! dummy loop arguments
405
406	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
407	zdiagi(ji,jj,jk) = zdiagi(ji,jj,jk) + e3t(ji,jj,jk,Kaa) * p2dt *edmfa(ji,jj,jk)
408	zdiags(ji,jj,jk) = zdiags(ji,jj,jk) + e3t(ji,jj,jk,Kaa) * p2dt *edmfc(ji,jj,jk)
409	zdiagd(ji,jj,jk) = zdiagd(ji,jj,jk) + e3t(ji,jj,jk,Kaa) * p2dt *edmfb(ji,jj,jk)
410	END_3D
411
412	END SUBROUTINE diag_mfc
413
414	SUBROUTINE rhs_mfc( zrhs, jjn )
415
416	REAL(dp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: zrhs ! inout: rhs trend
417	INTEGER , INTENT(in ) :: jjn ! tracer indices
418
419	INTEGER :: ji, jj, jk ! dummy loop arguments
420
421	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
422	zrhs(ji,jj,jk) = zrhs(ji,jj,jk) + edmftra(ji,jj,jk,jjn)
423	END_3D
424
425	END SUBROUTINE rhs_mfc
426
427
428
429	SUBROUTINE zdf_mfc_init
430	!!----------------------------------------------------------------------
431	!! * ROUTINE zdf_mfc_init *
432	!!
433	!! ** Purpose : Initialization of the vertical eddy diffivity and
434	!! mass flux
435	!!
436	!! ** Method : Read the namzdf_mfc namelist and check the parameters
437	!! called at the first timestep (nit000)
438	!!
439	!! ** input : Namlist namzdf_mfc
440	!!
441	!! ** Action : Increase by 1 the nstop flag is setting problem encounter
442	!!
443	!!----------------------------------------------------------------------
444	!
445	INTEGER :: jk ! dummy loop indices
446	INTEGER :: ios ! Local integer output status for namelist read
447	REAL(dp):: zcr ! local scalar
448	!!
449	NAMELIST/namzdf_mfc/ ln_edmfuv, rn_cemf, rn_cwmf, rn_cent, rn_cdet, rn_cap, App_max
450	!!----------------------------------------------------------
451	!
452	!
453	! REWIND( numnam_ref ) ! Namelist namzdf_mfc in reference namelist : Vertical eddy diffivity mass flux
454	READ ( numnam_ref, namzdf_mfc, IOSTAT = ios, ERR = 901)
455	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_edmf in reference namelist' )
456
457	! REWIND( numnam_cfg ) ! Namelist namzdf_mfc in configuration namelist : Vertical eddy diffivity mass flux
458	READ ( numnam_cfg, namzdf_mfc, IOSTAT = ios, ERR = 902 )
459	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_edmf in configuration namelist' )
460	IF(lwm) WRITE ( numond, namzdf_mfc )
461
462	IF(lwp) THEN !* Control print
463	WRITE(numout,*)
464	WRITE(numout,*) 'zdf_mfc_init'
465	WRITE(numout,*) '~~~~~~~~~~~~~'
466	WRITE(numout,*) ' Namelist namzdf_mfc : set eddy diffusivity Mass Flux Convection'
467	WRITE(numout,*) ' Apply mass flux on velocities (Not yet avail.) ln_edmfuv = ', ln_edmfuv
468	WRITE(numout,*) ' Coeff for entrain/detrain T/S of plume (Neg => cte) rn_cemf = ', rn_cemf
469	WRITE(numout,*) ' Coeff for entrain/detrain Wp of plume (Neg => cte) rn_cwmf = ', rn_cwmf
470	WRITE(numout,*) ' Coeff for entrain/detrain area of plume rn_cap = ', rn_cap
471	WRITE(numout,*) ' Coeff for entrain area of plume rn_cent = ', rn_cent
472	WRITE(numout,*) ' Coeff for detrain area of plume rn_cdet = ', rn_cdet
473	WRITE(numout,*) ' Max convective area App_max = ', App_max
474	ENDIF
475	!* allocate edmf arrays
476	IF( zdf_mfc_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_edmf_init : unable to allocate arrays' )
477	edmfa(:,:,:) = 0._wp
478	edmfb(:,:,:) = 0._wp
479	edmfc(:,:,:) = 0._wp
480	edmftra(:,:,:,:) = 0._wp
481	!
482	END SUBROUTINE zdf_mfc_init
483
484	!!======================================================================
485
486	!!======================================================================
487	END MODULE zdfmfc
488
489
490

Note: See TracBrowser for help on using the repository browser.

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

source: NEMO/branches/2021/dev_r14116_HPC-10_mcastril_Mixed_Precision_implementation/src/OCE/ZDF/zdfmfc.F90 @ 15540

Download in other formats: