Context Navigation

dynldf_iso.F90 @ 8568

Last change on this file since 8568 was 8568, checked in by gm, 7 years ago
#1911 (ENHANCE-09): PART I.2 - _NONE option + remove zts + see associated wiki page
Property svn:keywords set to `Id`
File size: 20.4 KB

Rev	Line
[3]	1	MODULE dynldf_iso
	2	!!======================================================================
	3	!! * MODULE dynldf_iso *
[5836]	4	!! Ocean dynamics: lateral viscosity trend (rotated laplacian operator)
[3]	5	!!======================================================================
[2715]	6	!! History : OPA ! 97-07 (G. Madec) Original code
	7	!! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module
	8	!! - ! 2004-08 (C. Talandier) New trends organization
	9	!! 2.0 ! 2005-11 (G. Madec) s-coordinate: horizontal diffusion
[5836]	10	!! 3.7 ! 2014-01 (F. Lemarie, G. Madec) restructuration/simplification of ahm specification,
	11	!! ! add velocity dependent coefficient and optional read in file
[2715]	12	!!----------------------------------------------------------------------
[5836]	13
[3]	14	!!----------------------------------------------------------------------
	15	!! dyn_ldf_iso : update the momentum trend with the horizontal part
	16	!! of the lateral diffusion using isopycnal or horizon-
	17	!! tal s-coordinate laplacian operator.
	18	!!----------------------------------------------------------------------
	19	USE oce ! ocean dynamics and tracers
	20	USE dom_oce ! ocean space and time domain
[5836]	21	USE ldfdyn ! lateral diffusion: eddy viscosity coef.
	22	USE ldftra ! lateral physics: eddy diffusivity
[3]	23	USE zdf_oce ! ocean vertical physics
	24	USE ldfslp ! iso-neutral slopes
[4990]	25	!
[3]	26	USE in_out_manager ! I/O manager
[2715]	27	USE lib_mpp ! MPP library
[5836]	28	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[258]	29	USE prtctl ! Print control
[3294]	30	USE timing ! Timing
[3]	31
	32	IMPLICIT NONE
	33	PRIVATE
	34
[2715]	35	PUBLIC dyn_ldf_iso ! called by step.F90
	36	PUBLIC dyn_ldf_iso_alloc ! called by nemogcm.F90
[3]	37
[8215]	38	REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: akzu, akzv !: vertical component of rotated lateral viscosity
	39
[2715]	40	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zfuw, zdiu, zdju, zdj1u ! 2D workspace (dyn_ldf_iso)
	41	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zfvw, zdiv, zdjv, zdj1v ! - -
	42
[3]	43	!! * Substitutions
	44	# include "vectopt_loop_substitute.h90"
	45	!!----------------------------------------------------------------------
[8568]	46	!! NEMO/OPA 4.0 , NEMO Consortium (2017)
[1156]	47	!! $Id$
[2715]	48	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	49	!!----------------------------------------------------------------------
	50	CONTAINS
	51
[2715]	52	INTEGER FUNCTION dyn_ldf_iso_alloc()
	53	!!----------------------------------------------------------------------
	54	!! * ROUTINE dyn_ldf_iso_alloc *
	55	!!----------------------------------------------------------------------
[8215]	56	ALLOCATE( akzu(jpi,jpj,jpk) , zfuw(jpi,jpk) , zdiu(jpi,jpk) , zdju(jpi,jpk) , zdj1u(jpi,jpk) , &
	57	& akzv(jpi,jpj,jpk) , zfvw(jpi,jpk) , zdiv(jpi,jpk) , zdjv(jpi,jpk) , zdj1v(jpi,jpk) , STAT=dyn_ldf_iso_alloc )
[2715]	58	!
	59	IF( dyn_ldf_iso_alloc /= 0 ) CALL ctl_warn('dyn_ldf_iso_alloc: array allocate failed.')
	60	END FUNCTION dyn_ldf_iso_alloc
	61
	62
[3]	63	SUBROUTINE dyn_ldf_iso( kt )
	64	!!----------------------------------------------------------------------
	65	!! * ROUTINE dyn_ldf_iso *
	66	!!
[455]	67	!! ** Purpose : Compute the before trend of the rotated laplacian
	68	!! operator of lateral momentum diffusion except the diagonal
	69	!! vertical term that will be computed in dynzdf module. Add it
	70	!! to the general trend of momentum equation.
[3]	71	!!
	72	!! ** Method :
[455]	73	!! The momentum lateral diffusive trend is provided by a 2nd
	74	!! order operator rotated along neutral or geopotential surfaces
	75	!! (in s-coordinates).
[3]	76	!! It is computed using before fields (forward in time) and isopyc-
[455]	77	!! nal or geopotential slopes computed in routine ldfslp.
[3]	78	!! Here, u and v components are considered as 2 independent scalar
	79	!! fields. Therefore, the property of splitting divergent and rota-
	80	!! tional part of the flow of the standard, z-coordinate laplacian
	81	!! momentum diffusion is lost.
	82	!! horizontal fluxes associated with the rotated lateral mixing:
	83	!! u-component:
[5836]	84	!! ziut = ( ahmt + rn_ahm_b ) e2t * e3t / e1t di[ ub ]
	85	!! - ahmt e2t * mi-1(uslp) dk[ mi(mk(ub)) ]
	86	!! zjuf = ( ahmf + rn_ahm_b ) e1f * e3f / e2f dj[ ub ]
	87	!! - ahmf e1f * mi(vslp) dk[ mj(mk(ub)) ]
[3]	88	!! v-component:
[5836]	89	!! zivf = ( ahmf + rn_ahm_b ) e2t * e3t / e1t di[ vb ]
	90	!! - ahmf e2t * mj(uslp) dk[ mi(mk(vb)) ]
	91	!! zjvt = ( ahmt + rn_ahm_b ) e1f * e3f / e2f dj[ ub ]
	92	!! - ahmt e1f * mj-1(vslp) dk[ mj(mk(vb)) ]
[3]	93	!! take the horizontal divergence of the fluxes:
	94	!! diffu = 1/(e1ue2ue3u) { di [ ziut ] + dj-1[ zjuf ] }
	95	!! diffv = 1/(e1ve2ve3v) { di-1[ zivf ] + dj [ zjvt ] }
	96	!! Add this trend to the general trend (ua,va):
	97	!! ua = ua + diffu
[455]	98	!! CAUTION: here the isopycnal part is with a coeff. of aht. This
	99	!! should be modified for applications others than orca_r2 (!!bug)
[3]	100	!!
	101	!! ** Action :
[8215]	102	!! -(ua,va) updated with the before geopotential harmonic mixing trend
	103	!! -(akzu,akzv) to accompt for the diagonal vertical component
	104	!! of the rotated operator in dynzdf module
[3]	105	!!----------------------------------------------------------------------
[2715]	106	INTEGER, INTENT( in ) :: kt ! ocean time-step index
	107	!
	108	INTEGER :: ji, jj, jk ! dummy loop indices
[8568]	109	REAL(wp) :: zabe1, zmskt, zmkt, zuav, zuwslpi, zuwslpj ! local scalars
	110	REAL(wp) :: zabe2, zmskf, zmkf, zvav, zvwslpi, zvwslpj ! - -
	111	REAL(wp) :: zcof0, zcof1, zcof2, zcof3, zcof4 ! - -
	112	REAL(wp), DIMENSION(jpi,jpj) :: ziut, zivf, zdku, zdk1u ! 2D workspace
	113	REAL(wp), DIMENSION(jpi,jpj) :: zjuf, zjvt, zdkv, zdk1v ! - -
[2715]	114	!!----------------------------------------------------------------------
[3294]	115	!
[8568]	116	IF( ln_timing ) CALL timing_start('dyn_ldf_iso')
[3294]	117	!
[3]	118	IF( kt == nit000 ) THEN
	119	IF(lwp) WRITE(numout,*)
	120	IF(lwp) WRITE(numout,*) 'dyn_ldf_iso : iso-neutral laplacian diffusive operator or '
	121	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ s-coordinate horizontal diffusive operator'
[2715]	122	! ! allocate dyn_ldf_bilap arrays
	123	IF( dyn_ldf_iso_alloc() /= 0 ) CALL ctl_stop('STOP', 'dyn_ldf_iso: failed to allocate arrays')
[3]	124	ENDIF
[216]	125
[5836]	126	!!gm bug is dyn_ldf_iso called before tra_ldf_iso .... <<<<<===== TO BE CHECKED
	127	! s-coordinate: Iso-level diffusion on momentum but not on tracer
[455]	128	IF( ln_dynldf_hor .AND. ln_traldf_iso ) THEN
[2715]	129	!
	130	DO jk = 1, jpk ! set the slopes of iso-level
[455]	131	DO jj = 2, jpjm1
[4488]	132	DO ji = 2, jpim1
[6140]	133	uslp (ji,jj,jk) = - ( gdept_b(ji+1,jj,jk) - gdept_b(ji ,jj ,jk) ) * r1_e1u(ji,jj) * umask(ji,jj,jk)
	134	vslp (ji,jj,jk) = - ( gdept_b(ji,jj+1,jk) - gdept_b(ji ,jj ,jk) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk)
	135	wslpi(ji,jj,jk) = - ( gdepw_b(ji+1,jj,jk) - gdepw_b(ji-1,jj,jk) ) * r1_e1t(ji,jj) * tmask(ji,jj,jk) * 0.5
	136	wslpj(ji,jj,jk) = - ( gdepw_b(ji,jj+1,jk) - gdepw_b(ji,jj-1,jk) ) * r1_e2t(ji,jj) * tmask(ji,jj,jk) * 0.5
[455]	137	END DO
	138	END DO
	139	END DO
	140	! Lateral boundary conditions on the slopes
	141	CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. )
	142	CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. )
[8215]	143	!
	144	ENDIF
[455]	145
[3]	146	! ! ===============
	147	DO jk = 1, jpkm1 ! Horizontal slab
	148	! ! ===============
	149
	150	! Vertical u- and v-shears at level jk and jk+1
	151	! ---------------------------------------------
	152	! surface boundary condition: zdku(jk=1)=zdku(jk=2)
	153	! zdkv(jk=1)=zdkv(jk=2)
	154
	155	zdk1u(:,:) = ( ub(:,:,jk) -ub(:,:,jk+1) ) * umask(:,:,jk+1)
	156	zdk1v(:,:) = ( vb(:,:,jk) -vb(:,:,jk+1) ) * vmask(:,:,jk+1)
	157
	158	IF( jk == 1 ) THEN
	159	zdku(:,:) = zdk1u(:,:)
	160	zdkv(:,:) = zdk1v(:,:)
	161	ELSE
	162	zdku(:,:) = ( ub(:,:,jk-1) - ub(:,:,jk) ) * umask(:,:,jk)
	163	zdkv(:,:) = ( vb(:,:,jk-1) - vb(:,:,jk) ) * vmask(:,:,jk)
	164	ENDIF
	165
	166	! -----f-----
	167	! Horizontal fluxes on U \|
	168	! --------------------=== t u t
	169	! \|
	170	! i-flux at t-point -----f-----
	171
[455]	172	IF( ln_zps ) THEN ! z-coordinate - partial steps : min(e3u)
	173	DO jj = 2, jpjm1
	174	DO ji = fs_2, jpi ! vector opt.
[6140]	175	zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * MIN( e3u_n(ji,jj,jk), e3u_n(ji-1,jj,jk) ) * r1_e1t(ji,jj)
[3]	176
[5836]	177	zmskt = 1._wp / MAX( umask(ji-1,jj,jk )+umask(ji,jj,jk+1) &
	178	& + umask(ji-1,jj,jk+1)+umask(ji,jj,jk ) , 1._wp )
[3]	179
[5836]	180	zcof1 = - rn_aht_0 * e2t(ji,jj) * zmskt * 0.5 * ( uslp(ji-1,jj,jk) + uslp(ji,jj,jk) )
[455]	181
[5836]	182	ziut(ji,jj) = ( zabe1 * ( ub(ji,jj,jk) - ub(ji-1,jj,jk) ) &
	183	& + zcof1 * ( zdku (ji,jj) + zdk1u(ji-1,jj) &
[455]	184	& +zdk1u(ji,jj) + zdku (ji-1,jj) ) ) * tmask(ji,jj,jk)
	185	END DO
	186	END DO
	187	ELSE ! other coordinate system (zco or sco) : e3t
	188	DO jj = 2, jpjm1
	189	DO ji = fs_2, jpi ! vector opt.
[6140]	190	zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * e3t_n(ji,jj,jk) * r1_e1t(ji,jj)
[3]	191
[5836]	192	zmskt = 1._wp / MAX( umask(ji-1,jj,jk ) + umask(ji,jj,jk+1) &
	193	& + umask(ji-1,jj,jk+1) + umask(ji,jj,jk ) , 1._wp )
[455]	194
[5836]	195	zcof1 = - rn_aht_0 * e2t(ji,jj) * zmskt * 0.5 * ( uslp(ji-1,jj,jk) + uslp(ji,jj,jk) )
[455]	196
	197	ziut(ji,jj) = ( zabe1 * ( ub(ji,jj,jk) - ub(ji-1,jj,jk) ) &
	198	& + zcof1 * ( zdku (ji,jj) + zdk1u(ji-1,jj) &
	199	& +zdk1u(ji,jj) + zdku (ji-1,jj) ) ) * tmask(ji,jj,jk)
	200	END DO
[3]	201	END DO
[455]	202	ENDIF
[3]	203
	204	! j-flux at f-point
	205	DO jj = 1, jpjm1
	206	DO ji = 1, fs_jpim1 ! vector opt.
[6140]	207	zabe2 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e1f(ji,jj) * e3f_n(ji,jj,jk) * r1_e2f(ji,jj)
[3]	208
[5836]	209	zmskf = 1._wp / MAX( umask(ji,jj+1,jk )+umask(ji,jj,jk+1) &
	210	& + umask(ji,jj+1,jk+1)+umask(ji,jj,jk ) , 1._wp )
[3]	211
[5836]	212	zcof2 = - rn_aht_0 * e1f(ji,jj) * zmskf * 0.5 * ( vslp(ji+1,jj,jk) + vslp(ji,jj,jk) )
[3]	213
[455]	214	zjuf(ji,jj) = ( zabe2 * ( ub(ji,jj+1,jk) - ub(ji,jj,jk) ) &
	215	& + zcof2 * ( zdku (ji,jj+1) + zdk1u(ji,jj) &
	216	& +zdk1u(ji,jj+1) + zdku (ji,jj) ) ) * fmask(ji,jj,jk)
[3]	217	END DO
	218	END DO
	219
	220	! \| t \|
	221	! Horizontal fluxes on V \| \|
	222	! --------------------=== f---v---f
	223	! \| \|
	224	! i-flux at f-point \| t \|
	225
	226	DO jj = 2, jpjm1
	227	DO ji = 1, fs_jpim1 ! vector opt.
[6140]	228	zabe1 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e2f(ji,jj) * e3f_n(ji,jj,jk) * r1_e1f(ji,jj)
[3]	229
[5836]	230	zmskf = 1._wp / MAX( vmask(ji+1,jj,jk )+vmask(ji,jj,jk+1) &
	231	& + vmask(ji+1,jj,jk+1)+vmask(ji,jj,jk ) , 1._wp )
[3]	232
[5836]	233	zcof1 = - rn_aht_0 * e2f(ji,jj) * zmskf * 0.5 * ( uslp(ji,jj+1,jk) + uslp(ji,jj,jk) )
[3]	234
[5836]	235	zivf(ji,jj) = ( zabe1 * ( vb(ji+1,jj,jk) - vb(ji,jj,jk) ) &
	236	& + zcof1 * ( zdkv (ji,jj) + zdk1v(ji+1,jj) &
	237	& + zdk1v(ji,jj) + zdkv (ji+1,jj) ) ) * fmask(ji,jj,jk)
[3]	238	END DO
	239	END DO
	240
	241	! j-flux at t-point
[455]	242	IF( ln_zps ) THEN ! z-coordinate - partial steps : min(e3u)
	243	DO jj = 2, jpj
	244	DO ji = 1, fs_jpim1 ! vector opt.
[6140]	245	zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * MIN( e3v_n(ji,jj,jk), e3v_n(ji,jj-1,jk) ) * r1_e2t(ji,jj)
[3]	246
[5836]	247	zmskt = 1._wp / MAX( vmask(ji,jj-1,jk )+vmask(ji,jj,jk+1) &
	248	& + vmask(ji,jj-1,jk+1)+vmask(ji,jj,jk ) , 1._wp )
[3]	249
[5836]	250	zcof2 = - rn_aht_0 * e1t(ji,jj) * zmskt * 0.5 * ( vslp(ji,jj-1,jk) + vslp(ji,jj,jk) )
[3]	251
[5836]	252	zjvt(ji,jj) = ( zabe2 * ( vb(ji,jj,jk) - vb(ji,jj-1,jk) ) &
	253	& + zcof2 * ( zdkv (ji,jj-1) + zdk1v(ji,jj) &
[455]	254	& +zdk1v(ji,jj-1) + zdkv (ji,jj) ) ) * tmask(ji,jj,jk)
	255	END DO
[3]	256	END DO
[455]	257	ELSE ! other coordinate system (zco or sco) : e3t
	258	DO jj = 2, jpj
	259	DO ji = 1, fs_jpim1 ! vector opt.
[6140]	260	zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * e3t_n(ji,jj,jk) * r1_e2t(ji,jj)
[3]	261
[455]	262	zmskt = 1./MAX( vmask(ji,jj-1,jk )+vmask(ji,jj,jk+1) &
	263	& + vmask(ji,jj-1,jk+1)+vmask(ji,jj,jk ), 1. )
[3]	264
[5836]	265	zcof2 = - rn_aht_0 * e1t(ji,jj) * zmskt * 0.5 * ( vslp(ji,jj-1,jk) + vslp(ji,jj,jk) )
[455]	266
	267	zjvt(ji,jj) = ( zabe2 * ( vb(ji,jj,jk) - vb(ji,jj-1,jk) ) &
	268	& + zcof2 * ( zdkv (ji,jj-1) + zdk1v(ji,jj) &
	269	& +zdk1v(ji,jj-1) + zdkv (ji,jj) ) ) * tmask(ji,jj,jk)
	270	END DO
	271	END DO
	272	ENDIF
	273
	274
[3]	275	! Second derivative (divergence) and add to the general trend
	276	! -----------------------------------------------------------
	277	DO jj = 2, jpjm1
[5836]	278	DO ji = 2, jpim1 !!gm Question vectop possible??? !!bug
[6140]	279	ua(ji,jj,jk) = ua(ji,jj,jk) + ( ziut(ji+1,jj) - ziut(ji,jj ) &
	280	& + zjuf(ji ,jj) - zjuf(ji,jj-1) ) * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk)
	281	va(ji,jj,jk) = va(ji,jj,jk) + ( zivf(ji,jj ) - zivf(ji-1,jj) &
	282	& + zjvt(ji,jj+1) - zjvt(ji,jj ) ) * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk)
[3]	283	END DO
	284	END DO
	285	! ! ===============
	286	END DO ! End of slab
	287	! ! ===============
[216]	288
[455]	289	! print sum trends (used for debugging)
	290	IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' ldfh - Ua: ', mask1=umask, &
	291	& tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[216]	292
	293
[455]	294	! ! ===============
	295	DO jj = 2, jpjm1 ! Vertical slab
	296	! ! ===============
	297
	298
	299	! I. vertical trends associated with the lateral mixing
	300	! =====================================================
	301	! (excluding the vertical flux proportional to dk[t]
	302
	303
	304	! I.1 horizontal momentum gradient
	305	! --------------------------------
	306
	307	DO jk = 1, jpk
	308	DO ji = 2, jpi
	309	! i-gradient of u at jj
	310	zdiu (ji,jk) = tmask(ji,jj ,jk) * ( ub(ji,jj ,jk) - ub(ji-1,jj ,jk) )
	311	! j-gradient of u and v at jj
	312	zdju (ji,jk) = fmask(ji,jj ,jk) * ( ub(ji,jj+1,jk) - ub(ji ,jj ,jk) )
	313	zdjv (ji,jk) = tmask(ji,jj ,jk) * ( vb(ji,jj ,jk) - vb(ji ,jj-1,jk) )
	314	! j-gradient of u and v at jj+1
	315	zdj1u(ji,jk) = fmask(ji,jj-1,jk) * ( ub(ji,jj ,jk) - ub(ji ,jj-1,jk) )
	316	zdj1v(ji,jk) = tmask(ji,jj+1,jk) * ( vb(ji,jj+1,jk) - vb(ji ,jj ,jk) )
	317	END DO
	318	END DO
	319	DO jk = 1, jpk
	320	DO ji = 1, jpim1
	321	! i-gradient of v at jj
	322	zdiv (ji,jk) = fmask(ji,jj ,jk) * ( vb(ji+1,jj,jk) - vb(ji ,jj ,jk) )
	323	END DO
	324	END DO
	325
	326
	327	! I.2 Vertical fluxes
	328	! -------------------
	329
	330	! Surface and bottom vertical fluxes set to zero
	331	DO ji = 1, jpi
	332	zfuw(ji, 1 ) = 0.e0
	333	zfvw(ji, 1 ) = 0.e0
	334	zfuw(ji,jpk) = 0.e0
	335	zfvw(ji,jpk) = 0.e0
	336	END DO
	337
	338	! interior (2=<jk=<jpk-1) on U field
	339	DO jk = 2, jpkm1
	340	DO ji = 2, jpim1
[8568]	341	zcof0 = 0.5_wp * rn_aht_0 * umask(ji,jj,jk)
[5836]	342	!
[8568]	343	zuwslpi = zcof0 * ( wslpi(ji+1,jj,jk) + wslpi(ji,jj,jk) )
	344	zuwslpj = zcof0 * ( wslpj(ji+1,jj,jk) + wslpj(ji,jj,jk) )
[5836]	345	!
[8568]	346	zmkt = 1./MAX( tmask(ji,jj,jk-1)+tmask(ji+1,jj,jk-1) &
	347	+ tmask(ji,jj,jk )+tmask(ji+1,jj,jk ) , 1. )
	348	zmkf = 1./MAX( fmask(ji,jj-1,jk-1) + fmask(ji,jj,jk-1) &
	349	+ fmask(ji,jj-1,jk ) + fmask(ji,jj,jk ) , 1. )
[455]	350
[8568]	351	zcof3 = - e2u(ji,jj) * zmkt * zuwslpi
	352	zcof4 = - e1u(ji,jj) * zmkf * zuwslpj
[455]	353	! vertical flux on u field
[8568]	354	zfuw(ji,jk) = zcof3 * ( zdiu (ji,jk-1) + zdiu (ji+1,jk-1) &
	355	& + zdiu (ji,jk ) + zdiu (ji+1,jk ) ) &
	356	& + zcof4 * ( zdj1u(ji,jk-1) + zdju (ji ,jk-1) &
	357	& + zdj1u(ji,jk ) + zdju (ji ,jk ) )
[8215]	358	! vertical mixing coefficient (akzu)
[8568]	359	! Note: zcof0 include rn_aht_0, so divided by rn_aht_0 to obtain slp^2 * rn_aht_0
[8215]	360	akzu(ji,jj,jk) = ( zuwslpi * zuwslpi + zuwslpj * zuwslpj ) / rn_aht_0
[455]	361	END DO
	362	END DO
	363
	364	! interior (2=<jk=<jpk-1) on V field
	365	DO jk = 2, jpkm1
	366	DO ji = 2, jpim1
[8568]	367	zcof0 = 0.5_wp * rn_aht_0 * vmask(ji,jj,jk)
	368	!
	369	zvwslpi = zcof0 * ( wslpi(ji,jj+1,jk) + wslpi(ji,jj,jk) )
	370	zvwslpj = zcof0 * ( wslpj(ji,jj+1,jk) + wslpj(ji,jj,jk) )
	371	!
	372	zmkf = 1./MAX( fmask(ji-1,jj,jk-1)+fmask(ji,jj,jk-1) &
	373	& + fmask(ji-1,jj,jk )+fmask(ji,jj,jk ) , 1. )
	374	zmkt = 1./MAX( tmask(ji,jj,jk-1)+tmask(ji,jj+1,jk-1) &
	375	& + tmask(ji,jj,jk )+tmask(ji,jj+1,jk ) , 1. )
[455]	376
[8568]	377	zcof3 = - e2v(ji,jj) * zmkf * zvwslpi
	378	zcof4 = - e1v(ji,jj) * zmkt * zvwslpj
[455]	379	! vertical flux on v field
[8568]	380	zfvw(ji,jk) = zcof3 * ( zdiv (ji,jk-1) + zdiv (ji-1,jk-1) &
	381	& + zdiv (ji,jk ) + zdiv (ji-1,jk ) ) &
	382	& + zcof4 * ( zdjv (ji,jk-1) + zdj1v(ji ,jk-1) &
	383	& + zdjv (ji,jk ) + zdj1v(ji ,jk ) )
[8215]	384	! vertical mixing coefficient (akzv)
[8568]	385	! Note: zcof0 include rn_aht_0, so divided by rn_aht_0 to obtain slp^2 * rn_aht_0
[8215]	386	akzv(ji,jj,jk) = ( zvwslpi * zvwslpi + zvwslpj * zvwslpj ) / rn_aht_0
[455]	387	END DO
	388	END DO
	389
	390
	391	! I.3 Divergence of vertical fluxes added to the general tracer trend
	392	! -------------------------------------------------------------------
	393	DO jk = 1, jpkm1
	394	DO ji = 2, jpim1
[6140]	395	ua(ji,jj,jk) = ua(ji,jj,jk) + ( zfuw(ji,jk) - zfuw(ji,jk+1) ) * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk)
	396	va(ji,jj,jk) = va(ji,jj,jk) + ( zfvw(ji,jk) - zfvw(ji,jk+1) ) * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk)
[455]	397	END DO
	398	END DO
	399	! ! ===============
	400	END DO ! End of slab
	401	! ! ===============
[2715]	402	!
[8568]	403	IF( ln_timing ) CALL timing_stop('dyn_ldf_iso')
[3294]	404	!
[3]	405	END SUBROUTINE dyn_ldf_iso
	406
	407	!!======================================================================
	408	END MODULE dynldf_iso

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

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

source: branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_iso.F90 @ 8568

Download in other formats: