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.

dynzdf_imp.F90 in NEMO/branches/2018/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/DYN – NEMO

Context Navigation

source: NEMO/branches/2018/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90 @ 10171

Last change on this file since 10171 was 10115, checked in by cbricaud, 6 years ago
phase 3.6 coarsening branch with nemo_3.6_rev9192
Property svn:keywords set to `Id`
File size: 19.0 KB

Rev	Line
[3]	1	MODULE dynzdf_imp
[2715]	2	!!======================================================================
[3]	3	!! * MODULE dynzdf_imp *
	4	!! Ocean dynamics: vertical component(s) of the momentum mixing trend
[2715]	5	!!======================================================================
[2528]	6	!! History : OPA ! 1990-10 (B. Blanke) Original code
	7	!! 8.0 ! 1997-05 (G. Madec) vertical component of isopycnal
[2715]	8	!! NEMO 0.5 ! 2002-08 (G. Madec) F90: Free form and module
[2528]	9	!! 3.3 ! 2010-04 (M. Leclair, G. Madec) Forcing averaged over 2 time steps
[3294]	10	!! 3.4 ! 2012-01 (H. Liu) Semi-implicit bottom friction
[503]	11	!!----------------------------------------------------------------------
[3]	12
	13	!!----------------------------------------------------------------------
[2715]	14	!! dyn_zdf_imp : update the momentum trend with the vertical diffusion using a implicit time-stepping
[3]	15	!!----------------------------------------------------------------------
	16	USE oce ! ocean dynamics and tracers
	17	USE dom_oce ! ocean space and time domain
[4292]	18	USE domvvl ! variable volume
[888]	19	USE sbc_oce ! surface boundary condition: ocean
	20	USE zdf_oce ! ocean vertical physics
[719]	21	USE phycst ! physical constants
[10115]	22	USE dynadv ! dynamics: vector invariant versus flux form
	23	USE dynspg_oce, ONLY: lk_dynspg_ts
	24	USE zdfbfr ! Bottom friction setup
	25	!
[3]	26	USE in_out_manager ! I/O manager
[10115]	27	USE iom ! I/O library
[2715]	28	USE lib_mpp ! MPP library
[3294]	29	USE wrk_nemo ! Memory Allocation
	30	USE timing ! Timing
[3]	31
	32	IMPLICIT NONE
	33	PRIVATE
	34
[2528]	35	PUBLIC dyn_zdf_imp ! called by step.F90
[3]	36
[4292]	37	REAL(wp) :: r_vvl ! variable volume indicator, =1 if lk_vvl=T, =0 otherwise
	38
[3]	39	!! * Substitutions
	40	# include "domzgr_substitute.h90"
	41	# include "vectopt_loop_substitute.h90"
	42	!!----------------------------------------------------------------------
[2528]	43	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[888]	44	!! $Id$
[2528]	45	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	46	!!----------------------------------------------------------------------
	47	CONTAINS
	48
[503]	49	SUBROUTINE dyn_zdf_imp( kt, p2dt )
[3]	50	!!----------------------------------------------------------------------
	51	!! * ROUTINE dyn_zdf_imp *
	52	!!
	53	!! ** Purpose : Compute the trend due to the vert. momentum diffusion
	54	!! and the surface forcing, and add it to the general trend of
	55	!! the momentum equations.
	56	!!
	57	!! ** Method : The vertical momentum mixing trend is given by :
	58	!! dz( avmu dz(u) ) = 1/e3u dk+1( avmu/e3uw dk(ua) )
	59	!! backward time stepping
[2528]	60	!! Surface boundary conditions: wind stress input (averaged over kt-1/2 & kt+1/2)
[3]	61	!! Bottom boundary conditions : bottom stress (cf zdfbfr.F)
	62	!! Add this trend to the general trend ua :
	63	!! ua = ua + dz( avmu dz(u) )
	64	!!
[2528]	65	!! ** Action : - Update (ua,va) arrays with the after vertical diffusive mixing trend.
[3]	66	!!---------------------------------------------------------------------
[3294]	67	INTEGER , INTENT(in) :: kt ! ocean time-step index
[2715]	68	REAL(wp), INTENT(in) :: p2dt ! vertical profile of tracer time-step
[2528]	69	!!
[2715]	70	INTEGER :: ji, jj, jk ! dummy loop indices
[3294]	71	INTEGER :: ikbu, ikbv ! local integers
[2715]	72	REAL(wp) :: z1_p2dt, zcoef, zzwi, zzws, zrhs ! local scalars
[10115]	73	REAL(wp) :: ze3ua, ze3va, zzz
	74	REAL(wp), POINTER, DIMENSION(:,:) :: z2d
[3294]	75	REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwd, zws
	76	!!----------------------------------------------------------------------
	77	!
	78	IF( nn_timing == 1 ) CALL timing_start('dyn_zdf_imp')
	79	!
	80	CALL wrk_alloc( jpi,jpj,jpk, zwi, zwd, zws )
	81	!
[3]	82	IF( kt == nit000 ) THEN
	83	IF(lwp) WRITE(numout,*)
	84	IF(lwp) WRITE(numout,*) 'dyn_zdf_imp : vertical momentum diffusion implicit operator'
	85	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
[4292]	86	!
	87	IF( lk_vvl ) THEN ; r_vvl = 1._wp ! Variable volume indicator
	88	ELSE ; r_vvl = 0._wp
	89	ENDIF
[3]	90	ENDIF
	91
	92	! 0. Local constant initialization
	93	! --------------------------------
[2528]	94	z1_p2dt = 1._wp / p2dt ! inverse of the timestep
[455]	95
[3294]	96	! 1. Apply semi-implicit bottom friction
	97	! --------------------------------------
	98	! Only needed for semi-implicit bottom friction setup. The explicit
	99	! bottom friction has been included in "u(v)a" which act as the R.H.S
	100	! column vector of the tri-diagonal matrix equation
	101	!
	102
	103	IF( ln_bfrimp ) THEN
[4292]	104	DO jj = 2, jpjm1
	105	DO ji = 2, jpim1
	106	ikbu = mbku(ji,jj) ! ocean bottom level at u- and v-points
	107	ikbv = mbkv(ji,jj) ! (deepest ocean u- and v-points)
	108	avmu(ji,jj,ikbu+1) = -bfrua(ji,jj) * fse3uw(ji,jj,ikbu+1)
	109	avmv(ji,jj,ikbv+1) = -bfrva(ji,jj) * fse3vw(ji,jj,ikbv+1)
	110	END DO
[3294]	111	END DO
[5602]	112	IF ( ln_isfcav ) THEN
	113	DO jj = 2, jpjm1
	114	DO ji = 2, jpim1
	115	ikbu = miku(ji,jj) ! ocean top level at u- and v-points
	116	ikbv = mikv(ji,jj) ! (first wet ocean u- and v-points)
	117	IF (ikbu .GE. 2) avmu(ji,jj,ikbu) = -tfrua(ji,jj) * fse3uw(ji,jj,ikbu)
	118	IF (ikbv .GE. 2) avmv(ji,jj,ikbv) = -tfrva(ji,jj) * fse3vw(ji,jj,ikbv)
	119	END DO
	120	END DO
	121	END IF
[3294]	122	ENDIF
	123
[4292]	124	#if defined key_dynspg_ts
	125	IF( ln_dynadv_vec .OR. .NOT. lk_vvl ) THEN ! applied on velocity
	126	DO jk = 1, jpkm1
	127	ua(:,:,jk) = ( ub(:,:,jk) + p2dt * ua(:,:,jk) ) * umask(:,:,jk)
	128	va(:,:,jk) = ( vb(:,:,jk) + p2dt * va(:,:,jk) ) * vmask(:,:,jk)
	129	END DO
	130	ELSE ! applied on thickness weighted velocity
	131	DO jk = 1, jpkm1
	132	ua(:,:,jk) = ( ub(:,:,jk) * fse3u_b(:,:,jk) &
	133	& + p2dt * ua(:,:,jk) * fse3u_n(:,:,jk) ) &
	134	& / fse3u_a(:,:,jk) * umask(:,:,jk)
	135	va(:,:,jk) = ( vb(:,:,jk) * fse3v_b(:,:,jk) &
	136	& + p2dt * va(:,:,jk) * fse3v_n(:,:,jk) ) &
	137	& / fse3v_a(:,:,jk) * vmask(:,:,jk)
	138	END DO
	139	ENDIF
	140
	141	IF ( ln_bfrimp .AND.lk_dynspg_ts ) THEN
	142	! remove barotropic velocities:
	143	DO jk = 1, jpkm1
	144	ua(:,:,jk) = (ua(:,:,jk) - ua_b(:,:)) * umask(:,:,jk)
	145	va(:,:,jk) = (va(:,:,jk) - va_b(:,:)) * vmask(:,:,jk)
[4990]	146	END DO
	147	! Add bottom/top stress due to barotropic component only:
[4292]	148	DO jj = 2, jpjm1
	149	DO ji = fs_2, fs_jpim1 ! vector opt.
	150	ikbu = mbku(ji,jj) ! ocean bottom level at u- and v-points
	151	ikbv = mbkv(ji,jj) ! (deepest ocean u- and v-points)
	152	ze3ua = ( 1._wp - r_vvl ) * fse3u_n(ji,jj,ikbu) + r_vvl * fse3u_a(ji,jj,ikbu)
	153	ze3va = ( 1._wp - r_vvl ) * fse3v_n(ji,jj,ikbv) + r_vvl * fse3v_a(ji,jj,ikbv)
	154	ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * bfrua(ji,jj) * ua_b(ji,jj) / ze3ua
	155	va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * bfrva(ji,jj) * va_b(ji,jj) / ze3va
	156	END DO
	157	END DO
[5602]	158	IF ( ln_isfcav ) THEN
	159	DO jj = 2, jpjm1
	160	DO ji = fs_2, fs_jpim1 ! vector opt.
	161	ikbu = miku(ji,jj) ! top ocean level at u- and v-points
	162	ikbv = mikv(ji,jj) ! (first wet ocean u- and v-points)
	163	ze3ua = ( 1._wp - r_vvl ) * fse3u_n(ji,jj,ikbu) + r_vvl * fse3u_a(ji,jj,ikbu)
	164	ze3va = ( 1._wp - r_vvl ) * fse3v_n(ji,jj,ikbv) + r_vvl * fse3v_a(ji,jj,ikbv)
	165	ua(ji,jj,ikbu) = ua(ji,jj,ikbu) + p2dt * tfrua(ji,jj) * ua_b(ji,jj) / ze3ua
	166	va(ji,jj,ikbv) = va(ji,jj,ikbv) + p2dt * tfrva(ji,jj) * va_b(ji,jj) / ze3va
	167	END DO
	168	END DO
	169	END IF
[4292]	170	ENDIF
	171	#endif
	172
[3294]	173	! 2. Vertical diffusion on u
[3]	174	! ---------------------------
	175	! Matrix and second member construction
[1662]	176	! bottom boundary condition: both zwi and zws must be masked as avmu can take
[3294]	177	! non zero value at the ocean bottom depending on the bottom friction used.
[2528]	178	!
	179	DO jk = 1, jpkm1 ! Matrix
[3]	180	DO jj = 2, jpjm1
	181	DO ji = fs_2, fs_jpim1 ! vector opt.
[4292]	182	ze3ua = ( 1._wp - r_vvl ) * fse3u_n(ji,jj,jk) + r_vvl * fse3u_a(ji,jj,jk) ! after scale factor at T-point
	183	zcoef = - p2dt / ze3ua
[5602]	184	zzwi = zcoef * avmu (ji,jj,jk ) / fse3uw(ji,jj,jk )
	185	zwi(ji,jj,jk) = zzwi * wumask(ji,jj,jk )
	186	zzws = zcoef * avmu (ji,jj,jk+1) / fse3uw(ji,jj,jk+1)
	187	zws(ji,jj,jk) = zzws * wumask(ji,jj,jk+1)
	188	zwd(ji,jj,jk) = 1._wp - zzwi - zzws
[3]	189	END DO
	190	END DO
	191	END DO
[4292]	192	DO jj = 2, jpjm1 ! Surface boundary conditions
[3]	193	DO ji = fs_2, fs_jpim1 ! vector opt.
[2528]	194	zwi(ji,jj,1) = 0._wp
	195	zwd(ji,jj,1) = 1._wp - zws(ji,jj,1)
[3]	196	END DO
	197	END DO
	198
	199	! Matrix inversion starting from the first level
	200	!-----------------------------------------------------------------------
	201	! solve m.x = y where m is a tri diagonal matrix ( jpk*jpk )
	202	!
	203	! ( zwd1 zws1 0 0 0 )( zwx1 ) ( zwy1 )
	204	! ( zwi2 zwd2 zws2 0 0 )( zwx2 ) ( zwy2 )
	205	! ( 0 zwi3 zwd3 zws3 0 )( zwx3 )=( zwy3 )
	206	! ( ... )( ... ) ( ... )
	207	! ( 0 0 0 zwik zwdk )( zwxk ) ( zwyk )
	208	!
	209	! m is decomposed in the product of an upper and a lower triangular matrix
	210	! The 3 diagonal terms are in 2d arrays: zwd, zws, zwi
	211	! The solution (the after velocity) is in ua
	212	!-----------------------------------------------------------------------
[2528]	213	!
[4990]	214	!== First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 (increasing k) ==
[5602]	215	DO jk = 2, jpkm1
	216	DO jj = 2, jpjm1
	217	DO ji = fs_2, fs_jpim1 ! vector opt.
[3]	218	zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1)
	219	END DO
	220	END DO
	221	END DO
[2528]	222	!
	223	DO jj = 2, jpjm1 !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==
[3]	224	DO ji = fs_2, fs_jpim1 ! vector opt.
[4292]	225	#if defined key_dynspg_ts
[5602]	226	ze3ua = ( 1._wp - r_vvl ) * fse3u_n(ji,jj,1) + r_vvl * fse3u_a(ji,jj,1)
	227	ua(ji,jj,1) = ua(ji,jj,1) + p2dt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) ) &
	228	& / ( ze3ua * rau0 ) * umask(ji,jj,1)
[4292]	229	#else
[5602]	230	ua(ji,jj,1) = ub(ji,jj,1) &
	231	& + p2dt (ua(ji,jj,1) + 0.5_wp ( utau_b(ji,jj) + utau(ji,jj) ) &
	232	& / ( fse3u(ji,jj,1) * rau0 ) * umask(ji,jj,1) )
[4292]	233	#endif
[5602]	234	END DO
	235	END DO
	236	DO jk = 2, jpkm1
	237	DO jj = 2, jpjm1
	238	DO ji = fs_2, fs_jpim1
[4292]	239	#if defined key_dynspg_ts
	240	zrhs = ua(ji,jj,jk) ! zrhs=right hand side
	241	#else
	242	zrhs = ub(ji,jj,jk) + p2dt * ua(ji,jj,jk)
	243	#endif
[3]	244	ua(ji,jj,jk) = zrhs - zwi(ji,jj,jk) / zwd(ji,jj,jk-1) * ua(ji,jj,jk-1)
	245	END DO
	246	END DO
	247	END DO
[2528]	248	!
	249	DO jj = 2, jpjm1 !== thrid recurrence : SOLk = ( Lk - Uk * Ek+1 ) / Dk ==
[3]	250	DO ji = fs_2, fs_jpim1 ! vector opt.
	251	ua(ji,jj,jpkm1) = ua(ji,jj,jpkm1) / zwd(ji,jj,jpkm1)
[5602]	252	END DO
	253	END DO
	254	DO jk = jpk-2, 1, -1
	255	DO jj = 2, jpjm1
	256	DO ji = fs_2, fs_jpim1
[2528]	257	ua(ji,jj,jk) = ( ua(ji,jj,jk) - zws(ji,jj,jk) * ua(ji,jj,jk+1) ) / zwd(ji,jj,jk)
[3]	258	END DO
	259	END DO
	260	END DO
	261
[3294]	262	! 3. Vertical diffusion on v
[3]	263	! ---------------------------
	264	! Matrix and second member construction
[1662]	265	! bottom boundary condition: both zwi and zws must be masked as avmv can take
[3294]	266	! non zero value at the ocean bottom depending on the bottom friction used
[2528]	267	!
	268	DO jk = 1, jpkm1 ! Matrix
[3]	269	DO jj = 2, jpjm1
	270	DO ji = fs_2, fs_jpim1 ! vector opt.
[4292]	271	ze3va = ( 1._wp - r_vvl ) * fse3v_n(ji,jj,jk) + r_vvl * fse3v_a(ji,jj,jk) ! after scale factor at T-point
	272	zcoef = - p2dt / ze3va
[2528]	273	zzwi = zcoef * avmv (ji,jj,jk ) / fse3vw(ji,jj,jk )
[5602]	274	zwi(ji,jj,jk) = zzwi * wvmask(ji,jj,jk)
[2528]	275	zzws = zcoef * avmv (ji,jj,jk+1) / fse3vw(ji,jj,jk+1)
[5602]	276	zws(ji,jj,jk) = zzws * wvmask(ji,jj,jk+1)
	277	zwd(ji,jj,jk) = 1._wp - zzwi - zzws
[3]	278	END DO
	279	END DO
	280	END DO
[4292]	281	DO jj = 2, jpjm1 ! Surface boundary conditions
[3]	282	DO ji = fs_2, fs_jpim1 ! vector opt.
[2528]	283	zwi(ji,jj,1) = 0._wp
	284	zwd(ji,jj,1) = 1._wp - zws(ji,jj,1)
[3]	285	END DO
	286	END DO
	287
	288	! Matrix inversion
	289	!-----------------------------------------------------------------------
	290	! solve m.x = y where m is a tri diagonal matrix ( jpk*jpk )
	291	!
	292	! ( zwd1 zws1 0 0 0 )( zwx1 ) ( zwy1 )
	293	! ( zwi2 zwd2 zws2 0 0 )( zwx2 ) ( zwy2 )
	294	! ( 0 zwi3 zwd3 zws3 0 )( zwx3 )=( zwy3 )
	295	! ( ... )( ... ) ( ... )
	296	! ( 0 0 0 zwik zwdk )( zwxk ) ( zwyk )
	297	!
[2528]	298	! m is decomposed in the product of an upper and lower triangular matrix
[3]	299	! The 3 diagonal terms are in 2d arrays: zwd, zws, zwi
	300	! The solution (after velocity) is in 2d array va
	301	!-----------------------------------------------------------------------
[2528]	302	!
[4990]	303	!== First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 (increasing k) ==
[5602]	304	DO jk = 2, jpkm1
	305	DO jj = 2, jpjm1
	306	DO ji = fs_2, fs_jpim1 ! vector opt.
[3]	307	zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1)
	308	END DO
	309	END DO
	310	END DO
[2528]	311	!
	312	DO jj = 2, jpjm1 !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==
[3]	313	DO ji = fs_2, fs_jpim1 ! vector opt.
[4292]	314	#if defined key_dynspg_ts
[5602]	315	ze3va = ( 1._wp - r_vvl ) * fse3v_n(ji,jj,1) + r_vvl * fse3v_a(ji,jj,1)
	316	va(ji,jj,1) = va(ji,jj,1) + p2dt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) ) &
[7256]	317	& / ( ze3va * rau0 ) * vmask(ji,jj,1)
[4292]	318	#else
[5602]	319	va(ji,jj,1) = vb(ji,jj,1) &
	320	& + p2dt (va(ji,jj,1) + 0.5_wp ( vtau_b(ji,jj) + vtau(ji,jj) ) &
[7256]	321	& / ( fse3v(ji,jj,1) * rau0 ) * vmask(ji,jj,1) )
[4292]	322	#endif
[5602]	323	END DO
	324	END DO
	325	DO jk = 2, jpkm1
	326	DO jj = 2, jpjm1
	327	DO ji = fs_2, fs_jpim1 ! vector opt.
[4292]	328	#if defined key_dynspg_ts
	329	zrhs = va(ji,jj,jk) ! zrhs=right hand side
	330	#else
	331	zrhs = vb(ji,jj,jk) + p2dt * va(ji,jj,jk)
	332	#endif
[3]	333	va(ji,jj,jk) = zrhs - zwi(ji,jj,jk) / zwd(ji,jj,jk-1) * va(ji,jj,jk-1)
	334	END DO
	335	END DO
	336	END DO
[2528]	337	!
[4292]	338	DO jj = 2, jpjm1 !== third recurrence : SOLk = ( Lk - Uk * SOLk+1 ) / Dk ==
[3]	339	DO ji = fs_2, fs_jpim1 ! vector opt.
	340	va(ji,jj,jpkm1) = va(ji,jj,jpkm1) / zwd(ji,jj,jpkm1)
[5602]	341	END DO
	342	END DO
	343	DO jk = jpk-2, 1, -1
	344	DO jj = 2, jpjm1
	345	DO ji = fs_2, fs_jpim1
[2528]	346	va(ji,jj,jk) = ( va(ji,jj,jk) - zws(ji,jj,jk) * va(ji,jj,jk+1) ) / zwd(ji,jj,jk)
[3]	347	END DO
	348	END DO
	349	END DO
	350
[10115]	351	IF( iom_use( 'dispkevfo' ) ) THEN ! ocean kinetic energy dissipation per unit area
	352	! ! due to v friction (v=vertical)
	353	! ! see NEMO_book appendix C, §C.8 (N.B. here averaged at t-points)
	354	! ! Note that formally, in a Leap-Frog environment, the shear**2 should be the product of
	355	! ! now by before shears, i.e. the source term of TKE (local positivity is not ensured).
	356	CALL wrk_alloc(jpi,jpj, z2d )
	357	z2d(:,:) = 0._wp
	358	DO jk = 1, jpkm1
	359	DO jj = 2, jpjm1
	360	DO ji = 2, jpim1
	361	z2d(ji,jj) = z2d(ji,jj) + ( &
	362	& avmu(ji ,jj,jk) * ( ua(ji ,jj,jk-1) - ua(ji ,jj,jk) )*2 / fse3uw(ji ,jj,jk) wumask(ji ,jj,jk) &
	363	& + avmu(ji-1,jj,jk) * ( ua(ji-1,jj,jk-1) - ua(ji-1,jj,jk) )*2 / fse3uw(ji-1,jj,jk) wumask(ji-1,jj,jk) &
	364	& + avmv(ji,jj ,jk) * ( va(ji,jj ,jk-1) - va(ji,jj ,jk) )*2 / fse3vw(ji,jj ,jk) wvmask(ji,jj ,jk) &
	365	& + avmv(ji,jj-1,jk) * ( va(ji,jj-1,jk-1) - va(ji,jj-1,jk) )*2 / fse3vw(ji,jj-1,jk) wvmask(ji,jj-1,jk) &
	366	& )
	367	END DO
	368	END DO
	369	END DO
	370	zzz= - 0.5_wp* rau0 ! caution sign minus here
	371	z2d(:,:) = zzz * z2d(:,:)
	372	CALL lbc_lnk( z2d,'T', 1. )
	373	CALL iom_put( 'dispkevfo', z2d )
	374	CALL wrk_dealloc(jpi,jpj, z2d )
	375	ENDIF
	376
[4292]	377	#if ! defined key_dynspg_ts
[10115]	378	!!gm this can be removed if tranxt is changed like in the trunk so that implicit outcome with
	379	!!gm the after velocity, not a trend
	380	! Normalization to obtain the general momentum trend ua
[3]	381	DO jk = 1, jpkm1
	382	DO jj = 2, jpjm1
	383	DO ji = fs_2, fs_jpim1 ! vector opt.
[10115]	384	ua(ji,jj,jk) = ( ua(ji,jj,jk) - ub(ji,jj,jk) ) * z1_p2dt
[2528]	385	va(ji,jj,jk) = ( va(ji,jj,jk) - vb(ji,jj,jk) ) * z1_p2dt
[3]	386	END DO
	387	END DO
	388	END DO
[4292]	389	#endif
[3294]	390
[4292]	391	! J. Chanut: Lines below are useless ?
[3294]	392	!! restore bottom layer avmu(v)
	393	IF( ln_bfrimp ) THEN
[4990]	394	DO jj = 2, jpjm1
	395	DO ji = 2, jpim1
	396	ikbu = mbku(ji,jj) ! ocean bottom level at u- and v-points
	397	ikbv = mbkv(ji,jj) ! (deepest ocean u- and v-points)
	398	avmu(ji,jj,ikbu+1) = 0.e0
	399	avmv(ji,jj,ikbv+1) = 0.e0
	400	END DO
	401	END DO
[5602]	402	IF (ln_isfcav) THEN
	403	DO jj = 2, jpjm1
	404	DO ji = 2, jpim1
	405	ikbu = miku(ji,jj) ! ocean top level at u- and v-points
	406	ikbv = mikv(ji,jj) ! (first wet ocean u- and v-points)
	407	IF (ikbu > 1) avmu(ji,jj,ikbu) = 0.e0
	408	IF (ikbv > 1) avmv(ji,jj,ikbv) = 0.e0
	409	END DO
	410	END DO
	411	END IF
[3294]	412	ENDIF
[2528]	413	!
[3294]	414	CALL wrk_dealloc( jpi,jpj,jpk, zwi, zwd, zws)
[2715]	415	!
[10115]	416	!
[3294]	417	IF( nn_timing == 1 ) CALL timing_stop('dyn_zdf_imp')
	418	!
[3]	419	END SUBROUTINE dyn_zdf_imp
	420
	421	!!==============================================================================
	422	END MODULE dynzdf_imp

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

Download in other formats: