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.

trdicp.F90 in trunk/NEMO/OPA_SRC/TRD – NEMO

Context Navigation

source: trunk/NEMO/OPA_SRC/TRD/trdicp.F90 @ 1887

Last change on this file since 1887 was 1708, checked in by rblod, 15 years ago
Stability for explicit bottom friction, see ticket #588
Property svn:eol-style set to `native` Property svn:keywords set to `Id`
File size: 33.6 KB

Rev	Line
[215]	1	MODULE trdicp
	2	!!======================================================================
	3	!! * MODULE trdicp *
	4	!! Ocean diagnostics: ocean tracers and dynamic trends
	5	!!=====================================================================
[503]	6	!! History : ! 91-12 (G. Madec)
	7	!! ! 92-06 (M. Imbard) add time step frequency
	8	!! ! 96-01 (G. Madec) terrain following coordinates
	9	!! 8.5 ! 02-06 (G. Madec) F90: Free form and module
	10	!! 9.0 ! 04-08 (C. Talandier) New trends organization
	11	!!----------------------------------------------------------------------
[215]	12	#if defined key_trdtra \|\| defined key_trddyn \|\| defined key_esopa
	13	!!----------------------------------------------------------------------
	14	!! 'key_trdtra' or active tracers trends diagnostics
	15	!! 'key_trddyn' momentum trends diagnostics
	16	!!----------------------------------------------------------------------
	17	!!----------------------------------------------------------------------
[503]	18	!! trd_icp : compute the basin averaged properties for tra/dyn
[215]	19	!! trd_dwr : print dynmaic trends in ocean.output file
	20	!! trd_twr : print tracers trends in ocean.output file
	21	!! trd_icp_init : initialization step
	22	!!----------------------------------------------------------------------
	23	USE oce ! ocean dynamics and tracers variables
	24	USE dom_oce ! ocean space and time domain variables
	25	USE trdmod_oce ! ocean variables trends
	26	USE ldftra_oce ! ocean active tracers: lateral physics
	27	USE ldfdyn_oce ! ocean dynamics: lateral physics
	28	USE zdf_oce ! ocean vertical physics
	29	USE in_out_manager ! I/O manager
	30	USE lib_mpp ! distibuted memory computing library
	31	USE eosbn2 ! equation of state
	32	USE phycst ! physical constants
	33
	34	IMPLICIT NONE
	35	PRIVATE
	36
[503]	37	INTERFACE trd_icp
[215]	38	MODULE PROCEDURE trd_2d, trd_3d
	39	END INTERFACE
	40
[503]	41	PUBLIC trd_icp ! called by trdmod.F90
	42	PUBLIC trd_dwr ! called by step.F90
	43	PUBLIC trd_twr ! called by step.F90
	44	PUBLIC trd_icp_init ! called by opa.F90
[215]	45
	46	!! * Substitutions
	47	# include "domzgr_substitute.h90"
	48	# include "vectopt_loop_substitute.h90"
	49	!!----------------------------------------------------------------------
[247]	50	!! OPA 9.0 , LOCEAN-IPSL (2005)
[1152]	51	!! $Id$
[503]	52	!! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
[215]	53	!!----------------------------------------------------------------------
	54
	55	CONTAINS
	56
[503]	57	SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype, clpas )
[215]	58	!!---------------------------------------------------------------------
	59	!! * ROUTINE trd_2d *
	60	!!
	61	!! ** Purpose : verify the basin averaged properties of tracers and/or
[1601]	62	!! momentum equations at every time step frequency nn_trd.
[503]	63	!!----------------------------------------------------------------------
	64	REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: ptrd2dx ! Temperature or U trend
	65	REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: ptrd2dy ! Salinity or V trend
	66	INTEGER , INTENT(in ) :: ktrd ! tracer trend index
	67	CHARACTER(len=3) , INTENT(in ) :: ctype ! momentum ('DYN') or tracers ('TRA') trends
	68	CHARACTER(len=3) , INTENT(in ), OPTIONAL :: clpas ! number of passage
[215]	69	!!
[503]	70	INTEGER :: ji, jj ! loop indices
	71	CHARACTER(len=3) :: cpas ! number of passage
	72	REAL(wp) :: zmsku, zbtu, zbt ! temporary scalars
	73	REAL(wp) :: zmskv, zbtv ! " "
[215]	74	!!----------------------------------------------------------------------
	75
[503]	76	! Control of optional arguments
	77	cpas = 'fst'
	78	IF( PRESENT(clpas) ) cpas = clpas
[215]	79
[503]	80	! 1. Mask trends
	81	! --------------
[215]	82
[503]	83	SELECT CASE( ctype )
	84	!
	85	CASE( 'DYN' ) ! Momentum
[215]	86	DO jj = 1, jpjm1
	87	DO ji = 1, jpim1
	88	zmsku = tmask_i(ji+1,jj ) * tmask_i(ji,jj) * umask(ji,jj,1)
	89	zmskv = tmask_i(ji ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,1)
	90	ptrd2dx(ji,jj) = ptrd2dx(ji,jj) * zmsku
	91	ptrd2dy(ji,jj) = ptrd2dy(ji,jj) * zmskv
	92	END DO
	93	END DO
	94	ptrd2dx(jpi, : ) = 0.e0 ; ptrd2dy(jpi, : ) = 0.e0
	95	ptrd2dx( : ,jpj) = 0.e0 ; ptrd2dy( : ,jpj) = 0.e0
[503]	96	!
	97	CASE( 'TRA' ) ! Tracers
[215]	98	ptrd2dx(:,:) = ptrd2dx(:,:) * tmask_i(:,:)
	99	ptrd2dy(:,:) = ptrd2dy(:,:) * tmask_i(:,:)
[503]	100	!
[215]	101	END SELECT
	102
[503]	103	! 2. Basin averaged tracer/momentum trends
	104	! ----------------------------------------
[215]	105
[503]	106	SELECT CASE( ctype )
	107	!
	108	CASE( 'DYN' ) ! Momentum
[215]	109	umo(ktrd) = 0.e0
	110	vmo(ktrd) = 0.e0
[503]	111	!
	112	SELECT CASE( ktrd )
	113	!
	114	CASE( jpdyn_trd_swf ) ! surface forcing
[215]	115	DO jj = 1, jpj
	116	DO ji = 1, jpi
	117	umo(ktrd) = umo(ktrd) + ptrd2dx(ji,jj) * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1)
	118	vmo(ktrd) = vmo(ktrd) + ptrd2dy(ji,jj) * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1)
	119	END DO
	120	END DO
[503]	121	!
[215]	122	END SELECT
[503]	123	!
	124	CASE( 'TRA' ) ! Tracers
	125	IF( cpas == 'fst' ) THEN
	126	tmo(ktrd) = 0.e0
	127	smo(ktrd) = 0.e0
	128	ENDIF
[215]	129	DO jj = 1, jpj
	130	DO ji = 1, jpi
	131	zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
	132	tmo(ktrd) = tmo(ktrd) + ptrd2dx(ji,jj) * zbt
	133	smo(ktrd) = smo(ktrd) + ptrd2dy(ji,jj) * zbt
	134	END DO
	135	END DO
[503]	136	!
[215]	137	END SELECT
	138
[503]	139	! 3. Basin averaged tracer/momentum square trends
	140	! ----------------------------------------------
[215]	141	! c a u t i o n: field now
	142
[503]	143	SELECT CASE( ctype )
	144	!
	145	CASE( 'DYN' ) ! Momentum
[215]	146	hke(ktrd) = 0.e0
	147	DO jj = 1, jpj
	148	DO ji = 1, jpi
	149	zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1)
	150	zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1)
	151	hke(ktrd) = hke(ktrd) &
	152	& + un(ji,jj,1) * ptrd2dx(ji,jj) * zbtu &
	153	& + vn(ji,jj,1) * ptrd2dy(ji,jj) * zbtv
	154	END DO
	155	END DO
[503]	156	!
	157	CASE( 'TRA' ) ! Tracers
	158	IF( cpas == 'fst' ) THEN
	159	t2(ktrd) = 0.e0
	160	s2(ktrd) = 0.e0
	161	ENDIF
[215]	162	DO jj = 1, jpj
	163	DO ji = 1, jpi
	164	zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
	165	t2(ktrd) = t2(ktrd) + ptrd2dx(ji,jj) * zbt * tn(ji,jj,1)
	166	s2(ktrd) = s2(ktrd) + ptrd2dy(ji,jj) * zbt * sn(ji,jj,1)
	167	END DO
	168	END DO
[503]	169	!
[215]	170	END SELECT
[503]	171	!
[215]	172	END SUBROUTINE trd_2d
	173
	174
[503]	175	SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd, ctype, clpas )
[215]	176	!!---------------------------------------------------------------------
	177	!! * ROUTINE trd_3d *
	178	!!
	179	!! ** Purpose : verify the basin averaged properties of tracers and/or
[1601]	180	!! momentum equations at every time step frequency nn_trd.
[503]	181	!!----------------------------------------------------------------------
	182	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: ptrd3dx ! Temperature or U trend
	183	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: ptrd3dy ! Salinity or V trend
	184	INTEGER, INTENT(in ) :: ktrd ! momentum or tracer trend index
	185	CHARACTER(len=3), INTENT(in ) :: ctype ! momentum ('DYN') or tracers ('TRA') trends
	186	CHARACTER(len=3), INTENT(in ), OPTIONAL :: clpas ! number of passage
[215]	187	!!
	188	INTEGER :: ji, jj, jk
[503]	189	CHARACTER(len=3) :: cpas ! number of passage
	190	REAL(wp) :: zbt, zbtu, zbtv, zmsku, zmskv ! temporary scalars
[215]	191	!!----------------------------------------------------------------------
	192
[503]	193	! Control of optional arguments
	194	cpas = 'fst'
	195	IF( PRESENT(clpas) ) cpas = clpas
[215]	196
[503]	197	! 1. Mask the trends
	198	! ------------------
[215]	199
[503]	200	SELECT CASE( ctype )
	201	!
	202	CASE( 'DYN' ) ! Momentum
[215]	203	DO jk = 1, jpk
	204	DO jj = 1, jpjm1
	205	DO ji = 1, jpim1
	206	zmsku = tmask_i(ji+1,jj ) * tmask_i(ji,jj) * umask(ji,jj,jk)
	207	zmskv = tmask_i(ji ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
	208	ptrd3dx(ji,jj,jk) = ptrd3dx(ji,jj,jk) * zmsku
	209	ptrd3dy(ji,jj,jk) = ptrd3dy(ji,jj,jk) * zmskv
[503]	210	END DO
	211	END DO
	212	END DO
[215]	213	ptrd3dx(jpi, : ,:) = 0.e0 ; ptrd3dy(jpi, : ,:) = 0.e0
	214	ptrd3dx( : ,jpj,:) = 0.e0 ; ptrd3dy( : ,jpj,:) = 0.e0
[503]	215	!
	216	CASE( 'TRA' ) ! Tracers
[215]	217	DO jk = 1, jpk
	218	ptrd3dx(:,:,jk) = ptrd3dx(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
	219	ptrd3dy(:,:,jk) = ptrd3dy(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
[503]	220	END DO
	221	!
[215]	222	END SELECT
	223
[503]	224	! 2. Basin averaged tracer/momentum trends
	225	! ----------------------------------------
[215]	226
[503]	227	SELECT CASE( ctype )
	228	!
	229	CASE( 'DYN' ) ! Momentum
[215]	230	umo(ktrd) = 0.e0
	231	vmo(ktrd) = 0.e0
	232	DO jk = 1, jpk
	233	DO jj = 1, jpj
	234	DO ji = 1, jpi
	235	zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
	236	zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
	237	umo(ktrd) = umo(ktrd) + ptrd3dx(ji,jj,jk) * zbtu
	238	vmo(ktrd) = vmo(ktrd) + ptrd3dy(ji,jj,jk) * zbtv
	239	END DO
	240	END DO
	241	END DO
[503]	242	!
	243	CASE( 'TRA' ) ! Tracers
	244	IF( cpas == 'fst' ) THEN
	245	tmo(ktrd) = 0.e0
	246	smo(ktrd) = 0.e0
	247	ENDIF
[215]	248	DO jk = 1, jpkm1
	249	DO jj = 1, jpj
	250	DO ji = 1, jpi
	251	zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	252	tmo(ktrd) = tmo(ktrd) + ptrd3dx(ji,jj,jk) * zbt
	253	smo(ktrd) = smo(ktrd) + ptrd3dy(ji,jj,jk) * zbt
	254	END DO
	255	END DO
	256	END DO
[503]	257	!
[215]	258	END SELECT
	259
[503]	260	! 3. Basin averaged tracer/momentum square trends
	261	! -----------------------------------------------
[215]	262	! c a u t i o n: field now
	263
[503]	264	SELECT CASE( ctype )
	265	!
	266	CASE( 'DYN' ) ! Momentum
[215]	267	hke(ktrd) = 0.e0
	268	DO jk = 1, jpk
	269	DO jj = 1, jpj
	270	DO ji = 1, jpi
	271	zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
	272	zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
	273	hke(ktrd) = hke(ktrd) &
	274	& + un(ji,jj,jk) * ptrd3dx(ji,jj,jk) * zbtu &
	275	& + vn(ji,jj,jk) * ptrd3dy(ji,jj,jk) * zbtv
	276	END DO
	277	END DO
	278	END DO
[503]	279	!
	280	CASE( 'TRA' ) ! Tracers
	281	IF( cpas == 'fst' ) THEN
	282	t2(ktrd) = 0.e0
	283	s2(ktrd) = 0.e0
	284	ENDIF
[215]	285	DO jk = 1, jpk
	286	DO jj = 1, jpj
	287	DO ji = 1, jpi
	288	zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	289	t2(ktrd) = t2(ktrd) + ptrd3dx(ji,jj,jk) * zbt * tn(ji,jj,jk)
	290	s2(ktrd) = s2(ktrd) + ptrd3dy(ji,jj,jk) * zbt * sn(ji,jj,jk)
	291	END DO
	292	END DO
	293	END DO
[503]	294	!
[215]	295	END SELECT
[503]	296	!
[215]	297	END SUBROUTINE trd_3d
	298
	299
	300
	301	SUBROUTINE trd_icp_init
	302	!!---------------------------------------------------------------------
	303	!! * ROUTINE trd_icp_init *
	304	!!
[503]	305	!! ** Purpose : Read the namtrd namelist
[215]	306	!!----------------------------------------------------------------------
[503]	307	INTEGER :: ji, jj, jk
[215]	308	REAL(wp) :: zmskt
	309	#if defined key_trddyn
[503]	310	REAL(wp) :: zmsku, zmskv
[215]	311	#endif
	312	!!----------------------------------------------------------------------
	313
	314	IF(lwp) THEN
	315	WRITE(numout,*)
	316	WRITE(numout,*) 'trd_icp_init : integral constraints properties trends'
	317	WRITE(numout,*) '~~~~~~~~~~~~~'
	318	ENDIF
	319
	320	! Total volume at t-points:
	321	tvolt = 0.e0
	322	DO jk = 1, jpkm1
	323	DO jj = 2, jpjm1
	324	DO ji = fs_2, fs_jpim1 ! vector opt.
	325	zmskt = tmask(ji,jj,jk) * tmask_i(ji,jj)
	326	tvolt = tvolt + zmskt * e1t(ji,jj) e2t(ji,jj) fse3t(ji,jj,jk)
	327	END DO
	328	END DO
	329	END DO
	330	IF( lk_mpp ) CALL mpp_sum( tvolt ) ! sum over the global domain
	331
[503]	332	IF(lwp) WRITE(numout,*) ' total ocean volume at T-point tvolt = ',tvolt
[215]	333
	334	#if defined key_trddyn
	335	! Initialization of potential to kinetic energy conversion
	336	rpktrd = 0.e0
	337
	338	! Total volume at u-, v- points:
	339	tvolu = 0.e0
	340	tvolv = 0.e0
	341
	342	DO jk = 1, jpk
	343	DO jj = 2, jpjm1
	344	DO ji = fs_2, fs_jpim1 ! vector opt.
	345	zmsku = tmask_i(ji+1,jj ) * tmask_i(ji,jj) * umask(ji,jj,jk)
	346	zmskv = tmask_i(ji ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
	347	tvolu = tvolu + zmsku * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
	348	tvolv = tvolv + zmskv * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
	349	END DO
	350	END DO
	351	END DO
	352	IF( lk_mpp ) CALL mpp_sum( tvolu ) ! sums over the global domain
	353	IF( lk_mpp ) CALL mpp_sum( tvolv )
	354
	355	IF(lwp) THEN
[503]	356	WRITE(numout,*) ' total ocean volume at U-point tvolu = ',tvolu
	357	WRITE(numout,*) ' total ocean volume at V-point tvolv = ',tvolv
[215]	358	ENDIF
	359	#endif
[503]	360	!
[215]	361	END SUBROUTINE trd_icp_init
	362
	363
	364	SUBROUTINE trd_dwr( kt )
	365	!!---------------------------------------------------------------------
	366	!! * ROUTINE trd_dwr *
	367	!!
	368	!! ** Purpose : write dynamic trends in ocean.output
[503]	369	!!----------------------------------------------------------------------
	370	INTEGER, INTENT(in) :: kt ! ocean time-step index
[215]	371	!!
[503]	372	INTEGER :: ji, jj, jk
	373	REAL(wp) :: ze1e2w, zcof, zbe1ru, zbe2rv, zbtr, ztz, zth ! " scalars
	374	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zkepe, zkx, zky, zkz ! temporary arrays
[215]	375	!!----------------------------------------------------------------------
	376
	377	! I. Momentum trends
	378	! -------------------
	379
[1601]	380	IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
[215]	381
	382	! I.1 Conversion potential energy - kinetic energy
	383	! --------------------------------------------------
	384	! c a u t i o n here, trends are computed at kt+1 (now , but after the swap)
	385
[503]	386	zkx(:,:,:) = 0.e0
	387	zky(:,:,:) = 0.e0
	388	zkz(:,:,:) = 0.e0
[215]	389	zkepe(:,:,:) = 0.e0
	390
	391	CALL eos( tn, sn, rhd, rhop ) ! now potential and in situ densities
	392
[503]	393	! Density flux at w-point
[215]	394	DO jk = 2, jpk
	395	DO jj = 1, jpj
	396	DO ji = 1, jpi
	397	ze1e2w = 0.5 * e1t(ji,jj) * e2t(ji,jj) * wn(ji,jj,jk) * tmask_i(ji,jj)
	398	zkz(ji,jj,jk) = ze1e2w / rau0 * ( rhop(ji,jj,jk) + rhop(ji,jj,jk-1) )
	399	END DO
	400	END DO
	401	END DO
[503]	402	zkz(:,:,1) = 0.e0
[215]	403
	404	! Density flux at u and v-points
	405	DO jk = 1, jpk
	406	DO jj = 1, jpjm1
	407	DO ji = 1, jpim1
	408	zcof = 0.5 / rau0
	409	zbe1ru = zcof * e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk)
	410	zbe2rv = zcof * e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk)
	411	zkx(ji,jj,jk) = zbe1ru * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) )
	412	zky(ji,jj,jk) = zbe2rv * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) )
	413	END DO
	414	END DO
	415	END DO
	416
	417	! Density flux divergence at t-point
	418	DO jk = 1, jpkm1
	419	DO jj = 2, jpjm1
	420	DO ji = 2, jpim1
	421	zbtr = 1. / ( e1t(ji,jj)e2t(ji,jj)fse3t(ji,jj,jk) )
	422	ztz = - zbtr * ( zkz(ji,jj,jk) - zkz(ji,jj,jk+1) )
	423	zth = - zbtr * ( ( zkx(ji,jj,jk) - zkx(ji-1,jj,jk) ) &
	424	& + ( zky(ji,jj,jk) - zky(ji,jj-1,jk) ) )
	425	zkepe(ji,jj,jk) = (zth + ztz) * tmask(ji,jj,jk) * tmask_i(ji,jj)
	426	END DO
	427	END DO
	428	END DO
	429	zkepe( : , : ,jpk) = 0.e0
	430	zkepe( : ,jpj, : ) = 0.e0
	431	zkepe(jpi, : , : ) = 0.e0
	432
	433	! I.2 Basin averaged kinetic energy trend
	434	! ----------------------------------------
	435	peke = 0.e0
	436	DO jk = 1,jpk
	437	DO jj = 1, jpj
	438	DO ji = 1, jpi
	439	peke = peke + zkepe(ji,jj,jk) * grav * fsdept(ji,jj,jk) &
	440	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	441	END DO
	442	END DO
	443	END DO
	444
	445	! I.3 Sums over the global domain
	446	! ---------------------------------
	447	IF( lk_mpp ) THEN
[503]	448	CALL mpp_sum( peke )
	449	CALL mpp_sum( umo , jptot_dyn )
	450	CALL mpp_sum( vmo , jptot_dyn )
	451	CALL mpp_sum( hke , jptot_dyn )
	452	ENDIF
[215]	453
	454	! I.2 Print dynamic trends in the ocean.output file
	455	! --------------------------------------------------
	456
	457	IF(lwp) THEN
	458	WRITE (numout,*)
	459	WRITE (numout,*)
	460	WRITE (numout,9500) kt
[503]	461	WRITE (numout,9501) umo(jpicpd_hpg) / tvolu, vmo(jpicpd_hpg) / tvolv
	462	WRITE (numout,9502) umo(jpicpd_keg) / tvolu, vmo(jpicpd_keg) / tvolv
	463	WRITE (numout,9503) umo(jpicpd_rvo) / tvolu, vmo(jpicpd_rvo) / tvolv
	464	WRITE (numout,9504) umo(jpicpd_pvo) / tvolu, vmo(jpicpd_pvo) / tvolv
	465	WRITE (numout,9505) umo(jpicpd_ldf) / tvolu, vmo(jpicpd_ldf) / tvolv
[1129]	466	WRITE (numout,9506) umo(jpicpd_had) / tvolu, vmo(jpicpd_had) / tvolv
	467	WRITE (numout,9507) umo(jpicpd_zad) / tvolu, vmo(jpicpd_zad) / tvolv
	468	WRITE (numout,9508) umo(jpicpd_zdf) / tvolu, vmo(jpicpd_zdf) / tvolv
	469	WRITE (numout,9509) umo(jpicpd_spg) / tvolu, vmo(jpicpd_spg) / tvolv
	470	WRITE (numout,9510) umo(jpicpd_swf) / tvolu, vmo(jpicpd_swf) / tvolv
	471	WRITE (numout,9511) umo(jpicpd_dat) / tvolu, vmo(jpicpd_dat) / tvolv
	472	WRITE (numout,9512) umo(jpicpd_bfr) / tvolu, vmo(jpicpd_bfr) / tvolv
	473	WRITE (numout,9513)
	474	WRITE (numout,9514) &
[503]	475	& ( umo(jpicpd_hpg) + umo(jpicpd_keg) + umo(jpicpd_rvo) + umo(jpicpd_pvo) + umo(jpicpd_ldf) &
[1129]	476	& + umo(jpicpd_had) + umo(jpicpd_zad) + umo(jpicpd_zdf) + umo(jpicpd_spg) + umo(jpicpd_dat) &
	477	& + umo(jpicpd_swf) + umo(jpicpd_bfr) ) / tvolu, &
[503]	478	& ( vmo(jpicpd_hpg) + vmo(jpicpd_keg) + vmo(jpicpd_rvo) + vmo(jpicpd_pvo) + vmo(jpicpd_ldf) &
[1129]	479	& + vmo(jpicpd_had) + vmo(jpicpd_zad) + vmo(jpicpd_zdf) + vmo(jpicpd_spg) + vmo(jpicpd_dat) &
	480	& + vmo(jpicpd_swf) + vmo(jpicpd_bfr) ) / tvolv
[215]	481	ENDIF
	482
	483	9500 FORMAT(' momentum trend at it= ', i6, ' :', /' ==============================')
	484	9501 FORMAT(' pressure gradient u= ', e20.13, ' v= ', e20.13)
	485	9502 FORMAT(' ke gradient u= ', e20.13, ' v= ', e20.13)
	486	9503 FORMAT(' relative vorticity term u= ', e20.13, ' v= ', e20.13)
	487	9504 FORMAT(' coriolis term u= ', e20.13, ' v= ', e20.13)
	488	9505 FORMAT(' horizontal diffusion u= ', e20.13, ' v= ', e20.13)
[1129]	489	9506 FORMAT(' horizontal advection u= ', e20.13, ' v= ', e20.13)
	490	9507 FORMAT(' vertical advection u= ', e20.13, ' v= ', e20.13)
	491	9508 FORMAT(' vertical diffusion u= ', e20.13, ' v= ', e20.13)
	492	9509 FORMAT(' surface pressure gradient u= ', e20.13, ' v= ', e20.13)
	493	9510 FORMAT(' surface wind forcing u= ', e20.13, ' v= ', e20.13)
	494	9511 FORMAT(' dampimg term u= ', e20.13, ' v= ', e20.13)
	495	9512 FORMAT(' bottom flux u= ', e20.13, ' v= ', e20.13)
	496	9513 FORMAT(' -----------------------------------------------------------------------------')
	497	9514 FORMAT(' total trend u= ', e20.13, ' v= ', e20.13)
[215]	498
	499	IF(lwp) THEN
	500	WRITE (numout,*)
	501	WRITE (numout,*)
	502	WRITE (numout,9520) kt
[503]	503	WRITE (numout,9521) hke(jpicpd_hpg) / tvolt
	504	WRITE (numout,9522) hke(jpicpd_keg) / tvolt
	505	WRITE (numout,9523) hke(jpicpd_rvo) / tvolt
	506	WRITE (numout,9524) hke(jpicpd_pvo) / tvolt
	507	WRITE (numout,9525) hke(jpicpd_ldf) / tvolt
[1129]	508	WRITE (numout,9526) hke(jpicpd_had) / tvolt
	509	WRITE (numout,9527) hke(jpicpd_zad) / tvolt
	510	WRITE (numout,9528) hke(jpicpd_zdf) / tvolt
	511	WRITE (numout,9529) hke(jpicpd_spg) / tvolt
	512	WRITE (numout,9530) hke(jpicpd_swf) / tvolt
	513	WRITE (numout,9531) hke(jpicpd_dat) / tvolt
[1708]	514	WRITE (numout,9532) hke(jpicpd_bfr) / tvolt
	515	WRITE (numout,9533)
	516	WRITE (numout,9534) &
[503]	517	& ( hke(jpicpd_hpg) + hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_pvo) + hke(jpicpd_ldf) &
[1129]	518	& + hke(jpicpd_had) + hke(jpicpd_zad) + hke(jpicpd_zdf) + hke(jpicpd_spg) + hke(jpicpd_dat) &
[1708]	519	& + hke(jpicpd_swf) + hke(jpicpd_bfr) ) / tvolt
[215]	520	ENDIF
	521
	522	9520 FORMAT(' kinetic energy trend at it= ', i6, ' :', /' ====================================')
	523	9521 FORMAT(' pressure gradient u2= ', e20.13)
	524	9522 FORMAT(' ke gradient u2= ', e20.13)
	525	9523 FORMAT(' relative vorticity term u2= ', e20.13)
	526	9524 FORMAT(' coriolis term u2= ', e20.13)
	527	9525 FORMAT(' horizontal diffusion u2= ', e20.13)
[1129]	528	9526 FORMAT(' horizontal advection u2= ', e20.13)
	529	9527 FORMAT(' vertical advection u2= ', e20.13)
	530	9528 FORMAT(' vertical diffusion u2= ', e20.13)
	531	9529 FORMAT(' surface pressure gradient u2= ', e20.13)
	532	9530 FORMAT(' surface wind forcing u2= ', e20.13)
	533	9531 FORMAT(' dampimg term u2= ', e20.13)
[1708]	534	9532 FORMAT(' bottom flux u2= ', e20.13)
	535	9533 FORMAT(' --------------------------------------------------')
	536	9534 FORMAT(' total trend u2= ', e20.13)
[215]	537
	538	IF(lwp) THEN
	539	WRITE (numout,*)
	540	WRITE (numout,*)
	541	WRITE (numout,9540) kt
[1129]	542	WRITE (numout,9541) ( hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
	543	WRITE (numout,9542) ( hke(jpicpd_keg) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
[503]	544	WRITE (numout,9543) ( hke(jpicpd_pvo) ) / tvolt
	545	WRITE (numout,9544) ( hke(jpicpd_rvo) ) / tvolt
	546	WRITE (numout,9545) ( hke(jpicpd_spg) ) / tvolt
	547	WRITE (numout,9546) ( hke(jpicpd_ldf) ) / tvolt
	548	WRITE (numout,9547) ( hke(jpicpd_zdf) ) / tvolt
	549	WRITE (numout,9548) ( hke(jpicpd_hpg) ) / tvolt, rpktrd / tvolt
	550	WRITE (numout,*)
	551	WRITE (numout,*)
[215]	552	ENDIF
	553
	554	9540 FORMAT(' energetic consistency at it= ', i6, ' :', /' =========================================')
	555	9541 FORMAT(' 0 = non linear term(true if key_vorenergy or key_combined): ', e20.13)
[1129]	556	9542 FORMAT(' 0 = ke gradient + horizontal + vertical advection : ', e20.13)
[215]	557	9543 FORMAT(' 0 = coriolis term (true if key_vorenergy or key_combined): ', e20.13)
[503]	558	9544 FORMAT(' 0 = uh.( rot(u) x uh ) (true if enstrophy conser.) : ', e20.13)
	559	9545 FORMAT(' 0 = surface pressure gradient : ', e20.13)
	560	9546 FORMAT(' 0 > horizontal diffusion : ', e20.13)
	561	9547 FORMAT(' 0 > vertical diffusion : ', e20.13)
	562	9548 FORMAT(' pressure gradient u2 = - 1/rau0 u.dz(rhop) : ', e20.13, ' u.dz(rhop) =', e20.13)
	563	!
[215]	564	! Save potential to kinetic energy conversion for next time step
	565	rpktrd = peke
[503]	566	!
[215]	567	ENDIF
[503]	568	!
[215]	569	END SUBROUTINE trd_dwr
	570
	571
	572	SUBROUTINE trd_twr( kt )
	573	!!---------------------------------------------------------------------
	574	!! * ROUTINE trd_twr *
	575	!!
	576	!! ** Purpose : write active tracers trends in ocean.output
	577	!!----------------------------------------------------------------------
[503]	578	INTEGER, INTENT(in) :: kt ! ocean time-step index
[215]	579	!!----------------------------------------------------------------------
	580
	581	! I. Tracers trends
	582	! -----------------
	583
[1601]	584	IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
[215]	585
	586	! I.1 Sums over the global domain
	587	! -------------------------------
	588	IF( lk_mpp ) THEN
[503]	589	CALL mpp_sum( tmo, jptot_tra )
	590	CALL mpp_sum( smo, jptot_tra )
	591	CALL mpp_sum( t2 , jptot_tra )
	592	CALL mpp_sum( s2 , jptot_tra )
[215]	593	ENDIF
	594
	595	! I.2 Print tracers trends in the ocean.output file
	596	! --------------------------------------------------
	597
	598	IF(lwp) THEN
	599	WRITE (numout,*)
	600	WRITE (numout,*)
	601	WRITE (numout,9400) kt
[503]	602	WRITE (numout,9401) (tmo(jpicpt_xad)+tmo(jpicpt_yad))/ tvolt, (smo(jpicpt_xad)+smo(jpicpt_yad))/ tvolt
	603	WRITE (numout,9402) tmo(jpicpt_zad) / tvolt, smo(jpicpt_zad) / tvolt
	604	WRITE (numout,9403) tmo(jpicpt_ldf) / tvolt, smo(jpicpt_ldf) / tvolt
	605	WRITE (numout,9404) tmo(jpicpt_zdf) / tvolt, smo(jpicpt_zdf) / tvolt
	606	WRITE (numout,9405) tmo(jpicpt_npc) / tvolt, smo(jpicpt_npc) / tvolt
	607	WRITE (numout,9406) tmo(jpicpt_dmp) / tvolt, smo(jpicpt_dmp) / tvolt
	608	WRITE (numout,9407) tmo(jpicpt_qsr) / tvolt
	609	WRITE (numout,9408) tmo(jpicpt_nsr) / tvolt, smo(jpicpt_nsr) / tvolt
	610	WRITE (numout,9409)
	611	WRITE (numout,9410) ( tmo(jpicpt_xad) + tmo(jpicpt_yad) + tmo(jpicpt_zad) + tmo(jpicpt_ldf) + tmo(jpicpt_zdf) &
	612	& + tmo(jpicpt_npc) + tmo(jpicpt_dmp) + tmo(jpicpt_qsr) + tmo(jpicpt_nsr) ) / tvolt, &
	613	& ( smo(jpicpt_xad) + smo(jpicpt_yad) + smo(jpicpt_zad) + smo(jpicpt_ldf) + smo(jpicpt_zdf) &
	614	& + smo(jpicpt_npc) + smo(jpicpt_dmp) + smo(jpicpt_nsr) ) / tvolt
[215]	615	ENDIF
	616
	617	9400 FORMAT(' tracer trend at it= ',i6,' : temperature', &
	618	' salinity',/' ============================')
	619	9401 FORMAT(' horizontal advection ',e20.13,' ',e20.13)
	620	9402 FORMAT(' vertical advection ',e20.13,' ',e20.13)
	621	9403 FORMAT(' horizontal diffusion ',e20.13,' ',e20.13)
	622	9404 FORMAT(' vertical diffusion ',e20.13,' ',e20.13)
[503]	623	9405 FORMAT(' static instability mixing ',e20.13,' ',e20.13)
[215]	624	9406 FORMAT(' damping term ',e20.13,' ',e20.13)
[503]	625	9407 FORMAT(' penetrative qsr ',e20.13)
	626	9408 FORMAT(' non solar radiation ',e20.13,' ',e20.13)
[215]	627	9409 FORMAT(' -------------------------------------------------------------------------')
	628	9410 FORMAT(' total trend ',e20.13,' ',e20.13)
	629
	630
	631	IF(lwp) THEN
	632	WRITE (numout,*)
	633	WRITE (numout,*)
	634	WRITE (numout,9420) kt
[503]	635	WRITE (numout,9421) ( t2(jpicpt_xad)+t2(jpicpt_yad) )/ tvolt, ( s2(jpicpt_xad)+s2(jpicpt_yad) )/ tvolt
	636	WRITE (numout,9422) t2(jpicpt_zad) / tvolt, s2(jpicpt_zad) / tvolt
	637	WRITE (numout,9423) t2(jpicpt_ldf) / tvolt, s2(jpicpt_ldf) / tvolt
	638	WRITE (numout,9424) t2(jpicpt_zdf) / tvolt, s2(jpicpt_zdf) / tvolt
	639	WRITE (numout,9425) t2(jpicpt_npc) / tvolt, s2(jpicpt_npc) / tvolt
	640	WRITE (numout,9426) t2(jpicpt_dmp) / tvolt, s2(jpicpt_dmp) / tvolt
	641	WRITE (numout,9427) t2(jpicpt_qsr) / tvolt
	642	WRITE (numout,9428) t2(jpicpt_nsr) / tvolt, s2(jpicpt_nsr) / tvolt
[215]	643	WRITE (numout,9429)
[503]	644	WRITE (numout,9430) ( t2(jpicpt_xad) + t2(jpicpt_yad) + t2(jpicpt_zad) + t2(jpicpt_ldf) + t2(jpicpt_zdf) &
	645	& + t2(jpicpt_npc) + t2(jpicpt_dmp) + t2(jpicpt_qsr) + t2(jpicpt_nsr) ) / tvolt, &
	646	& ( s2(jpicpt_xad) + s2(jpicpt_yad) + s2(jpicpt_zad) + s2(jpicpt_ldf) + s2(jpicpt_zdf) &
	647	& + s2(jpicpt_npc) + s2(jpicpt_dmp) + s2(jpicpt_nsr) ) / tvolt
[215]	648	ENDIF
	649
	650	9420 FORMAT(' tracer**2 trend at it= ', i6, ' : temperature', &
	651	' salinity', /, ' ===============================')
	652	9421 FORMAT(' horizontal advection * t ', e20.13, ' ', e20.13)
	653	9422 FORMAT(' vertical advection * t ', e20.13, ' ', e20.13)
	654	9423 FORMAT(' horizontal diffusion * t ', e20.13, ' ', e20.13)
	655	9424 FORMAT(' vertical diffusion * t ', e20.13, ' ', e20.13)
[503]	656	9425 FORMAT(' static instability mixing * t ', e20.13, ' ', e20.13)
[215]	657	9426 FORMAT(' damping term * t ', e20.13, ' ', e20.13)
[503]	658	9427 FORMAT(' penetrative qsr * t ', e20.13)
	659	9428 FORMAT(' non solar radiation * t ', e20.13, ' ', e20.13)
[215]	660	9429 FORMAT(' -----------------------------------------------------------------------------')
	661	9430 FORMAT(' total trend t = ', e20.13, ' s = ', e20.13)
	662
	663
	664	IF(lwp) THEN
	665	WRITE (numout,*)
	666	WRITE (numout,*)
	667	WRITE (numout,9440) kt
[503]	668	WRITE (numout,9441) ( tmo(jpicpt_xad)+tmo(jpicpt_yad)+tmo(jpicpt_zad) )/tvolt, &
	669	& ( smo(jpicpt_xad)+smo(jpicpt_yad)+smo(jpicpt_zad) )/tvolt
	670	WRITE (numout,9442) tmo(jpicpt_zl1)/tvolt, smo(jpicpt_zl1)/tvolt
	671	WRITE (numout,9443) tmo(jpicpt_ldf)/tvolt, smo(jpicpt_ldf)/tvolt
	672	WRITE (numout,9444) tmo(jpicpt_zdf)/tvolt, smo(jpicpt_zdf)/tvolt
	673	WRITE (numout,9445) tmo(jpicpt_npc)/tvolt, smo(jpicpt_npc)/tvolt
	674	WRITE (numout,9446) ( t2(jpicpt_xad)+t2(jpicpt_yad)+t2(jpicpt_zad) )/tvolt, &
	675	& ( s2(jpicpt_xad)+s2(jpicpt_yad)+s2(jpicpt_zad) )/tvolt
	676	WRITE (numout,9447) t2(jpicpt_ldf)/tvolt, s2(jpicpt_ldf)/tvolt
	677	WRITE (numout,9448) t2(jpicpt_zdf)/tvolt, s2(jpicpt_zdf)/tvolt
	678	WRITE (numout,9449) t2(jpicpt_npc)/tvolt, s2(jpicpt_npc)/tvolt
[215]	679	ENDIF
	680
	681	9440 FORMAT(' tracer consistency at it= ',i6, &
	682	' : temperature',' salinity',/, &
	683	' ==================================')
[503]	684	9441 FORMAT(' 0 = horizontal+vertical advection + ',e20.13,' ',e20.13)
	685	9442 FORMAT(' 1st lev vertical advection ',e20.13,' ',e20.13)
	686	9443 FORMAT(' 0 = horizontal diffusion ',e20.13,' ',e20.13)
	687	9444 FORMAT(' 0 = vertical diffusion ',e20.13,' ',e20.13)
	688	9445 FORMAT(' 0 = static instability mixing ',e20.13,' ',e20.13)
	689	9446 FORMAT(' 0 = horizontal+vertical advection * t ',e20.13,' ',e20.13)
	690	9447 FORMAT(' 0 > horizontal diffusion * t ',e20.13,' ',e20.13)
	691	9448 FORMAT(' 0 > vertical diffusion * t ',e20.13,' ',e20.13)
	692	9449 FORMAT(' 0 > static instability mixing * t ',e20.13,' ',e20.13)
	693	!
[215]	694	ENDIF
[503]	695	!
[215]	696	END SUBROUTINE trd_twr
	697
	698	# else
	699	!!----------------------------------------------------------------------
	700	!! Default case : Empty module
	701	!!----------------------------------------------------------------------
[601]	702	INTERFACE trd_icp
	703	MODULE PROCEDURE trd_2d, trd_3d
	704	END INTERFACE
	705
[215]	706	CONTAINS
[521]	707	SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype, clpas ) ! Empty routine
[503]	708	REAL, DIMENSION(:,:) :: ptrd2dx, ptrd2dy
[601]	709	INTEGER , INTENT(in ) :: ktrd ! tracer trend index
	710	CHARACTER(len=3) , INTENT(in ) :: ctype ! momentum ('DYN') or tracers ('TRA') trends
	711	CHARACTER(len=3), INTENT(in), OPTIONAL :: clpas ! number of passage
[521]	712	WRITE(,) 'trd_2d: You should not have seen this print! error ?', &
	713	& ptrd2dx(1,1), ptrd2dy(1,1), ktrd, ctype, clpas
[215]	714	END SUBROUTINE trd_2d
[521]	715	SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd , ctype, clpas ) ! Empty routine
[503]	716	REAL, DIMENSION(:,:,:) :: ptrd3dx, ptrd3dy
[601]	717	INTEGER , INTENT(in ) :: ktrd ! tracer trend index
	718	CHARACTER(len=3) , INTENT(in ) :: ctype ! momentum ('DYN') or tracers ('TRA') trends
	719	CHARACTER(len=3), INTENT(in), OPTIONAL :: clpas ! number of passage
[521]	720	WRITE(,) 'trd_3d: You should not have seen this print! error ?', &
	721	& ptrd3dx(1,1,1), ptrd3dy(1,1,1), ktrd, ctype, clpas
[215]	722	END SUBROUTINE trd_3d
	723	SUBROUTINE trd_icp_init ! Empty routine
	724	END SUBROUTINE trd_icp_init
	725	SUBROUTINE trd_dwr( kt ) ! Empty routine
	726	WRITE(,) 'trd_dwr: You should not have seen this print! error ?', kt
	727	END SUBROUTINE trd_dwr
	728	SUBROUTINE trd_twr( kt ) ! Empty routine
	729	WRITE(,) 'trd_twr: You should not have seen this print! error ?', kt
	730	END SUBROUTINE trd_twr
	731	#endif
	732
	733	!!======================================================================
	734	END MODULE trdicp

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

Download in other formats: