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 trunk/NEMOGCM/NEMO/OFF_SRC – NEMO

Context Navigation

source: trunk/NEMOGCM/NEMO/OFF_SRC/dtadyn.F90 @ 2528

Last change on this file since 2528 was 2528, checked in by rblod, 13 years ago
Update NEMOGCM from branch nemo_v3_3_beta
Property svn:keywords set to `Id`
File size: 34.0 KB

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

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

Download in other formats: