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.

zdftmx.F90 in branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

Context Navigation

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90 @ 4401

Last change on this file since 4401 was 3211, checked in by spickles2, 12 years ago

Stephen Pickles, 11 Dec 2011

Commit to bring the rest of the DCSE NEMO development branch
in line with the latest development version. This includes
array index re-ordering of all OPA_SRC/.

Property svn:keywords set to Id

File size: 32.2 KB

Line
1	MODULE zdftmx
2	!!========================================================================
3	!! * MODULE zdftmx *
4	!! Ocean physics: vertical tidal mixing coefficient
5	!!========================================================================
6	!! History : 1.0 ! 2004-04 (L. Bessieres, G. Madec) Original code
7	!! - ! 2006-08 (A. Koch-Larrouy) Indonesian strait
8	!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
9	!!----------------------------------------------------------------------
10	#if defined key_zdftmx \|\| defined key_esopa
11	!!----------------------------------------------------------------------
12	!! 'key_zdftmx' Tidal vertical mixing
13	!!----------------------------------------------------------------------
14	!! zdf_tmx : global momentum & tracer Kz with tidal induced Kz
15	!! tmx_itf : Indonesian momentum & tracer Kz with tidal induced Kz
16	!!----------------------------------------------------------------------
17	USE oce ! ocean dynamics and tracers variables
18	USE dom_oce ! ocean space and time domain variables
19	USE zdf_oce ! ocean vertical physics variables
20	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
21	USE eosbn2 ! ocean equation of state
22	USE phycst ! physical constants
23	USE prtctl ! Print control
24	USE in_out_manager ! I/O manager
25	USE iom ! I/O Manager
26	USE lib_mpp ! MPP library
27	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
28
29	IMPLICIT NONE
30	PRIVATE
31
32	PUBLIC zdf_tmx ! called in step module
33	PUBLIC zdf_tmx_init ! called in opa module
34	PUBLIC zdf_tmx_alloc ! called in nemogcm module
35
36	LOGICAL, PUBLIC, PARAMETER :: lk_zdftmx = .TRUE. !: tidal mixing flag
37
38	! !!* Namelist namzdf_tmx : tidal mixing *
39	REAL(wp) :: rn_htmx = 500. ! vertical decay scale for turbulence (meters)
40	REAL(wp) :: rn_n2min = 1.e-8 ! threshold of the Brunt-Vaisala frequency (s-1)
41	REAL(wp) :: rn_tfe = 1./3. ! tidal dissipation efficiency (St Laurent et al. 2002)
42	REAL(wp) :: rn_me = 0.2 ! mixing efficiency (Osborn 1980)
43	LOGICAL :: ln_tmx_itf = .TRUE. ! Indonesian Through Flow (ITF): Koch-Larrouy et al. (2007) parameterization
44	REAL(wp) :: rn_tfe_itf = 1. ! ITF tidal dissipation efficiency (St Laurent et al. 2002)
45
46	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: en_tmx ! energy available for tidal mixing (W/m2)
47	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: mask_itf ! mask to use over Indonesian area
48	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: az_tmx ! coefficient used to evaluate the tidal induced Kz
49
50	!! * Control permutation of array indices
51	# include "oce_ftrans.h90"
52	# include "dom_oce_ftrans.h90"
53	# include "zdf_oce_ftrans.h90"
54	!FTRANS az_tmx :I :I :z
55
56	!! * Substitutions
57	# include "domzgr_substitute.h90"
58	# include "vectopt_loop_substitute.h90"
59	!!----------------------------------------------------------------------
60	!! NEMO/OPA 4.0 , NEMO Consortium (2011)
61	!! $Id$
62	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
63	!!----------------------------------------------------------------------
64	CONTAINS
65
66	INTEGER FUNCTION zdf_tmx_alloc()
67	!!----------------------------------------------------------------------
68	!! * FUNCTION zdf_tmx_alloc *
69	!!----------------------------------------------------------------------
70	ALLOCATE(en_tmx(jpi,jpj), mask_itf(jpi,jpj), az_tmx(jpi,jpj,jpk), STAT=zdf_tmx_alloc )
71	!
72	IF( lk_mpp ) CALL mpp_sum ( zdf_tmx_alloc )
73	IF( zdf_tmx_alloc /= 0 ) CALL ctl_warn('zdf_tmx_alloc: failed to allocate arrays')
74	END FUNCTION zdf_tmx_alloc
75
76
77	SUBROUTINE zdf_tmx( kt )
78	!!----------------------------------------------------------------------
79	!! * ROUTINE zdf_tmx *
80	!!
81	!! ** Purpose : add to the vertical mixing coefficients the effect of
82	!! tidal mixing (Simmons et al 2004).
83	!!
84	!! ** Method : - tidal-induced vertical mixing is given by:
85	!! Kz_tides = az_tmx / max( rn_n2min, N^2 )
86	!! where az_tmx is a coefficient that specified the 3D space
87	!! distribution of the faction of tidal energy taht is used
88	!! for mixing. Its expression is set in zdf_tmx_init routine,
89	!! following Simmons et al. 2004.
90	!! NB: a specific bounding procedure is performed on av_tide
91	!! so that the input tidal energy is actually almost used. The
92	!! basic maximum value is 60 cm2/s, but values of 300 cm2/s
93	!! can be reached in area where bottom stratification is too
94	!! weak.
95	!!
96	!! - update av_tide in the Indonesian Through Flow area
97	!! following Koch-Larrouy et al. (2007) parameterisation
98	!! (see tmx_itf routine).
99	!!
100	!! - update the model vertical eddy viscosity and diffusivity:
101	!! avt = avt + av_tides
102	!! avm = avm + av_tides
103	!! avmu = avmu + mi(av_tides)
104	!! avmv = avmv + mj(av_tides)
105	!!
106	!! ** Action : avt, avm, avmu, avmv increased by tidal mixing
107	!!
108	!! References : Simmons et al. 2004, Ocean Modelling, 6, 3-4, 245-263.
109	!! Koch-Larrouy et al. 2007, GRL.
110	!!----------------------------------------------------------------------
111	USE oce, zav_tide => ua ! use ua as workspace
112	USE wrk_nemo, ONLY: zkz => wrk_2d_1
113	!!
114	INTEGER, INTENT(in) :: kt ! ocean time-step
115	!!
116	INTEGER :: ji, jj, jk ! dummy loop indices
117	REAL(wp) :: ztpc ! scalar workspace
118	#if defined key_z_first
119	REAL(wp) :: ztpc ! scalar workspace
120	#endif
121	!!----------------------------------------------------------------------
122
123	IF(wrk_in_use(2, 1))THEN
124	CALL ctl_stop('zdf_tmx : requested workspace array unavailable.') ; RETURN
125	END IF
126	! ! ----------------------- !
127	! ! Standard tidal mixing ! (compute zav_tide)
128	! ! ----------------------- !
129	! !* First estimation (with n2 bound by rn_n2min) bounded by 60 cm2/s
130	zav_tide(:,:,:) = MIN( 60.e-4, az_tmx(:,:,:) / MAX( rn_n2min, rn2(:,:,:) ) )
131
132	zkz(:,:) = 0.e0 !* Associated potential energy consummed over the whole water column
133	DO jk = 2, jpkm1
134	zkz(:,:) = zkz(:,:) + fse3w(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zav_tide(:,:,jk)* tmask(:,:,jk)
135	END DO
136
137	DO jj = 1, jpj !* Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx
138	DO ji = 1, jpi
139	IF( zkz(ji,jj) /= 0.e0 ) zkz(ji,jj) = en_tmx(ji,jj) / zkz(ji,jj)
140	END DO
141	END DO
142
143	#if defined key_z_first
144	DO jj = 1, jpj
145	DO ji = 1, jpi
146	zscal = MIN( zkz(ji,jj), 30./6. ) !kz max = 300 cm2/s
147	DO jk = 2, jpkm1 !* Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zav_tide bound by 300 cm2/s
148	zav_tide(ji,jj,jk) = zav_tide(ji,jj,jk) * zscal
149	END DO
150	END DO
151	END DO
152	#else
153	DO jk = 2, jpkm1 !* Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zav_tide bound by 300 cm2/s
154	zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) !kz max = 300 cm2/s
155	END DO
156	#endif
157
158	IF( kt == nit000 ) THEN !* check at first time-step: diagnose the energy consumed by zav_tide
159	ztpc = 0.e0
160	#if defined key_z_first
161	DO jj = 1, jpj
162	DO ji = 1, jpi
163	DO jk = 1, jpk
164	#else
165	DO jk= 1, jpk
166	DO jj= 1, jpj
167	DO ji= 1, jpi
168	#endif
169	ztpc = ztpc + fse3w(ji,jj,jk) * e1t(ji,jj) * e2t(ji,jj) &
170	& * MAX( 0.e0, rn2(ji,jj,jk) ) * zav_tide(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj)
171	END DO
172	END DO
173	END DO
174	ztpc= rau0 / ( rn_tfe * rn_me ) * ztpc
175	IF(lwp) WRITE(numout,*)
176	IF(lwp) WRITE(numout,) ' N Total power consumption by av_tide : ztpc = ', ztpc 1.e-12 ,'TW'
177	ENDIF
178
179	! ! ----------------------- !
180	! ! ITF tidal mixing ! (update zav_tide)
181	! ! ----------------------- !
182	IF( ln_tmx_itf ) CALL tmx_itf( kt, zav_tide )
183
184	! ! ----------------------- !
185	! ! Update mixing coefs !
186	! ! ----------------------- !
187	#if defined key_z_first
188	!* update momentum & tracer diffusivity with tidal mixing
189	DO jj = 1, jpj
190	DO ji = 1, jpi
191	DO jk = 2, jpkm1
192	avt(ji,jj,jk) = avt(ji,jj,jk) + zav_tide(ji,jj,jk)
193	avm(ji,jj,jk) = avm(ji,jj,jk) + zav_tide(ji,jj,jk)
194	END DO
195	END DO
196	END DO
197	DO jj = 2, jpjm1
198	DO ji = 2, fpim1
199	DO jk = 2, jpkm1
200	avmu(ji,jj,jk) = avmu(ji,jj,jk) + 0.5 * ( zav_tide(ji,jj,jk) + zav_tide(ji+1,jj ,jk) ) * umask(ji,jj,jk)
201	avmv(ji,jj,jk) = avmv(ji,jj,jk) + 0.5 * ( zav_tide(ji,jj,jk) + zav_tide(ji ,jj+1,jk) ) * vmask(ji,jj,jk)
202	END DO
203	END DO
204	END DO
205	#else
206	DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with tidal mixing
207	avt(:,:,jk) = avt(:,:,jk) + zav_tide(:,:,jk)
208	avm(:,:,jk) = avm(:,:,jk) + zav_tide(:,:,jk)
209	DO jj = 2, jpjm1
210	DO ji = fs_2, fs_jpim1 ! vector opt.
211	avmu(ji,jj,jk) = avmu(ji,jj,jk) + 0.5 * ( zav_tide(ji,jj,jk) + zav_tide(ji+1,jj ,jk) ) * umask(ji,jj,jk)
212	avmv(ji,jj,jk) = avmv(ji,jj,jk) + 0.5 * ( zav_tide(ji,jj,jk) + zav_tide(ji ,jj+1,jk) ) * vmask(ji,jj,jk)
213	END DO
214	END DO
215	END DO
216	#endif
217	CALL lbc_lnk( avmu, 'U', 1. ) ; CALL lbc_lnk( avmv, 'V', 1. ) ! lateral boundary condition
218
219	! !* output tidal mixing coefficient
220	CALL iom_put( "av_tide", zav_tide )
221
222	IF(ln_ctl) CALL prt_ctl(tab3d_1=zav_tide , clinfo1=' tmx - av_tide: ', tab3d_2=avt, clinfo2=' avt: ', ovlap=1, kdim=jpk)
223	!
224	IF(wrk_not_released(2, 1))THEN
225	CALL ctl_stop('zdf_tmx : failed to release workspace array.')
226	END IF
227	!
228	END SUBROUTINE zdf_tmx
229
230
231	SUBROUTINE tmx_itf( kt, pav )
232	!!----------------------------------------------------------------------
233	!! * ROUTINE tmx_itf *
234	!!
235	!! ** Purpose : modify the vertical eddy diffusivity coefficients
236	!! (pav) in the Indonesian Through Flow area (ITF).
237	!!
238	!! ** Method : - Following Koch-Larrouy et al. (2007), in the ITF defined
239	!! by msk_itf (read in a file, see tmx_init), the tidal
240	!! mixing coefficient is computed with :
241	!! * q=1 (i.e. all the tidal energy remains trapped in
242	!! the area and thus is used for mixing)
243	!! * the vertical distribution of the tifal energy is a
244	!! proportional to N above the thermocline (d(N^2)/dz > 0)
245	!! and to N^2 below the thermocline (d(N^2)/dz < 0)
246	!!
247	!! ** Action : av_tide updated in the ITF area (msk_itf)
248	!!
249	!! References : Koch-Larrouy et al. 2007, GRL
250	!!----------------------------------------------------------------------
251	USE wrk_nemo, ONLY: zkz => wrk_2d_5
252	USE wrk_nemo, ONLY: zsum1 => wrk_2d_2, zsum2 => wrk_2d_3, zsum => wrk_2d_4
253	USE wrk_nemo, ONLY: zempba_3d_1 => wrk_3d_1, zempba_3d_2 => wrk_3d_2
254	USE wrk_nemo, ONLY: zempba_3d => wrk_3d_3, zdn2dz => wrk_3d_4
255	USE wrk_nemo, ONLY: zavt_itf => wrk_3d_5
256	!! DCSE_NEMO: need additional directives for renamed module variables
257	!FTRANS zempba_3d_1 zempba_3d_2 zempba_3d zdn2dz zavt_itf :I :I :z
258	!!
259	INTEGER , INTENT(in ) :: kt ! ocean time-step
260
261	!! DCSE_NEMO: This style defeats ftrans
262	! REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pav ! Tidal mixing coef.
263	!FTRANS pav :I :I :z
264	REAL(wp), INTENT(inout) :: pav(jpi,jpj,jpk) ! Tidal mixing coef.
265	!!
266	INTEGER :: ji, jj, jk ! dummy loop indices
267	REAL(wp) :: zcoef, ztpc ! temporary scalar
268	!!----------------------------------------------------------------------
269	!
270	IF( wrk_in_use(2, 2,3,4,5) .OR. wrk_in_use(3, 1,2,3,4,5) )THEN
271	CALL ctl_stop('tmx_itf : requested workspace arrays unavailable.')
272	RETURN
273	END IF
274	! ! compute the form function using N2 at each time step
275	#if defined key_z_first
276	DO jj = 1, jpj
277	DO ji = 1, jpi
278	DO jk = 1, jpkm1
279	zdn2dz (ji,jj,jk) = rn2(ji,jj,jk) - rn2(ji,jj,jk+1) ! Vertical profile of dN2/dz
280	zempba_3d_1(ji,jj,jk) = SQRT( MAX( 0.e0, rn2(ji,jj,jk) ) ) ! - - of N
281	zempba_3d_2(ji,jj,jk) = MAX( 0.e0, rn2(ji,jj,jk) ) ! - - of N^2
282	END DO
283	zempba_3d_1(ji,jj,jpk) = 0.e0
284	zempba_3d_2(ji,jj,jpk) = 0.e0
285	END DO
286	END DO
287	#else
288	zempba_3d_1(:,:,jpk) = 0.e0
289	zempba_3d_2(:,:,jpk) = 0.e0
290	DO jk = 1, jpkm1
291	zdn2dz (:,:,jk) = rn2(:,:,jk) - rn2(:,:,jk+1) ! Vertical profile of dN2/dz
292	!CDIR NOVERRCHK
293	zempba_3d_1(:,:,jk) = SQRT( MAX( 0.e0, rn2(:,:,jk) ) ) ! - - of N
294	zempba_3d_2(:,:,jk) = MAX( 0.e0, rn2(:,:,jk) ) ! - - of N^2
295	END DO
296	#endif
297	!
298	#if defined key_z_first
299	DO jj = 1, jpj
300	DO ji = 1, jpj
301	zsum1(ji,jj) = 0.e0
302	zsum2(ji,jj) = 0.e0
303	DO jk= 2, jpk
304	zsum1(ji,jj) = zsum1(ji,jj) + zempba_3d_1(ji,jj,jk) * fse3w(ji,jj,jk)
305	zsum2(ji,jj) = zsum2(ji,jj) + zempba_3d_2(ji,jj,jk) * fse3w(ji,jj,jk)
306	END DO
307	IF( zsum1(ji,jj) /= 0.e0 ) zsum1(ji,jj) = 1.e0 / zsum1(ji,jj)
308	IF( zsum2(ji,jj) /= 0.e0 ) zsum2(ji,jj) = 1.e0 / zsum2(ji,jj)
309	END DO
310	END DO
311	#else
312	zsum1(:,:) = 0.e0
313	zsum2(:,:) = 0.e0
314	DO jk= 2, jpk
315	zsum1(:,:) = zsum1(:,:) + zempba_3d_1(:,:,jk) * fse3w(:,:,jk)
316	zsum2(:,:) = zsum2(:,:) + zempba_3d_2(:,:,jk) * fse3w(:,:,jk)
317	END DO
318	DO jj = 1, jpj
319	DO ji = 1, jpi
320	IF( zsum1(ji,jj) /= 0.e0 ) zsum1(ji,jj) = 1.e0 / zsum1(ji,jj)
321	IF( zsum2(ji,jj) /= 0.e0 ) zsum2(ji,jj) = 1.e0 / zsum2(ji,jj)
322	END DO
323	END DO
324	#endif
325
326	zsum (:,:) = 0.e0
327
328	#if defined key_z_first
329	DO jj = 1, jpj
330	DO ji = 1, jpi
331	DO jk = 1, jpk
332	#else
333	DO jk = 1, jpk
334	DO jj = 1, jpj
335	DO ji = 1, jpi
336	#endif
337	zcoef = 0.5 - SIGN( 0.5, zdn2dz(ji,jj,jk) ) ! =0 if dN2/dz > 0, =1 otherwise
338	ztpc = zempba_3d_1(ji,jj,jk) * zsum1(ji,jj) * zcoef &
339	& + zempba_3d_2(ji,jj,jk) * zsum2(ji,jj) * ( 1. - zcoef )
340	!
341	zempba_3d(ji,jj,jk) = ztpc
342	zsum (ji,jj) = zsum(ji,jj) + ztpc * fse3w(ji,jj,jk)
343	END DO
344	#if !defined key_z_first
345	END DO
346	END DO
347	DO jj = 1, jpj
348	DO ji = 1, jpi
349	#endif
350	IF( zsum(ji,jj) > 0.e0 ) zsum(ji,jj) = 1.e0 / zsum(ji,jj)
351	END DO
352	END DO
353
354	! ! first estimation bounded by 10 cm2/s (with n2 bounded by rn_n2min)
355	zcoef = rn_tfe_itf / ( rn_tfe * rau0 )
356	#if defined key_z_first
357	DO jj = 1, jpj
358	DO ji = 1, jpi
359	DO jk = 1, jpk
360	zavt_itf(ji,jj,jk) = MIN( 10.e-4, zcoef * en_tmx(ji,jj) * zsum(ji,jj) * zempba_3d(ji,jj,jk) &
361	& / MAX( rn_n2min, rn2(ji,jj,jk) ) * tmask(ji,jj,jk) )
362	END DO
363	END DO
364	END DO
365	#else
366	DO jk = 1, jpk
367	zavt_itf(:,:,jk) = MIN( 10.e-4, zcoef * en_tmx(:,:) * zsum(:,:) * zempba_3d(:,:,jk) &
368	& / MAX( rn_n2min, rn2(:,:,jk) ) * tmask(:,:,jk) )
369	END DO
370	#endif
371
372	#if defined key_z_first
373	DO jj = 1, jpj
374	DO ji = 1, jpi
375	zkz(ji,jj) = 0.e0 ! Associated potential energy consummed over the whole water column
376	DO jk = 2, jpkm1
377	zkz(ji,jj) = zkz(ji,jj) + fse3w(ji,jj,jk) &
378	& * MAX( 0.e0, rn2(ji,jj,jk) ) * rau0 * zavt_itf(ji,jj,jk) * tmask(ji,jj,jk)
379	END DO
380	END DO
381	END DO
382	#else
383	zkz(:,:) = 0.e0 ! Associated potential energy consummed over the whole water column
384	DO jk = 2, jpkm1
385	zkz(:,:) = zkz(:,:) + fse3w(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zavt_itf(:,:,jk) * tmask(:,:,jk)
386	END DO
387	#endif
388
389	DO jj = 1, jpj ! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx
390	DO ji = 1, jpi
391	IF( zkz(ji,jj) /= 0.e0 ) zkz(ji,jj) = en_tmx(ji,jj) * rn_tfe_itf / rn_tfe / zkz(ji,jj)
392	END DO
393	END DO
394
395	#if defined key_z_first
396	DO jj = 1, jpj
397	DO ji = 1, jpi
398	zcoef = MIN( zkz(:,:), 120./10. ) ! kz max = 120 cm2/s
399	DO jk = 2, jpkm1 ! Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zavt_itf bound by 300 cm2/s
400	zavt_itf(ji,jj,jk) = zavt_itf(ji,jj,jk) * zcoef
401	END DO
402	END DO
403	END DO
404	#else
405	DO jk = 2, jpkm1 ! Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zavt_itf bound by 300 cm2/s
406	zavt_itf(:,:,jk) = zavt_itf(:,:,jk) * MIN( zkz(:,:), 120./10. ) ! kz max = 120 cm2/s
407	END DO
408	#endif
409
410	IF( kt == nit000 ) THEN ! diagnose the energy consumed by zavt_itf
411	ztpc = 0.e0
412	#if defined key_z_first
413	DO jj = 1, jpj
414	DO ji = 1, jpi
415	DO jk = 1, jpk
416	#else
417	DO jk = 1, jpk
418	DO jj = 1, jpj
419	DO ji = 1, jpi
420	#endif
421	ztpc = ztpc + e1t(ji,jj) * e2t(ji,jj) * fse3w(ji,jj,jk) * MAX( 0.e0, rn2(ji,jj,jk) ) &
422	& * zavt_itf(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj)
423	END DO
424	END DO
425	END DO
426	ztpc= rau0 * ztpc / ( rn_me * rn_tfe_itf )
427	IF(lwp) WRITE(numout,) ' N Total power consumption by zavt_itf: ztpc = ', ztpc 1.e-12 ,'TW'
428	ENDIF
429
430	! ! Update pav with the ITF mixing coefficient
431	#if defined key_z_first
432	DO jj = 1, jpj
433	DO ji = 1, jpi
434	DO jk = 2, jpkm1
435	pav(ji,jj,jk) = pav (ji,jj,jk) * ( 1.e0 - mask_itf(ji,jj) ) &
436	& + zavt_itf(ji,jj,jk) * mask_itf(ji,jj)
437	END DO
438	END DO
439	END DO
440	#else
441	DO jk = 2, jpkm1
442	pav(:,:,jk) = pav (:,:,jk) * ( 1.e0 - mask_itf(:,:) ) &
443	& + zavt_itf(:,:,jk) * mask_itf(:,:)
444	END DO
445	#endif
446	!
447	IF( wrk_not_released(2, 2,3,4,5) .OR. &
448	wrk_not_released(3, 1,2,3,4,5) )THEN
449	CALL ctl_stop('tmx_itf : failed to release workspace arrays.')
450	END IF
451	!
452	END SUBROUTINE tmx_itf
453
454	!! * Reset control of array index permutation
455	# include "oce_ftrans.h90"
456	# include "dom_oce_ftrans.h90"
457	# include "zdf_oce_ftrans.h90"
458	!FTRANS az_tmx :I :I :z
459
460	SUBROUTINE zdf_tmx_init
461	!!----------------------------------------------------------------------
462	!! * ROUTINE zdf_tmx_init *
463	!!
464	!! ** Purpose : Initialization of the vertical tidal mixing, Reading
465	!! of M2 and K1 tidal energy in nc files
466	!!
467	!! ** Method : - Read the namtmx namelist and check the parameters
468	!!
469	!! - Read the input data in NetCDF files :
470	!! M2 and K1 tidal energy. The total tidal energy, en_tmx,
471	!! is the sum of M2, K1 and S2 energy where S2 is assumed
472	!! to be: S2=(1/2)^2 * M2
473	!! mask_itf, a mask array that determine where substituing
474	!! the standard Simmons et al. (2005) formulation with the
475	!! one of Koch_Larrouy et al. (2007).
476	!!
477	!! - Compute az_tmx, a 3D coefficient that allows to compute
478	!! the standard tidal-induced vertical mixing as follows:
479	!! Kz_tides = az_tmx / max( rn_n2min, N^2 )
480	!! with az_tmx a bottom intensified coefficient is given by:
481	!! az_tmx(z) = en_tmx / ( rau0 * rn_htmx ) * EXP( -(H-z)/rn_htmx )
482	!! / ( 1. - EXP( - H /rn_htmx ) )
483	!! where rn_htmx the characteristic length scale of the bottom
484	!! intensification, en_tmx the tidal energy, and H the ocean depth
485	!!
486	!! ** input : - Namlist namtmx
487	!! - NetCDF file : M2_ORCA2.nc, K1_ORCA2.nc, and mask_itf.nc
488	!!
489	!! ** Action : - Increase by 1 the nstop flag is setting problem encounter
490	!! - defined az_tmx used to compute tidal-induced mixing
491	!!
492	!! References : Simmons et al. 2004, Ocean Modelling, 6, 3-4, 245-263.
493	!! Koch-Larrouy et al. 2007, GRL.
494	!!----------------------------------------------------------------------
495	USE oce , zav_tide => ua ! ua used as workspace
496	USE wrk_nemo, ONLY: zem2 => wrk_2d_1 ! read M2 and
497	USE wrk_nemo, ONLY: zek1 => wrk_2d_2 ! K1 tidal energy
498	USE wrk_nemo, ONLY: zkz => wrk_2d_3 ! total M2, K1 and S2 tidal energy
499	USE wrk_nemo, ONLY: zfact => wrk_2d_4 ! used for vertical structure function
500	USE wrk_nemo, ONLY: zhdep => wrk_2d_5 ! Ocean depth
501	USE wrk_nemo, ONLY: zpc => wrk_3d_1 ! power consumption
502
503	!! DCSE_NEMO: need additional directives for renamed module variables
504	!FTRANS zpc :I :I :z
505
506	!!
507	INTEGER :: ji, jj, jk ! dummy loop indices
508	INTEGER :: inum ! local integer
509	REAL(wp) :: ztpc, ze_z ! local scalars
510	#if defined key_z_first
511	REAL(wp) :: zcoef ! local scalar
512	#endif
513
514	!!
515	NAMELIST/namzdf_tmx/ rn_htmx, rn_n2min, rn_tfe, rn_me, ln_tmx_itf, rn_tfe_itf
516	!!----------------------------------------------------------------------
517
518	IF( wrk_in_use(2, 1,2,3,4,5) .OR. wrk_in_use(3, 1) ) THEN
519	CALL ctl_stop('zdf_tmx_init : requested workspace arrays unavailable.') ; RETURN
520	END IF
521
522	REWIND( numnam ) ! Read Namelist namtmx : Tidal Mixing
523	READ ( numnam, namzdf_tmx )
524
525	IF(lwp) THEN ! Control print
526	WRITE(numout,*)
527	WRITE(numout,*) 'zdf_tmx_init : tidal mixing'
528	WRITE(numout,*) '~~~~~~~~~~~~'
529	WRITE(numout,*) ' Namelist namzdf_tmx : set tidal mixing parameters'
530	WRITE(numout,*) ' Vertical decay scale for turbulence = ', rn_htmx
531	WRITE(numout,*) ' Brunt-Vaisala frequency threshold = ', rn_n2min
532	WRITE(numout,*) ' Tidal dissipation efficiency = ', rn_tfe
533	WRITE(numout,*) ' Mixing efficiency = ', rn_me
534	WRITE(numout,*) ' ITF specific parameterisation = ', ln_tmx_itf
535	WRITE(numout,*) ' ITF tidal dissipation efficiency = ', rn_tfe_itf
536	ENDIF
537
538	! ! allocate tmx arrays
539	IF( zdf_tmx_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_tmx_init : unable to allocate tmx arrays' )
540
541	IF( ln_tmx_itf ) THEN ! read the Indonesian Through Flow mask
542	CALL iom_open('mask_itf',inum)
543	CALL iom_get (inum, jpdom_data, 'tmaskitf',mask_itf,1) !
544	CALL iom_close(inum)
545	ENDIF
546
547	! read M2 tidal energy flux : W/m2 ( zem2 < 0 )
548	CALL iom_open('M2rowdrg',inum)
549	CALL iom_get (inum, jpdom_data, 'field',zem2,1) !
550	CALL iom_close(inum)
551
552	! read K1 tidal energy flux : W/m2 ( zek1 < 0 )
553	CALL iom_open('K1rowdrg',inum)
554	CALL iom_get (inum, jpdom_data, 'field',zek1,1) !
555	CALL iom_close(inum)
556
557	! Total tidal energy ( M2, S2 and K1 with S2=(1/2)^2 * M2 )
558	! only the energy available for mixing is taken into account,
559	! (mixing efficiency tidal dissipation efficiency)
560	en_tmx(:,:) = - rn_tfe * rn_me * ( zem2(:,:) * 1.25 + zek1(:,:) ) * tmask(:,:,1)
561
562	! Vertical structure (az_tmx)
563	DO jj = 1, jpj ! part independent of the level
564	DO ji = 1, jpi
565	zhdep(ji,jj) = fsdepw(ji,jj,mbkt(ji,jj)+1) ! depth of the ocean
566	zfact(ji,jj) = rau0 * rn_htmx * ( 1. - EXP( -zhdep(ji,jj) / rn_htmx ) )
567	IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = en_tmx(ji,jj) / zfact(ji,jj)
568	END DO
569	END DO
570	#if defined key_z_first
571	DO jj = 1, jpj
572	DO ji = 1, jpi
573	DO jk= 1, jpk ! complete with the level-dependent part
574	#else
575	DO jk= 1, jpk ! complete with the level-dependent part
576	DO jj = 1, jpj
577	DO ji = 1, jpi
578	#endif
579	az_tmx(ji,jj,jk) = zfact(ji,jj) * EXP( -( zhdep(ji,jj)-fsdepw(ji,jj,jk) ) / rn_htmx ) * tmask(ji,jj,jk)
580	END DO
581	END DO
582	END DO
583
584	IF( nprint == 1 .AND. lwp ) THEN
585	! Control print
586	! Total power consumption due to vertical mixing
587	! zpc = rau0 * 1/rn_me * rn2 * zav_tide
588	#if defined key_z_first
589	DO jj = 1, jpj
590	DO ji = 1, jpi
591	zav_tide(ji,jj,1) = 0.e0
592	DO jk = 2, jpkm1
593	zav_tide(:,:,jk) = az_tmx(:,:,jk) / MAX( rn_n2min, rn2(:,:,jk) )
594	END DO
595	zav_tide(ji,jj,jpk) = 0.e0
596	END DO
597	END DO
598	#else
599	zav_tide(:,:,:) = 0.e0
600	DO jk = 2, jpkm1
601	zav_tide(:,:,jk) = az_tmx(:,:,jk) / MAX( rn_n2min, rn2(:,:,jk) )
602	END DO
603	#endif
604
605	ztpc = 0.e0
606	zpc(:,:,:) = MAX(rn_n2min,rn2(:,:,:)) * zav_tide(:,:,:)
607	#if defined key_z_first
608	DO jj = 1, jpj
609	DO ji = 1, jpi
610	DO jk= 2, jpkm1
611	#else
612	DO jk= 2, jpkm1
613	DO jj = 1, jpj
614	DO ji = 1, jpi
615	#endif
616	ztpc = ztpc + fse3w(ji,jj,jk) * e1t(ji,jj) * e2t(ji,jj) * zpc(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj)
617	END DO
618	END DO
619	END DO
620	ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc
621
622	WRITE(numout,*)
623	WRITE(numout,) ' Total power consumption of the tidally driven part of Kz : ztpc = ', ztpc 1.e-12 ,'TW'
624
625
626	! control print 2
627	zav_tide(:,:,:) = MIN( zav_tide(:,:,:), 60.e-4 )
628	#if defined key_z_first
629	DO jj = 1, jpj
630	DO ji = 1, jpi
631	zkz(ji,jj) = 0.e0
632	DO jk = 2, jpkm1
633	#else
634	zkz(:,:) = 0.e0
635	DO jk = 2, jpkm1
636	DO jj = 1, jpj
637	DO ji = 1, jpi
638	#endif
639	zkz(ji,jj) = zkz(ji,jj) + fse3w(ji,jj,jk) &
640	& * MAX( 0.e0, rn2(ji,jj,jk) ) * rau0 * zav_tide(ji,jj,jk)* tmask(ji,jj,jk)
641	END DO
642	END DO
643	END DO
644	! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz
645	DO jj = 1, jpj
646	DO ji = 1, jpi
647	IF( zkz(ji,jj) /= 0.e0 ) THEN
648	zkz(ji,jj) = en_tmx(ji,jj) / zkz(ji,jj)
649	ENDIF
650	END DO
651	END DO
652	ztpc = 1.e50
653	DO jj = 1, jpj
654	DO ji = 1, jpi
655	IF( zkz(ji,jj) /= 0.e0 ) THEN
656	ztpc = MIN( zkz(ji,jj), ztpc)
657	ENDIF
658	END DO
659	END DO
660	WRITE(numout,*) ' Min de zkz ', ztpc, ' Max = ', MAXVAL(zkz(:,:) )
661
662	#if defined key_z_first
663	DO jj = 1, jpj
664	DO ji = 1, jpi
665	zcoef = MIN( zkz(ji,jj), 30./6. ) !kz max = 300 cm2/s
666	DO jk = 2, jpkm1
667	zav_tide(ji,jj,jk) = zav_tide(ji,jj,jk) * zcoef
668	END DO
669	END DO
670	END DO
671	#else
672	DO jk = 2, jpkm1
673	zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) !kz max = 300 cm2/s
674	END DO
675	#endif
676	ztpc = 0.e0
677	zpc(:,:,:) = MAX(0.e0,rn2(:,:,:)) * zav_tide(:,:,:)
678	#if defined key_z_first
679	DO jj = 1, jpj
680	DO ji = 1, jpi
681	DO jk= 1, jpk
682	#else
683	DO jk= 1, jpk
684	DO jj = 1, jpj
685	DO ji = 1, jpi
686	#endif
687	ztpc = ztpc + fse3w(ji,jj,jk) * e1t(ji,jj) * e2t(ji,jj) * zpc(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj)
688	END DO
689	END DO
690	END DO
691	ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc
692	WRITE(numout,) ' 2 Total power consumption of the tidally driven part of Kz : ztpc = ', ztpc 1.e-12 ,'TW'
693
694	DO jk = 1, jpk
695	ze_z = SUM( e1t(:,:) * e2t(:,:) * zav_tide(:,:,jk) * tmask_i(:,:) ) &
696	& / MAX( 1.e-20, SUM( e1t(:,:) * e2t(:,:) * tmask (:,:,jk) * tmask_i(:,:) ) )
697	ztpc = 1.E50
698	DO jj = 1, jpj
699	DO ji = 1, jpi
700	IF( zav_tide(ji,jj,jk) /= 0.e0 ) ztpc =Min( ztpc, zav_tide(ji,jj,jk) )
701	END DO
702	END DO
703	WRITE(numout,) ' N2 min - jk= ', jk,' ', ze_z 1.e4,' cm2/s min= ',ztpc*1.e4, &
704	& 'max= ', MAXVAL(zav_tide(:,:,jk) )*1.e4, ' cm2/s'
705	END DO
706
707	WRITE(numout,) ' e_tide : ', SUM( e1te2ten_tmx ) / ( rn_tfe rn_me ) * 1.e-12, 'TW'
708	WRITE(numout,*)
709	WRITE(numout,*) ' Initial profile of tidal vertical mixing'
710	DO jk = 1, jpk
711	DO jj = 1,jpj
712	DO ji = 1,jpi
713	zkz(ji,jj) = az_tmx(ji,jj,jk) /MAX( rn_n2min, rn2(ji,jj,jk) )
714	END DO
715	END DO
716	ze_z = SUM( e1t(:,:) * e2t(:,:) * zkz(:,:) * tmask_i(:,:) ) &
717	& / MAX( 1.e-20, SUM( e1t(:,:) * e2t(:,:) * tmask (:,:,jk) * tmask_i(:,:) ) )
718	WRITE(numout,) ' jk= ', jk,' ', ze_z 1.e4,' cm2/s'
719	END DO
720	DO jk = 1, jpk
721	zkz(:,:) = az_tmx(:,:,jk) /rn_n2min
722	ze_z = SUM( e1t(:,:) * e2t(:,:) * zkz(:,:) * tmask_i(:,:) ) &
723	& / MAX( 1.e-20, SUM( e1t(:,:) * e2t(:,:) * tmask (:,:,jk) * tmask_i(:,:) ) )
724	WRITE(numout,*)
725	WRITE(numout,) ' N2 min - jk= ', jk,' ', ze_z 1.e4,' cm2/s min= ',MINVAL(zkz)*1.e4, &
726	& 'max= ', MAXVAL(zkz)*1.e4, ' cm2/s'
727	END DO
728	!
729	ENDIF
730	!
731	IF(wrk_not_released(2, 1,2,3,4,5) .OR. &
732	wrk_not_released(3, 1) ) CALL ctl_stop( 'zdf_tmx_init : failed to release workspace arrays' )
733	!
734	END SUBROUTINE zdf_tmx_init
735
736	#else
737	!!----------------------------------------------------------------------
738	!! Default option Dummy module NO Tidal MiXing
739	!!----------------------------------------------------------------------
740	LOGICAL, PUBLIC, PARAMETER :: lk_zdftmx = .FALSE. !: tidal mixing flag
741	CONTAINS
742	SUBROUTINE zdf_tmx_init ! Dummy routine
743	WRITE(,) 'zdf_tmx: You should not have seen this print! error?'
744	END SUBROUTINE zdf_tmx_init
745	SUBROUTINE zdf_tmx( kt ) ! Dummy routine
746	WRITE(,) 'zdf_tmx: You should not have seen this print! error?', kt
747	END SUBROUTINE zdf_tmx
748	#endif
749
750	!!======================================================================
751	END MODULE zdftmx

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

Download in other formats: