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.

dtadyn.F90 in branches/DEV_r2006_merge_TRA_TRC/NEMO/OFF_SRC – NEMO

Context Navigation

source: branches/DEV_r2006_merge_TRA_TRC/NEMO/OFF_SRC/dtadyn.F90 @ 2053

Last change on this file since 2053 was 2053, checked in by cetlod, 14 years ago
improve the offline part to take into account the merge of TRA-TRC, see ticket:702
Property svn:eol-style set to `native` Property svn:keywords set to `Id`
File size: 31.5 KB

Line
1	MODULE dtadyn
2	!!======================================================================
3	!! * MODULE dtadyn *
4	!! OFFLINE : interpolation of the physical fields
5	!!=====================================================================
6
7	!!----------------------------------------------------------------------
8	!! dta_dyn_init : initialization, namelist read, and parameters control
9	!! dta_dyn : Interpolation of the fields
10	!!----------------------------------------------------------------------
11	!! * Modules used
12	USE oce ! ocean dynamics and tracers variables
13	USE dom_oce ! ocean space and time domain variables
14	USE zdf_oce ! ocean vertical physics
15	USE in_out_manager ! I/O manager
16	USE phycst ! physical constants
17	USE sbc_oce
18	USE trabbl
19	USE ldfslp
20	USE ldfeiv ! eddy induced velocity coef.
21	USE ldftra_oce ! ocean tracer lateral physics
22	USE zdfmxl
23	USE eosbn2
24	USE zdfddm ! vertical physics: double diffusion
25	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
26	USE zpshde
27	USE lib_mpp ! distributed memory computing library
28
29	IMPLICIT NONE
30	PRIVATE
31
32	!! * Routine accessibility
33	PUBLIC dta_dyn_init ! called by opa.F90
34	PUBLIC dta_dyn ! called by step.F90
35
36	LOGICAL , PUBLIC :: &
37	lperdyn = .TRUE. , & ! boolean for periodic fields or not
38	lfirdyn = .TRUE. ! boolean for the first call or not
39
40	INTEGER , PUBLIC :: &
41	ndtadyn = 73 , & ! Number of dat in one year
42	ndtatot = 73 , & ! Number of data in the input field
43	nsptint = 1 ! type of spatial interpolation
44
45	CHARACTER(len=45) :: &
46	cfile_grid_T = 'dyna_grid_T.nc', & !: name of the grid_T file
47	cfile_grid_U = 'dyna_grid_U.nc', & !: name of the grid_U file
48	cfile_grid_V = 'dyna_grid_V.nc', & !: name of the grid_V file
49	cfile_grid_W = 'dyna_grid_W.nc' !: name of the grid_W file
50
51	REAL(wp) :: &
52	rnspdta , & !: number of time step per 2 consecutives data
53	rnspdta2 !: rnspdta * 0.5
54
55	INTEGER :: &
56	ndyn1, ndyn2 , &
57	nlecoff = 0 , & ! switch for the first read
58	numfl_t, numfl_u, &
59	numfl_v, numfl_w
60
61	REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: &
62	tdta , & ! temperature at two consecutive times
63	sdta , & ! salinity at two consecutive times
64	udta , & ! zonal velocity at two consecutive times
65	vdta , & ! meridional velocity at two consecutive times
66	wdta , & ! vertical velocity at two consecutive times
67	avtdta ! vertical diffusivity coefficient
68
69	REAL(wp), DIMENSION(jpi,jpj,2) :: &
70	hmlddta, & ! mixed layer depth at two consecutive times
71	wspddta, & ! wind speed at two consecutive times
72	frlddta, & ! sea-ice fraction at two consecutive times
73	empdta , & ! E-P at two consecutive times
74	qsrdta ! short wave heat flux at two consecutive times
75
76	#if defined key_ldfslp
77	REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: &
78	uslpdta , & ! zonal isopycnal slopes
79	vslpdta , & ! meridional isopycnal slopes
80	wslpidta , & ! zonal diapycnal slopes
81	wslpjdta ! meridional diapycnal slopes
82	#endif
83
84	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
85	REAL(wp), DIMENSION(jpi,jpj,2) :: &
86	aeiwdta ! G&M coefficient
87	#endif
88
89	#if defined key_degrad
90	REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: &
91	ahtudta, ahtvdta, ahtwdta ! Lateral diffusivity
92	# if defined key_traldf_eiv
93	REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: &
94	aeiudta, aeivdta, aeiwdta ! G&M coefficient
95	# endif
96
97	#endif
98
99	!! * Substitutions
100	# include "domzgr_substitute.h90"
101	# include "vectopt_loop_substitute.h90"
102	!!----------------------------------------------------------------------
103	!! OPA 9.0 , LOCEAN-IPSL (2005)
104	!! $Id$
105	!! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
106	!!----------------------------------------------------------------------
107
108	CONTAINS
109
110	SUBROUTINE dta_dyn( kt )
111	!!----------------------------------------------------------------------
112	!! * ROUTINE dta_dyn *
113	!!
114	!! ** Purpose : Prepares dynamics and physics fields from an
115	!! OPA9 simulation for an off-line simulation
116	!! for passive tracer
117	!!
118	!! ** Method : calculates the position of DATA to read READ DATA
119	!! (example month changement) computes slopes IF needed
120	!! interpolates DATA IF needed
121	!!
122	!! ** History :
123	!! ! original : 92-01 (M. Imbard: sub domain)
124	!! ! addition : 98-04 (L.Bopp MA Foujols: slopes for isopyc.)
125	!! ! addition : 98-05 (L. Bopp read output of coupled run)
126	!! ! addition : 05-03 (O. Aumont and A. El Moussaoui) F90
127	!! ! addition : 05-12 (C. Ethe) Adapted for DEGINT
128	!!----------------------------------------------------------------------
129	!! * Arguments
130	INTEGER, INTENT( in ) :: kt ! ocean time-step index
131
132	!! * Local declarations
133	INTEGER :: iper, iperm1, iswap, izt
134
135	REAL(wp) :: zt
136	REAL(wp) :: zweigh
137	!!----------------------------------------------------------------------
138
139	zt = ( FLOAT (kt) + rnspdta2 ) / rnspdta
140	izt = INT( zt )
141	zweigh = zt - FLOAT( INT(zt) )
142
143	IF( lperdyn ) THEN
144	iperm1 = MOD( izt, ndtadyn )
145	ELSE
146	iperm1 = MOD( izt, ndtatot - 1 ) + 1
147	ENDIF
148
149	iper = iperm1 + 1
150	IF( iperm1 == 0 ) THEN
151	IF( lperdyn ) THEN
152	iperm1 = ndtadyn
153	ELSE
154	IF( lfirdyn ) THEN
155	IF (lwp) WRITE (numout,*) &
156	& ' dynamic file is not periodic with or without interpolation &
157	& we take the first value for the previous period iperm1 = 0 '
158	END IF
159	END IF
160	END IF
161
162	iswap = 0
163
164	! 1. First call lfirdyn = true
165	! ----------------------------
166
167	IF( lfirdyn ) THEN
168	! store the information of the period read
169	ndyn1 = iperm1
170	ndyn2 = iper
171
172	IF (lwp) THEN
173	WRITE (numout,*) ' dynamics data read for the period ndyn1 =',ndyn1, &
174	& ' and for the period ndyn2 = ',ndyn2
175	WRITE (numout,*) ' time step is : ', kt
176	WRITE (numout,*) ' we have ndtadyn = ',ndtadyn,' records in the dynamic file for one year'
177	END IF
178	!
179	IF( iperm1 /= 0 ) THEN ! data read for the iperm1 period
180	CALL dynrea( kt, iperm1 )
181	ELSE
182	CALL dynrea( kt, 1 )
183	ENDIF
184
185	IF( lk_ldfslp ) THEN
186	! Computes slopes
187	! Caution : here tn, sn and avt are used as workspace
188	tn (:,:,:) = tdta (:,:,:,2)
189	sn (:,:,:) = sdta (:,:,:,2)
190	avt(:,:,:) = avtdta(:,:,:,2)
191
192	CALL eos( tn, sn, rhd, rhop ) ! Time-filtered in situ density
193	CALL bn2( tn, sn, rn2 ) ! before Brunt-Vaisala frequency
194	IF( ln_zps ) &
195	& CALL zps_hde( kt, tn , sn , rhd, & ! Partial steps: before Horizontal DErivative
196	& gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level
197	& gtv, gsv, grv )
198	CALL zdf_mxl( kt ) ! mixed layer depth
199	CALL ldf_slp( kt, rhd, rn2 )
200
201	uslpdta (:,:,:,2) = uslp (:,:,:)
202	vslpdta (:,:,:,2) = vslp (:,:,:)
203	wslpidta(:,:,:,2) = wslpi(:,:,:)
204	wslpjdta(:,:,:,2) = wslpj(:,:,:)
205	END IF
206
207	! swap from record 2 to 1
208	CALL swap_dyn_data
209
210	iswap = 1 ! indicates swap
211
212	CALL dynrea( kt, iper ) ! data read for the iper period
213
214	IF( lk_ldfslp ) THEN
215	! Computes slopes
216	! Caution : here tn, sn and avt are used as workspace
217	tn (:,:,:) = tdta (:,:,:,2)
218	sn (:,:,:) = sdta (:,:,:,2)
219	avt(:,:,:) = avtdta(:,:,:,2)
220
221	CALL eos( tn, sn, rhd, rhop ) ! Time-filtered in situ density
222	CALL bn2( tn, sn, rn2 ) ! before Brunt-Vaisala frequency
223	IF( ln_zps ) &
224	& CALL zps_hde( kt, tn , sn , rhd, & ! Partial steps: before Horizontal DErivative
225	& gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level
226	& gtv, gsv, grv )
227	CALL zdf_mxl( kt ) ! mixed layer depth
228	CALL ldf_slp( kt, rhd, rn2 )
229
230	uslpdta (:,:,:,2) = uslp (:,:,:)
231	vslpdta (:,:,:,2) = vslp (:,:,:)
232	wslpidta(:,:,:,2) = wslpi(:,:,:)
233	wslpjdta(:,:,:,2) = wslpj(:,:,:)
234	END IF
235	!
236	lfirdyn=.FALSE. ! trace the first call
237	ENDIF
238	!
239	! And now what we have to do at every time step
240	! check the validity of the period in memory
241	!
242	IF( iperm1 /= ndyn1 ) THEN
243
244	IF( iperm1 == 0. ) THEN
245	IF (lwp) THEN
246	WRITE (numout,*) ' dynamic file is not periodic with periodic interpolation'
247	WRITE (numout,*) ' we take the last value for the last period '
248	WRITE (numout,*) ' iperm1 = 12, iper = 13 '
249	ENDIF
250	iperm1 = 12
251	iper = 13
252	ENDIF
253	!
254	! We have to prepare a new read of data : swap from record 2 to 1
255	!
256	CALL swap_dyn_data
257
258	iswap = 1 ! indicates swap
259
260	CALL dynrea( kt, iper ) ! data read for the iper period
261
262	IF( lk_ldfslp ) THEN
263	! Computes slopes
264	! Caution : here tn, sn and avt are used as workspace
265	tn (:,:,:) = tdta (:,:,:,2)
266	sn (:,:,:) = sdta (:,:,:,2)
267	avt(:,:,:) = avtdta(:,:,:,2)
268
269	CALL eos( tn, sn, rhd, rhop ) ! Time-filtered in situ density
270	CALL bn2( tn, sn, rn2 ) ! before Brunt-Vaisala frequency
271	IF( ln_zps ) &
272	& CALL zps_hde( kt, tn , sn , rhd, & ! Partial steps: before Horizontal DErivative
273	& gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level
274	& gtv, gsv, grv )
275	CALL zdf_mxl( kt ) ! mixed layer depth
276	CALL ldf_slp( kt, rhd, rn2 )
277
278	uslpdta (:,:,:,2) = uslp (:,:,:)
279	vslpdta (:,:,:,2) = vslp (:,:,:)
280	wslpidta(:,:,:,2) = wslpi(:,:,:)
281	wslpjdta(:,:,:,2) = wslpj(:,:,:)
282	END IF
283
284	! store the information of the period read
285	ndyn1 = ndyn2
286	ndyn2 = iper
287
288	IF (lwp) THEN
289	WRITE (numout,*) ' dynamics data read for the period ndyn1 =',ndyn1, &
290	& ' and for the period ndyn2 = ',ndyn2
291	WRITE (numout,*) ' time step is : ', kt
292	END IF
293	!
294	END IF
295	!
296	! Compute the data at the given time step
297	!----------------------------------------
298
299	IF( nsptint == 0 ) THEN
300	! No spatial interpolation, data are probably correct
301	! We have to initialize data if we have changed the period
302	CALL assign_dyn_data
303	ELSE IF( nsptint == 1 ) THEN
304	! linear interpolation
305	CALL linear_interp_dyn_data( zweigh )
306	ELSE
307	! other interpolation
308	WRITE (numout,*) ' this kind of interpolation do not exist at the moment : we stop'
309	STOP 'dtadyn'
310	END IF
311
312	! In any case, we need rhop
313	CALL eos( tn, sn, rhd, rhop )
314
315	#if ! defined key_degrad && defined key_traldf_c2d
316	! In case of 2D varying coefficients, we need aeiv and aeiu
317	IF( lk_traldf_eiv ) CALL dta_eiv( kt ) ! eddy induced velocity coefficient
318	#endif
319
320	! Compute bbl coefficients if needed
321	IF( lk_trabbl ) THEN
322	tb(:,:,:) = tn(:,:,:)
323	sb(:,:,:) = sn(:,:,:)
324	CALL bbl( kt, 'TRC')
325	END IF
326
327	END SUBROUTINE dta_dyn
328
329	SUBROUTINE dynrea( kt, kenr )
330	!!----------------------------------------------------------------------
331	!! * ROUTINE dynrea *
332	!!
333	!! ** Purpose : READ dynamics fiels from OPA9 netcdf output
334	!!
335	!! ** Method : READ the kenr records of DATA and store in
336	!! in udta(...,2), ....
337	!!
338	!! ** History : additions : M. Levy et M. Benjelloul jan 2001
339	!! (netcdf FORMAT)
340	!! 05-03 (O. Aumont and A. El Moussaoui) F90
341	!! 06-07 : (C. Ethe) use of iom module
342	!!----------------------------------------------------------------------
343	!! * Modules used
344	USE iom
345
346	!! * Arguments
347	INTEGER, INTENT( in ) :: kt, kenr ! time index
348	!! * Local declarations
349	INTEGER :: jkenr
350
351	REAL(wp), DIMENSION(jpi,jpj,jpk) :: &
352	zu, zv, zw, zt, zs, zavt , & ! 3-D dynamical fields
353	zhdiv ! horizontal divergence
354
355	REAL(wp), DIMENSION(jpi,jpj) :: &
356	zemp, zqsr, zmld, zice, zwspd, &
357	ztaux, ztauy
358
359	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
360	REAL(wp), DIMENSION(jpi,jpj) :: zaeiw
361	#endif
362
363	#if defined key_degrad
364	REAL(wp), DIMENSION(jpi,jpj,jpk) :: &
365	zahtu, zahtv, zahtw ! Lateral diffusivity
366	# if defined key_traldf_eiv
367	REAL(wp), DIMENSION(jpi,jpj,jpk) :: &
368	zaeiu, zaeiv, zaeiw ! G&M coefficient
369	# endif
370	#endif
371
372	!---------------------------------------------------------------
373	! 0. Initialization
374
375	! cas d'un fichier non periodique : on utilise deux fois le premier et
376	! le dernier champ temporel
377
378	jkenr = kenr
379
380	IF(lwp) THEN
381	WRITE(numout,*)
382	WRITE(numout,*) 'Dynrea : reading dynamical fields, kenr = ', jkenr
383	WRITE(numout,*) ' ~~~~~~~'
384	#if defined key_degrad
385	WRITE(numout,*) ' Degraded fields'
386	#endif
387	WRITE(numout,*)
388	ENDIF
389
390
391	IF( kt == nit000 .AND. nlecoff == 0 ) THEN
392	nlecoff = 1
393	CALL iom_open ( cfile_grid_T, numfl_t )
394	CALL iom_open ( cfile_grid_U, numfl_u )
395	CALL iom_open ( cfile_grid_V, numfl_v )
396	CALL iom_open ( cfile_grid_W, numfl_w )
397	ENDIF
398
399	! file grid-T
400	!---------------
401	CALL iom_get( numfl_t, jpdom_data, 'votemper', zt (:,:,:), jkenr )
402	CALL iom_get( numfl_t, jpdom_data, 'vosaline', zs (:,:,:), jkenr )
403	CALL iom_get( numfl_t, jpdom_data, 'somixhgt', zmld (:,: ), jkenr )
404	CALL iom_get( numfl_t, jpdom_data, 'sowaflcd', zemp (:,: ), jkenr )
405	CALL iom_get( numfl_t, jpdom_data, 'soshfldo', zqsr (:,: ), jkenr )
406	CALL iom_get( numfl_t, jpdom_data, 'soicecov', zice (:,: ), jkenr )
407	IF( iom_varid( numfl_t, 'sowindsp', ldstop = .FALSE. ) > 0 ) THEN
408	CALL iom_get( numfl_t, jpdom_data, 'sowindsp', zwspd(:,:), jkenr )
409	ELSE
410	CALL iom_get( numfl_u, jpdom_data, 'sozotaux', ztaux(:,:), jkenr )
411	CALL iom_get( numfl_v, jpdom_data, 'sometauy', ztauy(:,:), jkenr )
412	CALL tau2wnd( ztaux, ztauy, zwspd )
413	ENDIF
414	! files grid-U / grid_V
415	CALL iom_get( numfl_u, jpdom_data, 'vozocrtx', zu (:,:,:), jkenr )
416	CALL iom_get( numfl_v, jpdom_data, 'vomecrty', zv (:,:,:), jkenr )
417
418	! file grid-W
419	!! CALL iom_get ( numfl_w, jpdom_data, 'vovecrtz', zw (:,:,:), jkenr )
420	! Computation of vertical velocity using horizontal divergence
421	CALL wzv( zu, zv, zw, zhdiv )
422
423	# if defined key_zdfddm
424	CALL iom_get( numfl_w, jpdom_data, 'voddmavs', zavt (:,:,:), jkenr )
425	#else
426	CALL iom_get( numfl_w, jpdom_data, 'votkeavt', zavt (:,:,:), jkenr )
427	#endif
428
429	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
430	CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw (:,: ), jkenr )
431	#endif
432
433	#if defined key_degrad
434	CALL iom_get( numfl_u, jpdom_data, 'vozoahtu', zahtu(:,:,:), jkenr )
435	CALL iom_get( numfl_v, jpdom_data, 'vomeahtv', zahtv(:,:,:), jkenr )
436	CALL iom_get( numfl_w, jpdom_data, 'voveahtw', zahtw(:,:,:), jkenr )
437	# if defined key_traldf_eiv
438	CALL iom_get( numfl_u, jpdom_data, 'vozoaeiu', zaeiu(:,:,:), jkenr )
439	CALL iom_get( numfl_v, jpdom_data, 'vomeaeiv', zaeiv(:,:,:), jkenr )
440	CALL iom_get( numfl_w, jpdom_data, 'voveaeiw', zaeiw(:,:,:), jkenr )
441	# endif
442	#endif
443
444	udta(:,:,:,2) = zu(:,:,:) * umask(:,:,:)
445	vdta(:,:,:,2) = zv(:,:,:) * vmask(:,:,:)
446	wdta(:,:,:,2) = zw(:,:,:) * tmask(:,:,:)
447
448	tdta(:,:,:,2) = zt (:,:,:) * tmask(:,:,:)
449	sdta(:,:,:,2) = zs (:,:,:) * tmask(:,:,:)
450	avtdta(:,:,:,2) = zavt(:,:,:) * tmask(:,:,:)
451
452	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
453	aeiwdta(:,:,2) = zaeiw(:,:) * tmask(:,:,1)
454	#endif
455
456	#if defined key_degrad
457	ahtudta(:,:,:,2) = zahtu(:,:,:) * umask(:,:,:)
458	ahtvdta(:,:,:,2) = zahtv(:,:,:) * vmask(:,:,:)
459	ahtwdta(:,:,:,2) = zahtw(:,:,:) * tmask(:,:,:)
460	# if defined key_traldf_eiv
461	aeiudta(:,:,:,2) = zaeiu(:,:,:) * umask(:,:,:)
462	aeivdta(:,:,:,2) = zaeiv(:,:,:) * vmask(:,:,:)
463	aeiwdta(:,:,:,2) = zaeiw(:,:,:) * tmask(:,:,:)
464	# endif
465	#endif
466
467	! fluxes
468	!
469	wspddta(:,:,2) = zwspd(:,:) * tmask(:,:,1)
470	frlddta(:,:,2) = MIN( 1., zice(:,:) ) * tmask(:,:,1)
471	empdta (:,:,2) = zemp(:,:) * tmask(:,:,1)
472	qsrdta (:,:,2) = zqsr(:,:) * tmask(:,:,1)
473	hmlddta(:,:,2) = zmld(:,:) * tmask(:,:,1)
474
475	IF( kt == nitend ) THEN
476	CALL iom_close ( numfl_t )
477	CALL iom_close ( numfl_u )
478	CALL iom_close ( numfl_v )
479	CALL iom_close ( numfl_w )
480	ENDIF
481
482	END SUBROUTINE dynrea
483
484	SUBROUTINE dta_dyn_init
485	!!----------------------------------------------------------------------
486	!! * ROUTINE dta_dyn_init *
487	!!
488	!! ** Purpose : initializations of parameters for the interpolation
489	!!
490	!! ** Method :
491	!!
492	!! History :
493	!! ! original : 92-01 (M. Imbard: sub domain)
494	!! ! 98-04 (L.Bopp MA Foujols: slopes for isopyc.)
495	!! ! 98-05 (L. Bopp read output of coupled run)
496	!! ! 05-03 (O. Aumont and A. El Moussaoui) F90
497	!!----------------------------------------------------------------------
498	!! * Modules used
499
500	!! * Local declarations
501
502	REAL(wp) :: znspyr !: number of time step per year
503
504	NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, lperdyn, &
505	& cfile_grid_T, cfile_grid_U, cfile_grid_V, cfile_grid_W
506	!!----------------------------------------------------------------------
507
508	! Define the dynamical input parameters
509	! ======================================
510
511	! Read Namelist namdyn : Lateral physics on tracers
512	REWIND( numnam )
513	READ ( numnam, namdyn )
514
515	IF(lwp) THEN
516	WRITE(numout,*)
517	WRITE(numout,*) 'namdyn : offline dynamical selection'
518	WRITE(numout,*) '~~~~~~~'
519	WRITE(numout,*) ' Namelist namdyn : set parameters for the lecture of the dynamical fields'
520	WRITE(numout,*)
521	WRITE(numout,*) ' number of elements in the FILE for a year ndtadyn = ' , ndtadyn
522	WRITE(numout,*) ' total number of elements in the FILE ndtatot = ' , ndtatot
523	WRITE(numout,*) ' type of interpolation nsptint = ' , nsptint
524	WRITE(numout,*) ' loop on the same FILE lperdyn = ' , lperdyn
525	WRITE(numout,*) ' '
526	WRITE(numout,*) ' name of grid_T file cfile_grid_T = ', TRIM(cfile_grid_T)
527	WRITE(numout,*) ' name of grid_U file cfile_grid_U = ', TRIM(cfile_grid_U)
528	WRITE(numout,*) ' name of grid_V file cfile_grid_V = ', TRIM(cfile_grid_V)
529	WRITE(numout,*) ' name of grid_W file cfile_grid_W = ', TRIM(cfile_grid_W)
530	WRITE(numout,*) ' '
531	ENDIF
532
533	znspyr = nyear_len(1) * rday / rdt
534	rnspdta = znspyr / FLOAT( ndtadyn )
535	rnspdta2 = rnspdta * 0.5
536
537	CALL dta_dyn( nit000 )
538
539	END SUBROUTINE dta_dyn_init
540
541	SUBROUTINE wzv( pu, pv, pw, phdiv )
542	!!----------------------------------------------------------------------
543	!! * ROUTINE wzv *
544	!!
545	!! ** Purpose : Compute the now vertical velocity after the array swap
546	!!
547	!! ** Method :
548	!! ** Method : - Divergence:
549	!! - compute the now divergence given by :
550	!! * z-coordinate
551	!! hdiv = 1/(e1t*e2t) [ di(e2u u) + dj(e1v v) ]
552	!! - Using the incompressibility hypothesis, the vertical
553	!! velocity is computed by integrating the horizontal divergence
554	!! from the bottom to the surface.
555	!! The boundary conditions are w=0 at the bottom (no flux) and,
556	!! in regid-lid case, w=0 at the sea surface.
557	!!
558	!!
559	!! History :
560	!! 9.0 ! 02-07 (G. Madec) Vector optimization
561	!!----------------------------------------------------------------------
562	!! * Arguments
563	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pu, pv !: horizontal velocities
564	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out ) :: pw !: verticla velocity
565	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: phdiv !: horizontal divergence
566
567	!! * Local declarations
568	INTEGER :: ji, jj, jk
569	REAL(wp) :: zu, zu1, zv, zv1, zet
570
571
572	! Computation of vertical velocity using horizontal divergence
573	phdiv(:,:,:) = 0.
574	DO jk = 1, jpkm1
575	DO jj = 2, jpjm1
576	DO ji = fs_2, fs_jpim1 ! vector opt.
577	#if defined key_zco
578	zu = pu(ji ,jj ,jk) * umask(ji ,jj ,jk) * e2u(ji ,jj )
579	zu1 = pu(ji-1,jj ,jk) * umask(ji-1,jj ,jk) * e2u(ji-1,jj )
580	zv = pv(ji ,jj ,jk) * vmask(ji ,jj ,jk) * e1v(ji ,jj )
581	zv1 = pv(ji ,jj-1,jk) * vmask(ji ,jj-1,jk) * e1v(ji ,jj-1)
582	zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) )
583	#else
584	zu = pu(ji ,jj ,jk) * umask(ji ,jj ,jk) * e2u(ji ,jj ) * fse3u(ji ,jj ,jk)
585	zu1 = pu(ji-1,jj ,jk) * umask(ji-1,jj ,jk) * e2u(ji-1,jj ) * fse3u(ji-1,jj ,jk)
586	zv = pv(ji ,jj ,jk) * vmask(ji ,jj ,jk) * e1v(ji ,jj ) * fse3v(ji ,jj ,jk)
587	zv1 = pv(ji ,jj-1,jk) * vmask(ji ,jj-1,jk) * e1v(ji ,jj-1) * fse3v(ji ,jj-1,jk)
588	zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
589	#endif
590	phdiv(ji,jj,jk) = ( zu - zu1 + zv - zv1 ) * zet
591	END DO
592	END DO
593	ENDDO
594
595	! Lateral boundary conditions on phdiv
596	CALL lbc_lnk( phdiv, 'T', 1. )
597
598
599	! computation of vertical velocity from the bottom
600	pw(:,:,jpk) = 0.
601	DO jk = jpkm1, 1, -1
602	pw(:,:,jk) = pw(:,:,jk+1) - fse3t(:,:,jk) * phdiv(:,:,jk)
603	END DO
604
605	END SUBROUTINE wzv
606
607	SUBROUTINE dta_eiv( kt )
608	!!----------------------------------------------------------------------
609	!! * ROUTINE dta_eiv *
610	!!
611	!! ** Purpose : Compute the eddy induced velocity coefficient from the
612	!! growth rate of baroclinic instability.
613	!!
614	!! ** Method : Specific to the offline model. Computes the horizontal
615	!! values from the vertical value
616	!!
617	!! History :
618	!! 9.0 ! 06-03 (O. Aumont) Free form, F90
619	!!----------------------------------------------------------------------
620	!! * Arguments
621	INTEGER, INTENT( in ) :: kt ! ocean time-step inedx
622
623	!! * Local declarations
624	INTEGER :: ji, jj ! dummy loop indices
625	!!----------------------------------------------------------------------
626
627	IF( kt == nit000 ) THEN
628	IF(lwp) WRITE(numout,*)
629	IF(lwp) WRITE(numout,*) 'dta_eiv : eddy induced velocity coefficients'
630	IF(lwp) WRITE(numout,*) '~~~~~~~'
631	ENDIF
632
633	! Average the diffusive coefficient at u- v- points
634	DO jj = 2, jpjm1
635	DO ji = fs_2, fs_jpim1 ! vector opt.
636	aeiu(ji,jj) = .5 * ( aeiw(ji,jj) + aeiw(ji+1,jj ) )
637	aeiv(ji,jj) = .5 * ( aeiw(ji,jj) + aeiw(ji ,jj+1) )
638	END DO
639	END DO
640
641	! lateral boundary condition on aeiu, aeiv
642	CALL lbc_lnk( aeiu, 'U', 1. )
643	CALL lbc_lnk( aeiv, 'V', 1. )
644
645	END SUBROUTINE dta_eiv
646
647	SUBROUTINE tau2wnd( ptaux, ptauy, pwspd )
648	!!---------------------------------------------------------------------
649	!! * ROUTINE sbc_tau2wnd *
650	!!
651	!! ** Purpose : Estimation of wind speed as a function of wind stress
652	!!
653	!! ** Method : \|tau\|=rhoaCd\|U\|^2
654	!!---------------------------------------------------------------------
655	!! * Arguments
656	REAL(wp), DIMENSION(jpi,jpj), INTENT( in ) :: &
657	ptaux, ptauy !: wind stress in i-j direction resp.
658	REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) :: &
659	pwspd !: wind speed
660	REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3
661	REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient
662	REAL(wp) :: ztx, zty, ztau, zcoef ! temporary variables
663	INTEGER :: ji, jj ! dummy indices
664	!!---------------------------------------------------------------------
665	zcoef = 1. / ( zrhoa * zcdrag )
666	!CDIR NOVERRCHK
667	DO jj = 2, jpjm1
668	!CDIR NOVERRCHK
669	DO ji = fs_2, fs_jpim1 ! vector opt.
670	ztx = ptaux(ji,jj) * umask(ji,jj,1) + ptaux(ji-1,jj ) * umask(ji-1,jj ,1)
671	zty = ptauy(ji,jj) * vmask(ji,jj,1) + ptauy(ji ,jj-1) * vmask(ji ,jj-1,1)
672	ztau = 0.5 * SQRT( ztx * ztx + zty * zty )
673	pwspd(ji,jj) = SQRT ( ztau * zcoef ) * tmask(ji,jj,1)
674	END DO
675	END DO
676	CALL lbc_lnk( pwspd(:,:), 'T', 1. )
677
678	END SUBROUTINE tau2wnd
679
680	SUBROUTINE swap_dyn_data
681	!!----------------------------------------------------------------------
682	!! * ROUTINE swap_dyn_data *
683	!!
684	!! ** Purpose : swap array data
685	!!
686	!! History :
687	!! 9.0 ! 07-09 (C. Ethe)
688	!!----------------------------------------------------------------------
689
690
691	! swap from record 2 to 1
692	tdta (:,:,:,1) = tdta (:,:,:,2)
693	sdta (:,:,:,1) = sdta (:,:,:,2)
694	avtdta (:,:,:,1) = avtdta (:,:,:,2)
695	udta (:,:,:,1) = udta (:,:,:,2)
696	vdta (:,:,:,1) = vdta (:,:,:,2)
697	wdta (:,:,:,1) = wdta (:,:,:,2)
698
699	#if defined key_ldfslp
700	uslpdta (:,:,:,1) = uslpdta (:,:,:,2)
701	vslpdta (:,:,:,1) = vslpdta (:,:,:,2)
702	wslpidta(:,:,:,1) = wslpidta(:,:,:,2)
703	wslpjdta(:,:,:,1) = wslpjdta(:,:,:,2)
704	#endif
705	hmlddta(:,:,1) = hmlddta(:,:,2)
706	wspddta(:,:,1) = wspddta(:,:,2)
707	frlddta(:,:,1) = frlddta(:,:,2)
708	empdta (:,:,1) = empdta (:,:,2)
709	qsrdta (:,:,1) = qsrdta (:,:,2)
710
711	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
712	aeiwdta(:,:,1) = aeiwdta(:,:,2)
713	#endif
714
715	#if defined key_degrad
716	ahtudta(:,:,:,1) = ahtudta(:,:,:,2)
717	ahtvdta(:,:,:,1) = ahtvdta(:,:,:,2)
718	ahtwdta(:,:,:,1) = ahtwdta(:,:,:,2)
719	# if defined key_traldf_eiv
720	aeiudta(:,:,:,1) = aeiudta(:,:,:,2)
721	aeivdta(:,:,:,1) = aeivdta(:,:,:,2)
722	aeiwdta(:,:,:,1) = aeiwdta(:,:,:,2)
723	# endif
724	#endif
725
726	END SUBROUTINE swap_dyn_data
727
728	SUBROUTINE assign_dyn_data
729	!!----------------------------------------------------------------------
730	!! * ROUTINE assign_dyn_data *
731	!!
732	!! ** Purpose : Assign dynamical data to the data that have been read
733	!! without time interpolation
734	!!
735	!!----------------------------------------------------------------------
736
737	tn (:,:,:) = tdta (:,:,:,2)
738	sn (:,:,:) = sdta (:,:,:,2)
739	avt(:,:,:) = avtdta(:,:,:,2)
740
741	un (:,:,:) = udta (:,:,:,2)
742	vn (:,:,:) = vdta (:,:,:,2)
743	wn (:,:,:) = wdta (:,:,:,2)
744
745	#if defined key_zdfddm
746	avs(:,:,:) = avtdta (:,:,:,2)
747	#endif
748
749
750	#if defined key_ldfslp
751	uslp (:,:,:) = uslpdta (:,:,:,2)
752	vslp (:,:,:) = vslpdta (:,:,:,2)
753	wslpi(:,:,:) = wslpidta(:,:,:,2)
754	wslpj(:,:,:) = wslpjdta(:,:,:,2)
755	#endif
756
757	hmld(:,:) = hmlddta(:,:,2)
758	wndm(:,:) = wspddta(:,:,2)
759	fr_i(:,:) = frlddta(:,:,2)
760	emp (:,:) = empdta (:,:,2)
761	emps(:,:) = emp(:,:)
762	qsr (:,:) = qsrdta (:,:,2)
763
764	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
765	aeiw(:,:) = aeiwdta(:,:,2)
766	#endif
767
768	#if defined key_degrad
769	ahtu(:,:,:) = ahtudta(:,:,:,2)
770	ahtv(:,:,:) = ahtvdta(:,:,:,2)
771	ahtw(:,:,:) = ahtwdta(:,:,:,2)
772	# if defined key_traldf_eiv
773	aeiu(:,:,:) = aeiudta(:,:,:,2)
774	aeiv(:,:,:) = aeivdta(:,:,:,2)
775	aeiw(:,:,:) = aeiwdta(:,:,:,2)
776	# endif
777
778	#endif
779
780	END SUBROUTINE assign_dyn_data
781
782	SUBROUTINE linear_interp_dyn_data( pweigh )
783	!!----------------------------------------------------------------------
784	!! * ROUTINE linear_interp_dyn_data *
785	!!
786	!! ** Purpose : linear interpolation of data
787	!!
788	!!----------------------------------------------------------------------
789	!! * Argument
790	REAL(wp), INTENT( in ) :: pweigh ! weigh
791
792	!! * Local declarations
793	REAL(wp) :: zweighm1
794	!!----------------------------------------------------------------------
795
796	zweighm1 = 1. - pweigh
797
798	tn (:,:,:) = zweighm1 * tdta (:,:,:,1) + pweigh * tdta (:,:,:,2)
799	sn (:,:,:) = zweighm1 * sdta (:,:,:,1) + pweigh * sdta (:,:,:,2)
800	avt(:,:,:) = zweighm1 * avtdta(:,:,:,1) + pweigh * avtdta(:,:,:,2)
801
802	un (:,:,:) = zweighm1 * udta (:,:,:,1) + pweigh * udta (:,:,:,2)
803	vn (:,:,:) = zweighm1 * vdta (:,:,:,1) + pweigh * vdta (:,:,:,2)
804	wn (:,:,:) = zweighm1 * wdta (:,:,:,1) + pweigh * wdta (:,:,:,2)
805
806	#if defined key_zdfddm
807	avs(:,:,:) = zweighm1 * avtdta (:,:,:,1) + pweigh * avtdta (:,:,:,2)
808	#endif
809
810
811	#if defined key_ldfslp
812	uslp (:,:,:) = zweighm1 * uslpdta (:,:,:,1) + pweigh * uslpdta (:,:,:,2)
813	vslp (:,:,:) = zweighm1 * vslpdta (:,:,:,1) + pweigh * vslpdta (:,:,:,2)
814	wslpi(:,:,:) = zweighm1 * wslpidta(:,:,:,1) + pweigh * wslpidta(:,:,:,2)
815	wslpj(:,:,:) = zweighm1 * wslpjdta(:,:,:,1) + pweigh * wslpjdta(:,:,:,2)
816	#endif
817
818	hmld(:,:) = zweighm1 * hmlddta(:,:,1) + pweigh * hmlddta(:,:,2)
819	wndm(:,:) = zweighm1 * wspddta(:,:,1) + pweigh * wspddta(:,:,2)
820	fr_i(:,:) = zweighm1 * frlddta(:,:,1) + pweigh * frlddta(:,:,2)
821	emp (:,:) = zweighm1 * empdta (:,:,1) + pweigh * empdta (:,:,2)
822	emps(:,:) = emp(:,:)
823	qsr (:,:) = zweighm1 * qsrdta (:,:,1) + pweigh * qsrdta (:,:,2)
824
825	#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv
826	aeiw(:,:) = zweighm1 * aeiwdta(:,:,1) + pweigh * aeiwdta(:,:,2)
827	#endif
828
829	#if defined key_degrad
830	ahtu(:,:,:) = zweighm1 * ahtudta(:,:,:,1) + pweigh * ahtudta(:,:,:,2)
831	ahtv(:,:,:) = zweighm1 * ahtvdta(:,:,:,1) + pweigh * ahtvdta(:,:,:,2)
832	ahtw(:,:,:) = zweighm1 * ahtwdta(:,:,:,1) + pweigh * ahtwdta(:,:,:,2)
833	# if defined key_traldf_eiv
834	aeiu(:,:,:) = zweighm1 * aeiudta(:,:,:,1) + pweigh * aeiudta(:,:,:,2)
835	aeiv(:,:,:) = zweighm1 * aeivdta(:,:,:,1) + pweigh * aeivdta(:,:,:,2)
836	aeiw(:,:,:) = zweighm1 * aeiwdta(:,:,:,1) + pweigh * aeiwdta(:,:,:,2)
837	# endif
838	#endif
839
840	END SUBROUTINE linear_interp_dyn_data
841
842	END MODULE dtadyn

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

Download in other formats: