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.

traldf_iso.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA – NEMO

Context Navigation

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf_iso.F90 @ 10806

Last change on this file since 10806 was 10806, checked in by davestorkey, 5 years ago
2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps branch: Latest updates. Make sure all time-dependent 3D variables are converted in previously modified modules.
Property svn:keywords set to `Id`
File size: 21.5 KB

Rev	Line
[3]	1	MODULE traldf_iso
[503]	2	!!======================================================================
[457]	3	!! * MODULE traldf_iso *
[2528]	4	!! Ocean tracers: horizontal component of the lateral tracer mixing trend
[503]	5	!!======================================================================
[5836]	6	!! History : OPA ! 1994-08 (G. Madec, M. Imbard)
	7	!! 8.0 ! 1997-05 (G. Madec) split into traldf and trazdf
	8	!! NEMO ! 2002-08 (G. Madec) Free form, F90
	9	!! 1.0 ! 2005-11 (G. Madec) merge traldf and trazdf :-)
	10	!! 3.3 ! 2010-09 (C. Ethe, G. Madec) Merge TRA-TRC
	11	!! 3.7 ! 2014-01 (G. Madec, S. Masson) restructuration/simplification of aht/aeiv specification
	12	!! - ! 2014-02 (F. Lemarie, G. Madec) triad operator (Griffies) + Method of Stabilizing Correction
[503]	13	!!----------------------------------------------------------------------
[5836]	14
[3]	15	!!----------------------------------------------------------------------
[6140]	16	!! tra_ldf_iso : update the tracer trend with the horizontal component of a iso-neutral laplacian operator
	17	!! and with the vertical part of the isopycnal or geopotential s-coord. operator
[3]	18	!!----------------------------------------------------------------------
[6140]	19	USE oce ! ocean dynamics and active tracers
	20	USE dom_oce ! ocean space and time domain
	21	USE trc_oce ! share passive tracers/Ocean variables
	22	USE zdf_oce ! ocean vertical physics
	23	USE ldftra ! lateral diffusion: tracer eddy coefficients
	24	USE ldfslp ! iso-neutral slopes
	25	USE diaptr ! poleward transport diagnostics
[7646]	26	USE diaar5 ! AR5 diagnostics
[5836]	27	!
[6140]	28	USE in_out_manager ! I/O manager
	29	USE iom ! I/O library
	30	USE phycst ! physical constants
	31	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[3]	32
	33	IMPLICIT NONE
	34	PRIVATE
	35
[503]	36	PUBLIC tra_ldf_iso ! routine called by step.F90
[3]	37
[7646]	38	LOGICAL :: l_ptr ! flag to compute poleward transport
	39	LOGICAL :: l_hst ! flag to compute heat transport
	40
[3]	41	!! * Substitutions
	42	# include "vectopt_loop_substitute.h90"
	43	!!----------------------------------------------------------------------
[9598]	44	!! NEMO/OCE 4.0 , NEMO Consortium (2018)
[2528]	45	!! $Id$
[10068]	46	!! Software governed by the CeCILL license (see ./LICENSE)
[247]	47	!!----------------------------------------------------------------------
[3]	48	CONTAINS
	49
[10806]	50	SUBROUTINE tra_ldf_iso( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu , pgv , &
[5836]	51	& pgui, pgvi, &
[10806]	52	& pt , pt_lev0, pt_rhs , kjpt, kpass )
[3]	53	!!----------------------------------------------------------------------
	54	!! * ROUTINE tra_ldf_iso *
[457]	55	!!
[3]	56	!! ** Purpose : Compute the before horizontal tracer (t & s) diffusive
[457]	57	!! trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and
	58	!! add it to the general trend of tracer equation.
[3]	59	!!
	60	!! ** Method : The horizontal component of the lateral diffusive trends
	61	!! is provided by a 2nd order operator rotated along neural or geopo-
	62	!! tential surfaces to which an eddy induced advection can be added
	63	!! It is computed using before fields (forward in time) and isopyc-
	64	!! nal or geopotential slopes computed in routine ldfslp.
	65	!!
[2528]	66	!! 1st part : masked horizontal derivative of T ( di[ t ] )
[5836]	67	!! ======== with partial cell update if ln_zps=T
	68	!! with top cell update if ln_isfcav
[457]	69	!!
	70	!! 2nd part : horizontal fluxes of the lateral mixing operator
	71	!! ========
[5836]	72	!! zftu = pahu e2u*e3u/e1u di[ tb ]
	73	!! - pahu e2u*uslp dk[ mi(mk(tb)) ]
	74	!! zftv = pahv e1v*e3v/e2v dj[ tb ]
	75	!! - pahv e2u*vslp dk[ mj(mk(tb)) ]
[3]	76	!! take the horizontal divergence of the fluxes:
[5836]	77	!! difft = 1/(e1e2t*e3t) { di-1[ zftu ] + dj-1[ zftv ] }
[3]	78	!! Add this trend to the general trend (ta,sa):
	79	!! ta = ta + difft
	80	!!
[457]	81	!! 3rd part: vertical trends of the lateral mixing operator
	82	!! ======== (excluding the vertical flux proportional to dk[t] )
	83	!! vertical fluxes associated with the rotated lateral mixing:
[5836]	84	!! zftw = - { mi(mk(pahu)) * e2t*wslpi di[ mi(mk(tb)) ]
	85	!! + mj(mk(pahv)) * e1t*wslpj dj[ mj(mk(tb)) ] }
[457]	86	!! take the horizontal divergence of the fluxes:
[5836]	87	!! difft = 1/(e1e2t*e3t) dk[ zftw ]
[457]	88	!! Add this trend to the general trend (ta,sa):
[10806]	89	!! pt_rhs = pt_rhs + difft
[3]	90	!!
[10806]	91	!! ** Action : Update pt_rhs arrays with the before rotated diffusion
[503]	92	!!----------------------------------------------------------------------
[2528]	93	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[5836]	94	INTEGER , INTENT(in ) :: kit000 ! first time step index
[10806]	95	INTEGER , INTENT(in ) :: ktlev ! time level index for e3t
	96	INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level thicknesses
[2528]	97	CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
	98	INTEGER , INTENT(in ) :: kjpt ! number of tracers
[5836]	99	INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage
	100	REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s]
	101	REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels
	102	REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels
[10806]	103	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2)
	104	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev0 ! tracer (only used in kpass=2)
	105	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend
[2715]	106	!
[2528]	107	INTEGER :: ji, jj, jk, jn ! dummy loop indices
[6140]	108	INTEGER :: ikt
[5836]	109	INTEGER :: ierr ! local integer
	110	REAL(wp) :: zmsku, zahu_w, zabe1, zcof1, zcoef3 ! local scalars
	111	REAL(wp) :: zmskv, zahv_w, zabe2, zcof2, zcoef4 ! - -
	112	REAL(wp) :: zcoef0, ze3w_2, zsign, z2dt, z1_2dt ! - -
[9019]	113	REAL(wp), DIMENSION(jpi,jpj) :: zdkt, zdk1t, z2d
	114	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdit, zdjt, zftu, zftv, ztfw
[3]	115	!!----------------------------------------------------------------------
[3294]	116	!
[9779]	117	IF( kpass == 1 .AND. kt == kit000 ) THEN
[3]	118	IF(lwp) WRITE(numout,*)
[2528]	119	IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype
[457]	120	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[5836]	121	!
[7753]	122	akz (:,:,:) = 0._wp
	123	ah_wslp2(:,:,:) = 0._wp
[3]	124	ENDIF
[7646]	125	!
	126	l_hst = .FALSE.
	127	l_ptr = .FALSE.
	128	IF( cdtype == 'TRA' .AND. ln_diaptr ) l_ptr = .TRUE.
	129	IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
	130	& iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.
	131	!
	132	! ! set time step size (Euler/Leapfrog)
[6140]	133	IF( neuler == 0 .AND. kt == nit000 ) THEN ; z2dt = rdt ! at nit000 (Euler)
	134	ELSE ; z2dt = 2.* rdt ! (Leapfrog)
[5836]	135	ENDIF
	136	z1_2dt = 1._wp / z2dt
	137	!
	138	IF( kpass == 1 ) THEN ; zsign = 1._wp ! bilaplacian operator require a minus sign (eddy diffusivity >0)
	139	ELSE ; zsign = -1._wp
	140	ENDIF
	141
	142	!!----------------------------------------------------------------------
	143	!! 0 - calculate ah_wslp2 and akz
	144	!!----------------------------------------------------------------------
	145	!
	146	IF( kpass == 1 ) THEN !== first pass only ==!
	147	!
	148	DO jk = 2, jpkm1
	149	DO jj = 2, jpjm1
	150	DO ji = fs_2, fs_jpim1 ! vector opt.
	151	!
[6140]	152	zmsku = wmask(ji,jj,jk) / MAX( umask(ji ,jj,jk-1) + umask(ji-1,jj,jk) &
[5836]	153	& + umask(ji-1,jj,jk-1) + umask(ji ,jj,jk) , 1._wp )
[6140]	154	zmskv = wmask(ji,jj,jk) / MAX( vmask(ji,jj ,jk-1) + vmask(ji,jj-1,jk) &
[5836]	155	& + vmask(ji,jj-1,jk-1) + vmask(ji,jj ,jk) , 1._wp )
	156	!
	157	zahu_w = ( pahu(ji ,jj,jk-1) + pahu(ji-1,jj,jk) &
	158	& + pahu(ji-1,jj,jk-1) + pahu(ji ,jj,jk) ) * zmsku
	159	zahv_w = ( pahv(ji,jj ,jk-1) + pahv(ji,jj-1,jk) &
	160	& + pahv(ji,jj-1,jk-1) + pahv(ji,jj ,jk) ) * zmskv
	161	!
	162	ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk) &
	163	& + zahv_w * wslpj(ji,jj,jk) * wslpj(ji,jj,jk)
	164	END DO
	165	END DO
	166	END DO
	167	!
	168	IF( ln_traldf_msc ) THEN ! stabilizing vertical diffusivity coefficient
	169	DO jk = 2, jpkm1
	170	DO jj = 2, jpjm1
	171	DO ji = fs_2, fs_jpim1
	172	akz(ji,jj,jk) = 0.25_wp * ( &
	173	& ( pahu(ji ,jj,jk) + pahu(ji ,jj,jk-1) ) / ( e1u(ji ,jj) * e1u(ji ,jj) ) &
	174	& + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) ) &
	175	& + ( pahv(ji,jj ,jk) + pahv(ji,jj ,jk-1) ) / ( e2v(ji,jj ) * e2v(ji,jj ) ) &
	176	& + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) ) )
	177	END DO
	178	END DO
	179	END DO
	180	!
	181	IF( ln_traldf_blp ) THEN ! bilaplacian operator
	182	DO jk = 2, jpkm1
	183	DO jj = 1, jpjm1
	184	DO ji = 1, fs_jpim1
	185	akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk) &
[10806]	186	& * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,kt2lev) * e3w(ji,jj,jk,kt2lev) ) )
[5836]	187	END DO
	188	END DO
	189	END DO
	190	ELSEIF( ln_traldf_lap ) THEN ! laplacian operator
	191	DO jk = 2, jpkm1
	192	DO jj = 1, jpjm1
	193	DO ji = 1, fs_jpim1
[10806]	194	ze3w_2 = e3w(ji,jj,jk,kt2lev) * e3w(ji,jj,jk,kt2lev)
[5836]	195	zcoef0 = z2dt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 )
	196	akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt
	197	END DO
	198	END DO
	199	END DO
	200	ENDIF
	201	!
	202	ELSE ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit
[7753]	203	akz(:,:,:) = ah_wslp2(:,:,:)
[5836]	204	ENDIF
	205	ENDIF
	206	!
[2528]	207	! ! ===========
	208	DO jn = 1, kjpt ! tracer loop
	209	! ! ===========
	210	!
	211	!!----------------------------------------------------------------------
	212	!! I - masked horizontal derivative
	213	!!----------------------------------------------------------------------
[5836]	214	!!gm : bug.... why (x,:,:)? (1,jpj,:) and (jpi,1,:) should be sufficient....
[7753]	215	zdit (1,:,:) = 0._wp ; zdit (jpi,:,:) = 0._wp
	216	zdjt (1,:,:) = 0._wp ; zdjt (jpi,:,:) = 0._wp
[2528]	217	!!end
[3]	218
[2528]	219	! Horizontal tracer gradient
	220	DO jk = 1, jpkm1
	221	DO jj = 1, jpjm1
	222	DO ji = 1, fs_jpim1 ! vector opt.
[10806]	223	zdit(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk)
	224	zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk)
[2528]	225	END DO
[457]	226	END DO
	227	END DO
[5836]	228	IF( ln_zps ) THEN ! botton and surface ocean correction of the horizontal gradient
	229	DO jj = 1, jpjm1 ! bottom correction (partial bottom cell)
[2528]	230	DO ji = 1, fs_jpim1 ! vector opt.
[5120]	231	zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)
	232	zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn)
	233	END DO
	234	END DO
[5836]	235	IF( ln_isfcav ) THEN ! first wet level beneath a cavity
	236	DO jj = 1, jpjm1
	237	DO ji = 1, fs_jpim1 ! vector opt.
	238	IF( miku(ji,jj) > 1 ) zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn)
	239	IF( mikv(ji,jj) > 1 ) zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn)
	240	END DO
	241	END DO
	242	ENDIF
[5120]	243	ENDIF
[6140]	244	!
[2528]	245	!!----------------------------------------------------------------------
	246	!! II - horizontal trend (full)
	247	!!----------------------------------------------------------------------
[5836]	248	!
	249	DO jk = 1, jpkm1 ! Horizontal slab
	250	!
	251	! !== Vertical tracer gradient
[10806]	252	zdk1t(:,:) = ( pt(:,:,jk,jn) - pt(:,:,jk+1,jn) ) * wmask(:,:,jk+1) ! level jk+1
[5836]	253	!
[7753]	254	IF( jk == 1 ) THEN ; zdkt(:,:) = zdk1t(:,:) ! surface: zdkt(jk=1)=zdkt(jk=2)
[10806]	255	ELSE ; zdkt(:,:) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * wmask(:,:,jk)
[5836]	256	ENDIF
	257	DO jj = 1 , jpjm1 !== Horizontal fluxes
[5120]	258	DO ji = 1, fs_jpim1 ! vector opt.
[10806]	259	zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,ktlev)
	260	zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,ktlev)
[2528]	261	!
[6140]	262	zmsku = 1. / MAX( wmask(ji+1,jj,jk ) + wmask(ji,jj,jk+1) &
	263	& + wmask(ji+1,jj,jk+1) + wmask(ji,jj,jk ), 1. )
[2528]	264	!
[6140]	265	zmskv = 1. / MAX( wmask(ji,jj+1,jk ) + wmask(ji,jj,jk+1) &
	266	& + wmask(ji,jj+1,jk+1) + wmask(ji,jj,jk ), 1. )
[2528]	267	!
[5836]	268	zcof1 = - pahu(ji,jj,jk) * e2u(ji,jj) * uslp(ji,jj,jk) * zmsku
	269	zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv
[2528]	270	!
	271	zftu(ji,jj,jk ) = ( zabe1 * zdit(ji,jj,jk) &
[5836]	272	& + zcof1 * ( zdkt (ji+1,jj) + zdk1t(ji,jj) &
	273	& + zdk1t(ji+1,jj) + zdkt (ji,jj) ) ) * umask(ji,jj,jk)
[2528]	274	zftv(ji,jj,jk) = ( zabe2 * zdjt(ji,jj,jk) &
[5836]	275	& + zcof2 * ( zdkt (ji,jj+1) + zdk1t(ji,jj) &
	276	& + zdk1t(ji,jj+1) + zdkt (ji,jj) ) ) * vmask(ji,jj,jk)
[2528]	277	END DO
[3]	278	END DO
[5836]	279	!
[10806]	280	DO jj = 2 , jpjm1 !== horizontal divergence and add to pt_rhs
[5120]	281	DO ji = fs_2, fs_jpim1 ! vector opt.
[10806]	282	pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk) &
[5836]	283	& + zftv(ji,jj,jk) - zftv(ji,jj-1,jk) ) &
[10806]	284	& * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev)
[2528]	285	END DO
[3]	286	END DO
[2528]	287	END DO ! End of slab
[3]	288
[2528]	289	!!----------------------------------------------------------------------
[5836]	290	!! III - vertical trend (full)
[2528]	291	!!----------------------------------------------------------------------
[6140]	292	!
[7753]	293	ztfw(1,:,:) = 0._wp ; ztfw(jpi,:,:) = 0._wp
[6140]	294	!
[2528]	295	! Vertical fluxes
	296	! ---------------
[6140]	297	! ! Surface and bottom vertical fluxes set to zero
[7753]	298	ztfw(:,:, 1 ) = 0._wp ; ztfw(:,:,jpk) = 0._wp
[2528]	299
[6140]	300	DO jk = 2, jpkm1 ! interior (2=<jk=<jpk-1)
[2528]	301	DO jj = 2, jpjm1
	302	DO ji = fs_2, fs_jpim1 ! vector opt.
	303	!
[6140]	304	zmsku = wmask(ji,jj,jk) / MAX( umask(ji ,jj,jk-1) + umask(ji-1,jj,jk) &
[5836]	305	& + umask(ji-1,jj,jk-1) + umask(ji ,jj,jk) , 1._wp )
[6140]	306	zmskv = wmask(ji,jj,jk) / MAX( vmask(ji,jj ,jk-1) + vmask(ji,jj-1,jk) &
[5836]	307	& + vmask(ji,jj-1,jk-1) + vmask(ji,jj ,jk) , 1._wp )
	308	!
	309	zahu_w = ( pahu(ji ,jj,jk-1) + pahu(ji-1,jj,jk) &
	310	& + pahu(ji-1,jj,jk-1) + pahu(ji ,jj,jk) ) * zmsku
	311	zahv_w = ( pahv(ji,jj ,jk-1) + pahv(ji,jj-1,jk) &
	312	& + pahv(ji,jj-1,jk-1) + pahv(ji,jj ,jk) ) * zmskv
	313	!
	314	zcoef3 = - zahu_w * e2t(ji,jj) * zmsku * wslpi (ji,jj,jk) !wslpi & j are already w-masked
	315	zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk)
[2528]	316	!
	317	ztfw(ji,jj,jk) = zcoef3 * ( zdit(ji ,jj ,jk-1) + zdit(ji-1,jj ,jk) &
	318	& + zdit(ji-1,jj ,jk-1) + zdit(ji ,jj ,jk) ) &
	319	& + zcoef4 * ( zdjt(ji ,jj ,jk-1) + zdjt(ji ,jj-1,jk) &
	320	& + zdjt(ji ,jj-1,jk-1) + zdjt(ji ,jj ,jk) )
	321	END DO
[457]	322	END DO
	323	END DO
[5836]	324	! !== add the vertical 33 flux ==!
	325	IF( ln_traldf_lap ) THEN ! laplacian case: eddy coef = ah_wslp2 - akz
	326	DO jk = 2, jpkm1
	327	DO jj = 1, jpjm1
	328	DO ji = fs_2, fs_jpim1
[10806]	329	ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * wmask(ji,jj,jk) &
[5836]	330	& * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) ) &
[10806]	331	& * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) )
[5836]	332	END DO
	333	END DO
	334	END DO
	335	!
	336	ELSE ! bilaplacian
	337	SELECT CASE( kpass )
	338	CASE( 1 ) ! 1st pass : eddy coef = ah_wslp2
	339	DO jk = 2, jpkm1
	340	DO jj = 1, jpjm1
	341	DO ji = fs_2, fs_jpim1
	342	ztfw(ji,jj,jk) = ztfw(ji,jj,jk) &
	343	& + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj) &
[10806]	344	& * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) / e3w(ji,jj,jk,kt2lev) * wmask(ji,jj,jk)
[5836]	345	END DO
	346	END DO
	347	END DO
[10806]	348	CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt and pt_lev0 gradients, resp.
[5836]	349	DO jk = 2, jpkm1
	350	DO jj = 1, jpjm1
	351	DO ji = fs_2, fs_jpim1
[10806]	352	ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * wmask(ji,jj,jk) &
	353	& * ( ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) ) &
	354	& + akz (ji,jj,jk) * ( pt_lev0(ji,jj,jk-1,jn) - pt_lev0(ji,jj,jk,jn) ) )
[5836]	355	END DO
	356	END DO
	357	END DO
	358	END SELECT
	359	ENDIF
	360	!
[10806]	361	DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to pt_rhs ==!
[2528]	362	DO jj = 2, jpjm1
	363	DO ji = fs_2, fs_jpim1 ! vector opt.
[10806]	364	pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1) ) &
	365	& * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev)
[2528]	366	END DO
[457]	367	END DO
	368	END DO
[2528]	369	!
[5836]	370	IF( ( kpass == 1 .AND. ln_traldf_lap ) .OR. & !== first pass only ( laplacian) ==!
	371	( kpass == 2 .AND. ln_traldf_blp ) ) THEN !== 2nd pass (bilaplacian) ==!
	372	!
	373	! ! "Poleward" diffusive heat or salt transports (T-S case only)
[7646]	374	! note sign is reversed to give down-gradient diffusive transports )
	375	IF( l_ptr ) CALL dia_ptr_hst( jn, 'ldf', -zftv(:,:,:) )
	376	! ! Diffusive heat transports
	377	IF( l_hst ) CALL dia_ar5_hst( jn, 'ldf', -zftu(:,:,:), -zftv(:,:,:) )
[5836]	378	!
	379	ENDIF !== end pass selection ==!
	380	!
	381	! ! ===============
	382	END DO ! end tracer loop
[503]	383	!
[3]	384	END SUBROUTINE tra_ldf_iso
	385
	386	!!==============================================================================
	387	END MODULE traldf_iso

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

Download in other formats: