Context Navigation

agrif_opa_interp.F90 @ 9005

Last change on this file since 9005 was 9005, checked in by timgraham, 6 years ago
Add TOP update code
Property svn:keywords set to `Id`
File size: 55.4 KB

Rev	Line
[636]	1	MODULE agrif_opa_interp
[1605]	2	!!======================================================================
	3	!! * MODULE agrif_opa_interp *
	4	!! AGRIF: interpolation package
	5	!!======================================================================
	6	!! History : 2.0 ! 2002-06 (XXX) Original cade
	7	!! - ! 2005-11 (XXX)
	8	!! 3.2 ! 2009-04 (R. Benshila)
[5656]	9	!! 3.6 ! 2014-09 (R. Benshila)
[1605]	10	!!----------------------------------------------------------------------
[7646]	11	#if defined key_agrif
[1605]	12	!!----------------------------------------------------------------------
	13	!! 'key_agrif' AGRIF zoom
	14	!!----------------------------------------------------------------------
	15	!! Agrif_tra :
	16	!! Agrif_dyn :
	17	!! interpu :
	18	!! interpv :
	19	!!----------------------------------------------------------------------
[636]	20	USE par_oce
	21	USE oce
	22	USE dom_oce
[6140]	23	USE zdf_oce
[782]	24	USE agrif_oce
[1605]	25	USE phycst
[8993]	26	USE dynspg_ts, ONLY: un_adv, vn_adv
[6140]	27	!
[1605]	28	USE in_out_manager
[2715]	29	USE agrif_opa_sponge
	30	USE lib_mpp
[4292]	31	USE wrk_nemo
[5656]	32
[636]	33	IMPLICIT NONE
	34	PRIVATE
[4292]	35
[4486]	36	PUBLIC Agrif_tra, Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_ssh_ts, Agrif_dta_ts
[5930]	37	PUBLIC interpun, interpvn
[5656]	38	PUBLIC interptsn, interpsshn
	39	PUBLIC interpunb, interpvnb, interpub2b, interpvb2b
	40	PUBLIC interpe3t, interpumsk, interpvmsk
[8993]	41	# if defined key_zdftke \|\| defined key_zdfgls
	42	PUBLIC Agrif_avm, interpavm
[5656]	43	# endif
[390]	44
[6140]	45	INTEGER :: bdy_tinterp = 0
	46
[1605]	47	# include "vectopt_loop_substitute.h90"
[1156]	48	!!----------------------------------------------------------------------
[6140]	49	!! NEMO/NST 3.7 , NEMO Consortium (2015)
[1156]	50	!! $Id$
[2528]	51	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1156]	52	!!----------------------------------------------------------------------
[5656]	53	CONTAINS
	54
[782]	55	SUBROUTINE Agrif_tra
[1605]	56	!!----------------------------------------------------------------------
[5656]	57	!! * ROUTINE Agrif_tra *
[1605]	58	!!----------------------------------------------------------------------
[636]	59	!
[1605]	60	IF( Agrif_Root() ) RETURN
[6140]	61	!
	62	Agrif_SpecialValue = 0._wp
[636]	63	Agrif_UseSpecialValue = .TRUE.
[6140]	64	!
[5656]	65	CALL Agrif_Bc_variable( tsn_id, procname=interptsn )
[6140]	66	!
[636]	67	Agrif_UseSpecialValue = .FALSE.
[1605]	68	!
[636]	69	END SUBROUTINE Agrif_tra
	70
[1605]	71
[636]	72	SUBROUTINE Agrif_dyn( kt )
[1605]	73	!!----------------------------------------------------------------------
	74	!! * ROUTINE Agrif_DYN *
	75	!!----------------------------------------------------------------------
	76	INTEGER, INTENT(in) :: kt
[6140]	77	!
	78	INTEGER :: ji, jj, jk ! dummy loop indices
	79	INTEGER :: j1, j2, i1, i2
	80	REAL(wp), POINTER, DIMENSION(:,:) :: zub, zvb
[1605]	81	!!----------------------------------------------------------------------
[6140]	82	!
[1605]	83	IF( Agrif_Root() ) RETURN
[6140]	84	!
	85	CALL wrk_alloc( jpi,jpj, zub, zvb )
	86	!
	87	Agrif_SpecialValue = 0._wp
[5656]	88	Agrif_UseSpecialValue = ln_spc_dyn
[6140]	89	!
	90	CALL Agrif_Bc_variable( un_interp_id, procname=interpun )
	91	CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn )
	92	!
[5656]	93	Agrif_UseSpecialValue = .FALSE.
[6140]	94	!
[5656]	95	! prevent smoothing in ghost cells
[6140]	96	i1 = 1 ; i2 = jpi
	97	j1 = 1 ; j2 = jpj
	98	IF( nbondj == -1 .OR. nbondj == 2 ) j1 = 3
	99	IF( nbondj == +1 .OR. nbondj == 2 ) j2 = nlcj-2
	100	IF( nbondi == -1 .OR. nbondi == 2 ) i1 = 3
	101	IF( nbondi == +1 .OR. nbondi == 2 ) i2 = nlci-2
[782]	102
[6140]	103	IF( nbondi == -1 .OR. nbondi == 2 ) THEN
	104	!
[5930]	105	! Smoothing
	106	! ---------
[6140]	107	IF( .NOT.ln_dynspg_ts ) THEN ! Store transport
	108	ua_b(2,:) = 0._wp
	109	DO jk = 1, jpkm1
	110	DO jj = 1, jpj
	111	ua_b(2,jj) = ua_b(2,jj) + e3u_a(2,jj,jk) * ua(2,jj,jk)
[5930]	112	END DO
[636]	113	END DO
[6140]	114	DO jj = 1, jpj
	115	ua_b(2,jj) = ua_b(2,jj) * r1_hu_a(2,jj)
[636]	116	END DO
[5930]	117	ENDIF
[6140]	118	!
[8993]	119	IF (.NOT.lk_agrif_clp) THEN
[8998]	120	DO jk=1,jpkm1 ! Smooth
[8993]	121	DO jj=j1,j2
	122	ua(2,jj,jk) = 0.25_wp(ua(1,jj,jk)+2._wpua(2,jj,jk)+ua(3,jj,jk))
	123	ua(2,jj,jk) = ua(2,jj,jk) * umask(2,jj,jk)
	124	END DO
[636]	125	END DO
[8998]	126	ENDIF
[6140]	127	!
	128	zub(2,:) = 0._wp ! Correct transport
	129	DO jk = 1, jpkm1
	130	DO jj = 1, jpj
	131	zub(2,jj) = zub(2,jj) + e3u_a(2,jj,jk) * ua(2,jj,jk)
[636]	132	END DO
	133	END DO
	134	DO jj=1,jpj
[6140]	135	zub(2,jj) = zub(2,jj) * r1_hu_a(2,jj)
[636]	136	END DO
[390]	137
[636]	138	DO jk=1,jpkm1
[5930]	139	DO jj=1,jpj
	140	ua(2,jj,jk) = (ua(2,jj,jk)+ua_b(2,jj)-zub(2,jj))*umask(2,jj,jk)
[636]	141	END DO
	142	END DO
[390]	143
[4486]	144	! Set tangential velocities to time splitting estimate
[5930]	145	!-----------------------------------------------------
[6140]	146	IF( ln_dynspg_ts ) THEN
	147	zvb(2,:) = 0._wp
	148	DO jk = 1, jpkm1
	149	DO jj = 1, jpj
	150	zvb(2,jj) = zvb(2,jj) + e3v_a(2,jj,jk) * va(2,jj,jk)
[5930]	151	END DO
	152	END DO
[6140]	153	DO jj = 1, jpj
	154	zvb(2,jj) = zvb(2,jj) * r1_hv_a(2,jj)
[4486]	155	END DO
[6140]	156	DO jk = 1, jpkm1
	157	DO jj = 1, jpj
	158	va(2,jj,jk) = (va(2,jj,jk)+va_b(2,jj)-zvb(2,jj)) * vmask(2,jj,jk)
[5930]	159	END DO
	160	END DO
	161	ENDIF
[6140]	162	!
[5930]	163	! Mask domain edges:
	164	!-------------------
[6140]	165	DO jk = 1, jpkm1
	166	DO jj = 1, jpj
[5930]	167	ua(1,jj,jk) = 0._wp
	168	va(1,jj,jk) = 0._wp
[4486]	169	END DO
[5930]	170	END DO
[6140]	171	!
[636]	172	ENDIF
[390]	173
[6140]	174	IF( nbondi == 1 .OR. nbondi == 2 ) THEN
[5930]	175
	176	! Smoothing
	177	! ---------
[6140]	178	IF( .NOT.ln_dynspg_ts ) THEN ! Store transport
	179	ua_b(nlci-2,:) = 0._wp
[5930]	180	DO jk=1,jpkm1
	181	DO jj=1,jpj
[6140]	182	ua_b(nlci-2,jj) = ua_b(nlci-2,jj) + e3u_a(nlci-2,jj,jk) * ua(nlci-2,jj,jk)
[5930]	183	END DO
[636]	184	END DO
	185	DO jj=1,jpj
[6140]	186	ua_b(nlci-2,jj) = ua_b(nlci-2,jj) * r1_hu_a(nlci-2,jj)
[5930]	187	END DO
	188	ENDIF
	189
[8993]	190	IF (.NOT.lk_agrif_clp) THEN
[8998]	191	DO jk = 1, jpkm1 ! Smooth
[8993]	192	DO jj = j1, j2
	193	ua(nlci-2,jj,jk) = 0.25_wp * umask(nlci-2,jj,jk) &
	194	& * ( ua(nlci-3,jj,jk) + 2._wp*ua(nlci-2,jj,jk) + ua(nlci-1,jj,jk) )
	195	END DO
[636]	196	END DO
[8993]	197	ENDIF
[6140]	198	zub(nlci-2,:) = 0._wp ! Correct transport
	199	DO jk = 1, jpkm1
	200	DO jj = 1, jpj
	201	zub(nlci-2,jj) = zub(nlci-2,jj) + e3u_a(nlci-2,jj,jk) * ua(nlci-2,jj,jk)
[636]	202	END DO
	203	END DO
[6140]	204	DO jj = 1, jpj
	205	zub(nlci-2,jj) = zub(nlci-2,jj) * r1_hu_a(nlci-2,jj)
[636]	206	END DO
[5930]	207
[6140]	208	DO jk = 1, jpkm1
	209	DO jj = 1, jpj
	210	ua(nlci-2,jj,jk) = ( ua(nlci-2,jj,jk) + ua_b(nlci-2,jj) - zub(nlci-2,jj) ) * umask(nlci-2,jj,jk)
[636]	211	END DO
	212	END DO
[6140]	213	!
[4486]	214	! Set tangential velocities to time splitting estimate
[5930]	215	!-----------------------------------------------------
[6140]	216	IF( ln_dynspg_ts ) THEN
	217	zvb(nlci-1,:) = 0._wp
	218	DO jk = 1, jpkm1
	219	DO jj = 1, jpj
	220	zvb(nlci-1,jj) = zvb(nlci-1,jj) + e3v_a(nlci-1,jj,jk) * va(nlci-1,jj,jk)
[5930]	221	END DO
	222	END DO
[4486]	223	DO jj=1,jpj
[6140]	224	zvb(nlci-1,jj) = zvb(nlci-1,jj) * r1_hv_a(nlci-1,jj)
[4486]	225	END DO
[6140]	226	DO jk = 1, jpkm1
	227	DO jj = 1, jpj
	228	va(nlci-1,jj,jk) = ( va(nlci-1,jj,jk) + va_b(nlci-1,jj) - zvb(nlci-1,jj) ) * vmask(nlci-1,jj,jk)
[5930]	229	END DO
	230	END DO
	231	ENDIF
[6140]	232	!
[5930]	233	! Mask domain edges:
	234	!-------------------
[6140]	235	DO jk = 1, jpkm1
	236	DO jj = 1, jpj
[5930]	237	ua(nlci-1,jj,jk) = 0._wp
	238	va(nlci ,jj,jk) = 0._wp
[4486]	239	END DO
[5930]	240	END DO
[6140]	241	!
[636]	242	ENDIF
[390]	243
[6140]	244	IF( nbondj == -1 .OR. nbondj == 2 ) THEN
[390]	245
[5930]	246	! Smoothing
	247	! ---------
[6140]	248	IF( .NOT.ln_dynspg_ts ) THEN ! Store transport
	249	va_b(:,2) = 0._wp
	250	DO jk = 1, jpkm1
	251	DO ji = 1, jpi
	252	va_b(ji,2) = va_b(ji,2) + e3v_a(ji,2,jk) * va(ji,2,jk)
[5930]	253	END DO
[636]	254	END DO
	255	DO ji=1,jpi
[6140]	256	va_b(ji,2) = va_b(ji,2) * r1_hv_a(ji,2)
[636]	257	END DO
[5930]	258	ENDIF
[6140]	259	!
[8993]	260	IF (.NOT.lk_agrif_clp) THEN
	261	DO jk = 1, jpkm1 ! Smooth
	262	DO ji = i1, i2
	263	va(ji,2,jk) = 0.25_wp * vmask(ji,2,jk) &
	264	& * ( va(ji,1,jk) + 2._wp*va(ji,2,jk) + va(ji,3,jk) )
	265	END DO
[636]	266	END DO
[8993]	267	ENDIF
[6140]	268	!
	269	zvb(:,2) = 0._wp ! Correct transport
[636]	270	DO jk=1,jpkm1
	271	DO ji=1,jpi
[6140]	272	zvb(ji,2) = zvb(ji,2) + e3v_a(ji,2,jk) * va(ji,2,jk) * vmask(ji,2,jk)
[636]	273	END DO
	274	END DO
[6140]	275	DO ji = 1, jpi
	276	zvb(ji,2) = zvb(ji,2) * r1_hv_a(ji,2)
[636]	277	END DO
[6140]	278	DO jk = 1, jpkm1
	279	DO ji = 1, jpi
	280	va(ji,2,jk) = ( va(ji,2,jk) + va_b(ji,2) - zvb(ji,2) ) * vmask(ji,2,jk)
[636]	281	END DO
	282	END DO
[390]	283
[4486]	284	! Set tangential velocities to time splitting estimate
[5930]	285	!-----------------------------------------------------
[6140]	286	IF( ln_dynspg_ts ) THEN
	287	zub(:,2) = 0._wp
	288	DO jk = 1, jpkm1
	289	DO ji = 1, jpi
	290	zub(ji,2) = zub(ji,2) + e3u_a(ji,2,jk) * ua(ji,2,jk) * umask(ji,2,jk)
[5930]	291	END DO
	292	END DO
[6140]	293	DO ji = 1, jpi
	294	zub(ji,2) = zub(ji,2) * r1_hu_a(ji,2)
[4486]	295	END DO
	296
[6140]	297	DO jk = 1, jpkm1
	298	DO ji = 1, jpi
	299	ua(ji,2,jk) = ( ua(ji,2,jk) + ua_b(ji,2) - zub(ji,2) ) * umask(ji,2,jk)
[5930]	300	END DO
	301	END DO
	302	ENDIF
[4486]	303
[5930]	304	! Mask domain edges:
	305	!-------------------
[6140]	306	DO jk = 1, jpkm1
	307	DO ji = 1, jpi
[5930]	308	ua(ji,1,jk) = 0._wp
	309	va(ji,1,jk) = 0._wp
[4486]	310	END DO
[5930]	311	END DO
	312
[636]	313	ENDIF
[390]	314
[6140]	315	IF( nbondj == 1 .OR. nbondj == 2 ) THEN
	316	!
[5930]	317	! Smoothing
	318	! ---------
[6140]	319	IF( .NOT.ln_dynspg_ts ) THEN ! Store transport
	320	va_b(:,nlcj-2) = 0._wp
	321	DO jk = 1, jpkm1
	322	DO ji = 1, jpi
	323	va_b(ji,nlcj-2) = va_b(ji,nlcj-2) + e3v_a(ji,nlcj-2,jk) * va(ji,nlcj-2,jk)
[5930]	324	END DO
[636]	325	END DO
[6140]	326	DO ji = 1, jpi
	327	va_b(ji,nlcj-2) = va_b(ji,nlcj-2) * r1_hv_a(ji,nlcj-2)
[636]	328	END DO
[5930]	329	ENDIF
[6140]	330	!
[8993]	331	IF (.NOT.lk_agrif_clp) THEN
[8998]	332	DO jk = 1, jpkm1 ! Smooth
[8993]	333	DO ji = i1, i2
	334	va(ji,nlcj-2,jk) = 0.25_wp * vmask(ji,nlcj-2,jk) &
	335	& * ( va(ji,nlcj-3,jk) + 2._wp * va(ji,nlcj-2,jk) + va(ji,nlcj-1,jk) )
	336	END DO
[636]	337	END DO
[8998]	338	END IF
[6140]	339	!
	340	zvb(:,nlcj-2) = 0._wp ! Correct transport
	341	DO jk = 1, jpkm1
	342	DO ji = 1, jpi
	343	zvb(ji,nlcj-2) = zvb(ji,nlcj-2) + e3v_a(ji,nlcj-2,jk) * va(ji,nlcj-2,jk) * vmask(ji,nlcj-2,jk)
[636]	344	END DO
	345	END DO
[6140]	346	DO ji = 1, jpi
	347	zvb(ji,nlcj-2) = zvb(ji,nlcj-2) * r1_hv_a(ji,nlcj-2)
[636]	348	END DO
[6140]	349	DO jk = 1, jpkm1
	350	DO ji = 1, jpi
	351	va(ji,nlcj-2,jk) = ( va(ji,nlcj-2,jk) + va_b(ji,nlcj-2) - zvb(ji,nlcj-2) ) * vmask(ji,nlcj-2,jk)
[636]	352	END DO
	353	END DO
[6140]	354	!
[4486]	355	! Set tangential velocities to time splitting estimate
[5930]	356	!-----------------------------------------------------
[6140]	357	IF( ln_dynspg_ts ) THEN
	358	zub(:,nlcj-1) = 0._wp
	359	DO jk = 1, jpkm1
	360	DO ji = 1, jpi
	361	zub(ji,nlcj-1) = zub(ji,nlcj-1) + e3u_a(ji,nlcj-1,jk) * ua(ji,nlcj-1,jk) * umask(ji,nlcj-1,jk)
[5930]	362	END DO
	363	END DO
[6140]	364	DO ji = 1, jpi
	365	zub(ji,nlcj-1) = zub(ji,nlcj-1) * r1_hu_a(ji,nlcj-1)
[4486]	366	END DO
[6140]	367	!
	368	DO jk = 1, jpkm1
	369	DO ji = 1, jpi
	370	ua(ji,nlcj-1,jk) = ( ua(ji,nlcj-1,jk) + ua_b(ji,nlcj-1) - zub(ji,nlcj-1) ) * umask(ji,nlcj-1,jk)
[5930]	371	END DO
	372	END DO
	373	ENDIF
[6140]	374	!
[5930]	375	! Mask domain edges:
	376	!-------------------
[6140]	377	DO jk = 1, jpkm1
	378	DO ji = 1, jpi
[5930]	379	ua(ji,nlcj ,jk) = 0._wp
	380	va(ji,nlcj-1,jk) = 0._wp
[4486]	381	END DO
[5930]	382	END DO
[6140]	383	!
[636]	384	ENDIF
[2715]	385	!
[6140]	386	CALL wrk_dealloc( jpi,jpj, zub, zvb )
[2715]	387	!
[636]	388	END SUBROUTINE Agrif_dyn
[390]	389
[6140]	390
[4486]	391	SUBROUTINE Agrif_dyn_ts( jn )
[4292]	392	!!----------------------------------------------------------------------
	393	!! * ROUTINE Agrif_dyn_ts *
	394	!!----------------------------------------------------------------------
	395	!!
[4486]	396	INTEGER, INTENT(in) :: jn
[4292]	397	!!
	398	INTEGER :: ji, jj
[4486]	399	!!----------------------------------------------------------------------
[6140]	400	!
[4486]	401	IF( Agrif_Root() ) RETURN
[6140]	402	!
[4486]	403	IF((nbondi == -1).OR.(nbondi == 2)) THEN
	404	DO jj=1,jpj
	405	va_e(2,jj) = vbdy_w(jj) * hvr_e(2,jj)
[5656]	406	! Specified fluxes:
[4486]	407	ua_e(2,jj) = ubdy_w(jj) * hur_e(2,jj)
[5656]	408	! Characteristics method:
	409	!alt ua_e(2,jj) = 0.5_wp * ( ubdy_w(jj) * hur_e(2,jj) + ua_e(3,jj) &
	410	!alt & - sqrt(grav * hur_e(2,jj)) * (sshn_e(3,jj) - hbdy_w(jj)) )
[4486]	411	END DO
	412	ENDIF
[6140]	413	!
[4486]	414	IF((nbondi == 1).OR.(nbondi == 2)) THEN
	415	DO jj=1,jpj
	416	va_e(nlci-1,jj) = vbdy_e(jj) * hvr_e(nlci-1,jj)
[5656]	417	! Specified fluxes:
[4486]	418	ua_e(nlci-2,jj) = ubdy_e(jj) * hur_e(nlci-2,jj)
[5656]	419	! Characteristics method:
	420	!alt ua_e(nlci-2,jj) = 0.5_wp * ( ubdy_e(jj) * hur_e(nlci-2,jj) + ua_e(nlci-3,jj) &
	421	!alt & + sqrt(grav * hur_e(nlci-2,jj)) * (sshn_e(nlci-2,jj) - hbdy_e(jj)) )
[4486]	422	END DO
	423	ENDIF
[6140]	424	!
[4486]	425	IF((nbondj == -1).OR.(nbondj == 2)) THEN
	426	DO ji=1,jpi
	427	ua_e(ji,2) = ubdy_s(ji) * hur_e(ji,2)
[5656]	428	! Specified fluxes:
[4486]	429	va_e(ji,2) = vbdy_s(ji) * hvr_e(ji,2)
[5656]	430	! Characteristics method:
	431	!alt va_e(ji,2) = 0.5_wp * ( vbdy_s(ji) * hvr_e(ji,2) + va_e(ji,3) &
	432	!alt & - sqrt(grav * hvr_e(ji,2)) * (sshn_e(ji,3) - hbdy_s(ji)) )
[4486]	433	END DO
	434	ENDIF
[6140]	435	!
[4486]	436	IF((nbondj == 1).OR.(nbondj == 2)) THEN
	437	DO ji=1,jpi
	438	ua_e(ji,nlcj-1) = ubdy_n(ji) * hur_e(ji,nlcj-1)
[5656]	439	! Specified fluxes:
[4486]	440	va_e(ji,nlcj-2) = vbdy_n(ji) * hvr_e(ji,nlcj-2)
[5656]	441	! Characteristics method:
	442	!alt va_e(ji,nlcj-2) = 0.5_wp * ( vbdy_n(ji) * hvr_e(ji,nlcj-2) + va_e(ji,nlcj-3) &
	443	!alt & + sqrt(grav * hvr_e(ji,nlcj-2)) * (sshn_e(ji,nlcj-2) - hbdy_n(ji)) )
[4486]	444	END DO
	445	ENDIF
	446	!
	447	END SUBROUTINE Agrif_dyn_ts
	448
[6140]	449
[4486]	450	SUBROUTINE Agrif_dta_ts( kt )
	451	!!----------------------------------------------------------------------
	452	!! * ROUTINE Agrif_dta_ts *
	453	!!----------------------------------------------------------------------
	454	!!
	455	INTEGER, INTENT(in) :: kt
	456	!!
	457	INTEGER :: ji, jj
	458	LOGICAL :: ll_int_cons
[4292]	459	!!----------------------------------------------------------------------
[6140]	460	!
[4292]	461	IF( Agrif_Root() ) RETURN
[6140]	462	!
	463	ll_int_cons = ln_bt_fw ! Assume conservative temporal integration in the forward case only
	464	!
[8993]	465	! Enforce volume conservation if no time refinement:
	466	IF ( Agrif_rhot()==1 ) ll_int_cons=.TRUE.
[6140]	467	!
[4486]	468	! Interpolate barotropic fluxes
[8993]	469	Agrif_SpecialValue=0._wp
[4486]	470	Agrif_UseSpecialValue = ln_spc_dyn
[6140]	471	!
	472	IF( ll_int_cons ) THEN ! Conservative interpolation
[5656]	473	! orders matters here !!!!!!
[6140]	474	CALL Agrif_Bc_variable( ub2b_interp_id, calledweight=1._wp, procname=interpub2b ) ! Time integrated
	475	CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b )
[5656]	476	bdy_tinterp = 1
[6140]	477	CALL Agrif_Bc_variable( unb_id , calledweight=1._wp, procname=interpunb ) ! After
	478	CALL Agrif_Bc_variable( vnb_id , calledweight=1._wp, procname=interpvnb )
[5656]	479	bdy_tinterp = 2
[6140]	480	CALL Agrif_Bc_variable( unb_id , calledweight=0._wp, procname=interpunb ) ! Before
	481	CALL Agrif_Bc_variable( vnb_id , calledweight=0._wp, procname=interpvnb )
[4486]	482	ELSE ! Linear interpolation
[5656]	483	bdy_tinterp = 0
[6140]	484	ubdy_w(:) = 0._wp ; vbdy_w(:) = 0._wp
	485	ubdy_e(:) = 0._wp ; vbdy_e(:) = 0._wp
	486	ubdy_n(:) = 0._wp ; vbdy_n(:) = 0._wp
	487	ubdy_s(:) = 0._wp ; vbdy_s(:) = 0._wp
[8993]	488	CALL Agrif_Bc_variable( unb_id, procname=interpunb )
	489	CALL Agrif_Bc_variable( vnb_id, procname=interpvnb )
[4486]	490	ENDIF
	491	Agrif_UseSpecialValue = .FALSE.
[5656]	492	!
[4486]	493	END SUBROUTINE Agrif_dta_ts
	494
[6140]	495
[2486]	496	SUBROUTINE Agrif_ssh( kt )
	497	!!----------------------------------------------------------------------
[8993]	498	!! * ROUTINE Agrif_ssh *
[2486]	499	!!----------------------------------------------------------------------
	500	INTEGER, INTENT(in) :: kt
	501	!!
[8993]	502	INTEGER :: ji, jj
[2486]	503	!!----------------------------------------------------------------------
[6140]	504	!
[2486]	505	IF( Agrif_Root() ) RETURN
[8993]	506	!
	507	! Linear interpolation in time of sea level
[6140]	508	!
[8993]	509	Agrif_SpecialValue = 0._wp
	510	Agrif_UseSpecialValue = .TRUE.
	511	CALL Agrif_Bc_variable(sshn_id, procname=interpsshn )
	512	Agrif_UseSpecialValue = .FALSE.
	513	!
[2486]	514	IF((nbondi == -1).OR.(nbondi == 2)) THEN
[8993]	515	DO jj=1,jpj
	516	ssha(2,jj) = hbdy_w(jj)
	517	END DO
[2486]	518	ENDIF
[6140]	519	!
[2486]	520	IF((nbondi == 1).OR.(nbondi == 2)) THEN
[8993]	521	DO jj=1,jpj
	522	ssha(nlci-1,jj) = hbdy_e(jj)
	523	END DO
[2486]	524	ENDIF
[6140]	525	!
[2486]	526	IF((nbondj == -1).OR.(nbondj == 2)) THEN
[8993]	527	DO ji=1,jpi
	528	ssha(ji,2) = hbdy_s(ji)
	529	END DO
[2486]	530	ENDIF
[6140]	531	!
[2486]	532	IF((nbondj == 1).OR.(nbondj == 2)) THEN
[8993]	533	DO ji=1,jpi
	534	ssha(ji,nlcj-1) = hbdy_n(ji)
	535	END DO
[2486]	536	ENDIF
[6140]	537	!
[2486]	538	END SUBROUTINE Agrif_ssh
	539
[6140]	540
[4486]	541	SUBROUTINE Agrif_ssh_ts( jn )
[4292]	542	!!----------------------------------------------------------------------
	543	!! * ROUTINE Agrif_ssh_ts *
	544	!!----------------------------------------------------------------------
[4486]	545	INTEGER, INTENT(in) :: jn
[4292]	546	!!
[4486]	547	INTEGER :: ji,jj
[4292]	548	!!----------------------------------------------------------------------
[6140]	549	!
[8993]	550	!
	551	IF( Agrif_Root() ) RETURN
	552	!
[4292]	553	IF((nbondi == -1).OR.(nbondi == 2)) THEN
[6140]	554	DO jj = 1, jpj
[4486]	555	ssha_e(2,jj) = hbdy_w(jj)
	556	END DO
[4292]	557	ENDIF
[6140]	558	!
[4292]	559	IF((nbondi == 1).OR.(nbondi == 2)) THEN
[6140]	560	DO jj = 1, jpj
[4486]	561	ssha_e(nlci-1,jj) = hbdy_e(jj)
	562	END DO
[4292]	563	ENDIF
[6140]	564	!
[4292]	565	IF((nbondj == -1).OR.(nbondj == 2)) THEN
[6140]	566	DO ji = 1, jpi
[4486]	567	ssha_e(ji,2) = hbdy_s(ji)
	568	END DO
[4292]	569	ENDIF
[6140]	570	!
[4292]	571	IF((nbondj == 1).OR.(nbondj == 2)) THEN
[6140]	572	DO ji = 1, jpi
[4486]	573	ssha_e(ji,nlcj-1) = hbdy_n(ji)
	574	END DO
[4292]	575	ENDIF
[6140]	576	!
[4292]	577	END SUBROUTINE Agrif_ssh_ts
	578
[8993]	579	# if defined key_zdftke \|\| defined key_zdfgls
[6140]	580
[8993]	581	SUBROUTINE Agrif_avm
[4292]	582	!!----------------------------------------------------------------------
[8993]	583	!! * ROUTINE Agrif_avm *
[5656]	584	!!----------------------------------------------------------------------
	585	REAL(wp) :: zalpha
[6140]	586	!!----------------------------------------------------------------------
[5656]	587	!
[8993]	588	IF( Agrif_Root() ) RETURN
[6140]	589	!
[8993]	590	! zalpha = REAL( Agrif_NbStepint() + Agrif_IRhot() - 1, wp ) / REAL( Agrif_IRhot(), wp )
	591	! IF( zalpha > 1. ) zalpha = 1.
	592	zalpha = 1._wp ! JC: proper time interpolation impossible
	593	! => use last available value from parent
	594	!
[5656]	595	Agrif_SpecialValue = 0.e0
	596	Agrif_UseSpecialValue = .TRUE.
[6140]	597	!
[5656]	598	CALL Agrif_Bc_variable(avm_id ,calledweight=zalpha, procname=interpavm)
[6140]	599	!
[5656]	600	Agrif_UseSpecialValue = .FALSE.
	601	!
[8993]	602	END SUBROUTINE Agrif_avm
[6140]	603
[5656]	604	# endif
	605
[6140]	606	SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
	607	!!----------------------------------------------------------------------
[5656]	608	!! * ROUTINE interptsn *
[6140]	609	!!----------------------------------------------------------------------
	610	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: ptab
	611	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
	612	LOGICAL , INTENT(in ) :: before
	613	INTEGER , INTENT(in ) :: nb , ndir
[5656]	614	!
[8998]	615	INTEGER :: ji, jj, jk, jn, iref, jref ! dummy loop indices
	616	INTEGER :: imin, imax, jmin, jmax, N_in, N_out
[5656]	617	REAL(wp) :: zrhox , zalpha1, zalpha2, zalpha3
	618	REAL(wp) :: zalpha4, zalpha5, zalpha6, zalpha7
[8998]	619	LOGICAL :: western_side, eastern_side,northern_side,southern_side
	620	! vertical interpolation:
	621	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: ptab_child
	622	REAL(wp), DIMENSION(k1:k2,n1:n2-1) :: tabin
	623	REAL(wp), DIMENSION(k1:k2) :: h_in
	624	REAL(wp), DIMENSION(1:jpk) :: h_out(1:jpk)
	625	REAL(wp) :: h_diff, zrhoxy
	626
	627	zrhoxy = Agrif_rhox()*Agrif_rhoy()
[5656]	628	IF (before) THEN
[8998]	629	DO jn = 1,jpts
	630	DO jk=k1,k2
	631	DO jj=j1,j2
	632	DO ji=i1,i2
	633	ptab(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)
	634	END DO
	635	END DO
	636	END DO
	637	END DO
	638
	639	# if defined key_vertical
	640	DO jk=k1,k2
	641	DO jj=j1,j2
	642	DO ji=i1,i2
	643	ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk)
	644	END DO
	645	END DO
	646	END DO
	647	# endif
	648	ELSE
	649
	650	western_side = (nb == 1).AND.(ndir == 1)
	651	eastern_side = (nb == 1).AND.(ndir == 2)
	652	southern_side = (nb == 2).AND.(ndir == 1)
	653	northern_side = (nb == 2).AND.(ndir == 2)
	654
	655	# if defined key_vertical
	656	DO jj=j1,j2
	657	DO ji=i1,i2
	658	iref = ji
	659	jref = jj
	660	if(western_side) iref=MAX(2,ji)
	661	if(eastern_side) iref=MIN(nlci-1,ji)
	662	if(southern_side) jref=MAX(2,jj)
	663	if(northern_side) jref=MIN(nlcj-1,jj)
	664	N_in = 0
	665	DO jk=k1,k2 !k2 = jpk of parent grid
	666	IF (ptab(ji,jj,jk,n2) == 0) EXIT
	667	N_in = N_in + 1
	668	tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)
	669	h_in(N_in) = ptab(ji,jj,jk,n2)
	670	END DO
	671	N_out = 0
	672	DO jk=1,jpk ! jpk of child grid
	673	IF (tmask(iref,jref,jk) == 0) EXIT
	674	N_out = N_out + 1
	675	h_out(jk) = e3t_n(iref,jref,jk)
	676	ENDDO
	677	IF (N_in > 0) THEN
	678	DO jn=1,jpts
	679	call reconstructandremap(tabin(1:N_in,jn),h_in,ptab_child(ji,jj,1:N_out,jn),h_out,N_in,N_out)
	680	ENDDO
	681	ENDIF
	682	ENDDO
	683	ENDDO
	684	# else
	685	ptab_child(i1:i2,j1:j2,1:jpk,1:jpts) = ptab(i1:i2,j1:j2,1:jpk,1:jpts)
	686	# endif
	687	!
[8993]	688	IF (lk_agrif_clp) THEN
	689	DO jn = 1, jpts
	690	DO jk = 1, jpkm1
	691	DO ji = i1,i2
	692	DO jj = j1,j2
[8998]	693	tsa(ji,jj,jk,jn) = ptab_child(ji,jj,jk,jn)
[8993]	694	END DO
	695	END DO
	696	END DO
	697	END DO
	698	return
	699	ENDIF
[5656]	700	!
	701	zrhox = Agrif_Rhox()
	702	!
	703	zalpha1 = ( zrhox - 1. ) * 0.5
	704	zalpha2 = 1. - zalpha1
	705	!
	706	zalpha3 = ( zrhox - 1. ) / ( zrhox + 1. )
	707	zalpha4 = 1. - zalpha3
	708	!
	709	zalpha6 = 2. * ( zrhox - 1. ) / ( zrhox + 1. )
	710	zalpha7 = - ( zrhox - 1. ) / ( zrhox + 3. )
	711	zalpha5 = 1. - zalpha6 - zalpha7
	712	!
	713	imin = i1
	714	imax = i2
	715	jmin = j1
	716	jmax = j2
	717	!
	718	! Remove CORNERS
	719	IF((nbondj == -1).OR.(nbondj == 2)) jmin = 3
	720	IF((nbondj == +1).OR.(nbondj == 2)) jmax = nlcj-2
	721	IF((nbondi == -1).OR.(nbondi == 2)) imin = 3
	722	IF((nbondi == +1).OR.(nbondi == 2)) imax = nlci-2
	723	!
[6140]	724	IF( eastern_side ) THEN
[5656]	725	DO jn = 1, jpts
[8998]	726	tsa(nlci,j1:j2,1:jpk,jn) = zalpha1 * ptab_child(nlci,j1:j2,1:jpk,jn) + zalpha2 * ptab_child(nlci-1,j1:j2,1:jpk,jn)
[5656]	727	DO jk = 1, jpkm1
	728	DO jj = jmin,jmax
[6140]	729	IF( umask(nlci-2,jj,jk) == 0._wp ) THEN
[5656]	730	tsa(nlci-1,jj,jk,jn) = tsa(nlci,jj,jk,jn) * tmask(nlci-1,jj,jk)
	731	ELSE
	732	tsa(nlci-1,jj,jk,jn)=(zalpha4tsa(nlci,jj,jk,jn)+zalpha3tsa(nlci-2,jj,jk,jn))*tmask(nlci-1,jj,jk)
[6140]	733	IF( un(nlci-2,jj,jk) > 0._wp ) THEN
[5656]	734	tsa(nlci-1,jj,jk,jn)=( zalpha6tsa(nlci-2,jj,jk,jn)+zalpha5tsa(nlci,jj,jk,jn) &
	735	+ zalpha7tsa(nlci-3,jj,jk,jn) ) tmask(nlci-1,jj,jk)
	736	ENDIF
	737	ENDIF
	738	END DO
	739	END DO
[5930]	740	tsa(nlci,j1:j2,k1:k2,jn) = 0._wp
[6140]	741	END DO
[5656]	742	ENDIF
	743	!
	744	IF( northern_side ) THEN
	745	DO jn = 1, jpts
[8998]	746	tsa(i1:i2,nlcj,1:jpk,jn) = zalpha1 * ptab_child(i1:i2,nlcj,1:jpk,jn) + zalpha2 * ptab_child(i1:i2,nlcj-1,1:jpk,jn)
[5656]	747	DO jk = 1, jpkm1
	748	DO ji = imin,imax
[6140]	749	IF( vmask(ji,nlcj-2,jk) == 0._wp ) THEN
[5656]	750	tsa(ji,nlcj-1,jk,jn) = tsa(ji,nlcj,jk,jn) * tmask(ji,nlcj-1,jk)
	751	ELSE
	752	tsa(ji,nlcj-1,jk,jn)=(zalpha4tsa(ji,nlcj,jk,jn)+zalpha3tsa(ji,nlcj-2,jk,jn))*tmask(ji,nlcj-1,jk)
[6140]	753	IF (vn(ji,nlcj-2,jk) > 0._wp ) THEN
[5656]	754	tsa(ji,nlcj-1,jk,jn)=( zalpha6tsa(ji,nlcj-2,jk,jn)+zalpha5tsa(ji,nlcj,jk,jn) &
	755	+ zalpha7tsa(ji,nlcj-3,jk,jn) ) tmask(ji,nlcj-1,jk)
	756	ENDIF
	757	ENDIF
	758	END DO
	759	END DO
[5930]	760	tsa(i1:i2,nlcj,k1:k2,jn) = 0._wp
[6140]	761	END DO
[5656]	762	ENDIF
	763	!
[6140]	764	IF( western_side ) THEN
[5656]	765	DO jn = 1, jpts
[8998]	766	tsa(1,j1:j2,1:jpk,jn) = zalpha1 * ptab_child(1,j1:j2,1:jpk,jn) + zalpha2 * ptab_child(2,j1:j2,1:jpk,jn)
[5656]	767	DO jk = 1, jpkm1
	768	DO jj = jmin,jmax
[6140]	769	IF( umask(2,jj,jk) == 0._wp ) THEN
[5656]	770	tsa(2,jj,jk,jn) = tsa(1,jj,jk,jn) * tmask(2,jj,jk)
	771	ELSE
	772	tsa(2,jj,jk,jn)=(zalpha4tsa(1,jj,jk,jn)+zalpha3tsa(3,jj,jk,jn))*tmask(2,jj,jk)
[6140]	773	IF( un(2,jj,jk) < 0._wp ) THEN
[5656]	774	tsa(2,jj,jk,jn)=(zalpha6tsa(3,jj,jk,jn)+zalpha5tsa(1,jj,jk,jn)+zalpha7tsa(4,jj,jk,jn))tmask(2,jj,jk)
	775	ENDIF
	776	ENDIF
	777	END DO
	778	END DO
[5930]	779	tsa(1,j1:j2,k1:k2,jn) = 0._wp
[5656]	780	END DO
	781	ENDIF
	782	!
	783	IF( southern_side ) THEN
	784	DO jn = 1, jpts
[8998]	785	tsa(i1:i2,1,1:jpk,jn) = zalpha1 * ptab_child(i1:i2,1,1:jpk,jn) + zalpha2 * ptab_child(i1:i2,2,1:jpk,jn)
[6140]	786	DO jk = 1, jpk
[5656]	787	DO ji=imin,imax
[6140]	788	IF( vmask(ji,2,jk) == 0._wp ) THEN
[5656]	789	tsa(ji,2,jk,jn)=tsa(ji,1,jk,jn) * tmask(ji,2,jk)
	790	ELSE
	791	tsa(ji,2,jk,jn)=(zalpha4tsa(ji,1,jk,jn)+zalpha3tsa(ji,3,jk,jn))*tmask(ji,2,jk)
[6140]	792	IF( vn(ji,2,jk) < 0._wp ) THEN
[5656]	793	tsa(ji,2,jk,jn)=(zalpha6tsa(ji,3,jk,jn)+zalpha5tsa(ji,1,jk,jn)+zalpha7tsa(ji,4,jk,jn))tmask(ji,2,jk)
	794	ENDIF
	795	ENDIF
	796	END DO
	797	END DO
[5930]	798	tsa(i1:i2,1,k1:k2,jn) = 0._wp
[6140]	799	END DO
[5656]	800	ENDIF
	801	!
	802	! Treatment of corners
	803	!
	804	! East south
	805	IF ((eastern_side).AND.((nbondj == -1).OR.(nbondj == 2))) THEN
[8998]	806	tsa(nlci-1,2,:,:) = ptab_child(nlci-1,2,:,1:jpts)
[5656]	807	ENDIF
	808	! East north
	809	IF ((eastern_side).AND.((nbondj == 1).OR.(nbondj == 2))) THEN
[8998]	810	tsa(nlci-1,nlcj-1,:,:) = ptab_child(nlci-1,nlcj-1,:,1:jpts)
[5656]	811	ENDIF
	812	! West south
	813	IF ((western_side).AND.((nbondj == -1).OR.(nbondj == 2))) THEN
[8998]	814	tsa(2,2,:,:) = ptab_child(2,2,:,1:jpts)
[5656]	815	ENDIF
	816	! West north
	817	IF ((western_side).AND.((nbondj == 1).OR.(nbondj == 2))) THEN
[8998]	818	tsa(2,nlcj-1,:,:) = ptab_child(2,nlcj-1,:,1:jpts)
[5656]	819	ENDIF
	820	!
	821	ENDIF
	822	!
	823	END SUBROUTINE interptsn
	824
[6140]	825	SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before, nb, ndir )
[5656]	826	!!----------------------------------------------------------------------
[4292]	827	!! * ROUTINE interpsshn *
	828	!!----------------------------------------------------------------------
[6140]	829	INTEGER , INTENT(in ) :: i1, i2, j1, j2
	830	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
	831	LOGICAL , INTENT(in ) :: before
	832	INTEGER , INTENT(in ) :: nb , ndir
	833	!
[5656]	834	LOGICAL :: western_side, eastern_side,northern_side,southern_side
	835	!!----------------------------------------------------------------------
	836	!
	837	IF( before) THEN
	838	ptab(i1:i2,j1:j2) = sshn(i1:i2,j1:j2)
	839	ELSE
	840	western_side = (nb == 1).AND.(ndir == 1)
	841	eastern_side = (nb == 1).AND.(ndir == 2)
	842	southern_side = (nb == 2).AND.(ndir == 1)
	843	northern_side = (nb == 2).AND.(ndir == 2)
	844	IF(western_side) hbdy_w(j1:j2) = ptab(i1,j1:j2) * tmask(i1,j1:j2,1)
	845	IF(eastern_side) hbdy_e(j1:j2) = ptab(i1,j1:j2) * tmask(i1,j1:j2,1)
	846	IF(southern_side) hbdy_s(i1:i2) = ptab(i1:i2,j1) * tmask(i1:i2,j1,1)
	847	IF(northern_side) hbdy_n(i1:i2) = ptab(i1:i2,j1) * tmask(i1:i2,j1,1)
	848	ENDIF
	849	!
	850	END SUBROUTINE interpsshn
	851
[8998]	852	SUBROUTINE interpun( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before, nb, ndir )
[6140]	853	!!----------------------------------------------------------------------
[5656]	854	!! * ROUTINE interpun *
[8998]	855	!!---------------------------------------------
	856	!!
	857	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,m1,m2
	858	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
	859	LOGICAL, INTENT(in) :: before
	860	INTEGER, INTENT(in) :: nb , ndir
	861	!!
	862	INTEGER :: ji,jj,jk
	863	REAL(wp) :: zrhoy
	864	! vertical interpolation:
	865	REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
	866	REAL(wp), DIMENSION(1:jpk) :: h_out
	867	INTEGER :: N_in, N_out, iref
	868	REAL(wp) :: h_diff
	869	LOGICAL :: western_side, eastern_side
	870	!!---------------------------------------------
[5656]	871	!
[8998]	872	zrhoy = Agrif_rhoy()
	873	IF (before) THEN
	874	DO jk=1,jpk
	875	DO jj=j1,j2
	876	DO ji=i1,i2
[8999]	877	ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk)*umask(ji,jj,jk))
[8998]	878	# if defined key_vertical
	879	ptab(ji,jj,jk,2) = (umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk))
	880	# endif
	881	END DO
	882	END DO
[5656]	883	END DO
	884	ELSE
[8998]	885	zrhoy = Agrif_rhoy()
	886	# if defined key_vertical
	887	! VERTICAL REFINEMENT BEGIN
	888	western_side = (nb == 1).AND.(ndir == 1)
	889	eastern_side = (nb == 1).AND.(ndir == 2)
	890
	891	DO ji=i1,i2
	892	iref = ji
	893	IF (western_side) iref = MAX(2,ji)
	894	IF (eastern_side) iref = MIN(nlci-2,ji)
	895	DO jj=j1,j2
	896	N_in = 0
	897	DO jk=k1,k2
	898	IF (ptab(ji,jj,jk,2) == 0) EXIT
	899	N_in = N_in + 1
	900	tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2)
	901	h_in(N_in) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy)
	902	ENDDO
	903
	904	IF (N_in == 0) THEN
	905	ua(ji,jj,:) = 0._wp
	906	CYCLE
	907	ENDIF
	908
	909	N_out = 0
	910	DO jk=1,jpk
	911	if (umask(iref,jj,jk) == 0) EXIT
	912	N_out = N_out + 1
	913	h_out(N_out) = e3u_a(iref,jj,jk)
	914	ENDDO
	915
	916	IF (N_out == 0) THEN
	917	ua(ji,jj,:) = 0._wp
	918	CYCLE
	919	ENDIF
	920
	921	IF (N_in * N_out > 0) THEN
	922	h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
	923	! Should be able to remove the next IF/ELSEIF statement once scale factors are dealt with properly
	924	if (h_diff < -1.e4) then
	925	print *,'CHECK YOUR BATHY ...', h_diff, sum(h_out(1:N_out)), sum(h_in(1:N_in))
	926	! stop
	927	endif
	928	ENDIF
	929	call reconstructandremap(tabin(1:N_in),h_in(1:N_in),ua(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)
	930	ENDDO
	931	ENDDO
	932
	933	# else
[6140]	934	DO jk = 1, jpkm1
[5656]	935	DO jj=j1,j2
[8998]	936	ua(i1:i2,jj,jk) = ptab(i1:i2,jj,jk,1) / ( zrhoy * e2u(i1:i2,jj) * e3u_a(i1:i2,jj,jk) )
[5656]	937	END DO
	938	END DO
[8998]	939	# endif
	940
[5656]	941	ENDIF
	942	!
	943	END SUBROUTINE interpun
	944
[8998]	945	SUBROUTINE interpvn( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before, nb, ndir )
[6140]	946	!!----------------------------------------------------------------------
[5656]	947	!! * ROUTINE interpvn *
[6140]	948	!!----------------------------------------------------------------------
[5656]	949	!
[8998]	950	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,m1,m2
	951	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
	952	LOGICAL, INTENT(in) :: before
	953	INTEGER, INTENT(in) :: nb , ndir
	954	!
	955	INTEGER :: ji,jj,jk
	956	REAL(wp) :: zrhox
	957	! vertical interpolation:
	958	REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
	959	REAL(wp), DIMENSION(1:jpk) :: h_out
	960	INTEGER :: N_in, N_out, jref
	961	REAL(wp) :: h_diff
	962	LOGICAL :: northern_side,southern_side
	963	!!---------------------------------------------
[5656]	964	!
[8998]	965	IF (before) THEN
	966	DO jk=k1,k2
	967	DO jj=j1,j2
	968	DO ji=i1,i2
[8999]	969	ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)*vmask(ji,jj,jk))
[8998]	970	# if defined key_vertical
	971	ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v_n(ji,jj,jk)
	972	# endif
	973	END DO
	974	END DO
[5656]	975	END DO
[8998]	976	ELSE
	977	zrhox = Agrif_rhox()
	978	# if defined key_vertical
	979
	980	southern_side = (nb == 2).AND.(ndir == 1)
	981	northern_side = (nb == 2).AND.(ndir == 2)
	982
	983	DO jj=j1,j2
	984	jref = jj
	985	IF (southern_side) jref = MAX(2,jj)
	986	IF (northern_side) jref = MIN(nlcj-2,jj)
	987	DO ji=i1,i2
	988	N_in = 0
	989	DO jk=k1,k2
	990	if (ptab(ji,jj,jk,2) == 0) EXIT
	991	N_in = N_in + 1
	992	tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2)
	993	h_in(N_in) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox)
	994	END DO
	995	IF (N_in == 0) THEN
	996	va(ji,jj,:) = 0._wp
	997	CYCLE
	998	ENDIF
	999
	1000	N_out = 0
	1001	DO jk=1,jpk
	1002	if (vmask(ji,jref,jk) == 0) EXIT
	1003	N_out = N_out + 1
	1004	h_out(N_out) = e3v_a(ji,jref,jk)
	1005	END DO
	1006	IF (N_out == 0) THEN
	1007	va(ji,jj,:) = 0._wp
	1008	CYCLE
	1009	ENDIF
	1010	call reconstructandremap(tabin(1:N_in),h_in(1:N_in),va(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)
	1011	END DO
	1012	END DO
	1013	# else
[6140]	1014	DO jk = 1, jpkm1
[8998]	1015	va(i1:i2,j1:j2,jk) = ptab(i1:i2,j1:j2,jk,1) / ( zrhox * e1v(i1:i2,j1:j2) * e3v_a(i1:i2,j1:j2,jk) )
[5656]	1016	END DO
[8998]	1017	# endif
[5656]	1018	ENDIF
	1019	!
	1020	END SUBROUTINE interpvn
[636]	1021
[6140]	1022	SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before, nb, ndir )
[1605]	1023	!!----------------------------------------------------------------------
[5656]	1024	!! * ROUTINE interpunb *
[1605]	1025	!!----------------------------------------------------------------------
[6140]	1026	INTEGER , INTENT(in ) :: i1, i2, j1, j2
	1027	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
	1028	LOGICAL , INTENT(in ) :: before
	1029	INTEGER , INTENT(in ) :: nb , ndir
	1030	!
	1031	INTEGER :: ji, jj
	1032	REAL(wp) :: zrhoy, zrhot, zt0, zt1, ztcoeff
	1033	LOGICAL :: western_side, eastern_side,northern_side,southern_side
[1605]	1034	!!----------------------------------------------------------------------
[5656]	1035	!
[6140]	1036	IF( before ) THEN
	1037	ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu_n(i1:i2,j1:j2) * un_b(i1:i2,j1:j2)
[5656]	1038	ELSE
	1039	western_side = (nb == 1).AND.(ndir == 1)
	1040	eastern_side = (nb == 1).AND.(ndir == 2)
	1041	southern_side = (nb == 2).AND.(ndir == 1)
	1042	northern_side = (nb == 2).AND.(ndir == 2)
	1043	zrhoy = Agrif_Rhoy()
	1044	zrhot = Agrif_rhot()
	1045	! Time indexes bounds for integration
	1046	zt0 = REAL(Agrif_NbStepint() , wp) / zrhot
	1047	zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
	1048	! Polynomial interpolation coefficients:
	1049	IF( bdy_tinterp == 1 ) THEN
	1050	ztcoeff = zrhot * ( zt1*2._wp ( zt1 - 1._wp) &
[6140]	1051	& - zt0*2._wp ( zt0 - 1._wp) )
[5656]	1052	ELSEIF( bdy_tinterp == 2 ) THEN
	1053	ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp &
[6140]	1054	& - zt0 * ( zt0 - 1._wp)**2._wp )
[636]	1055
[5656]	1056	ELSE
	1057	ztcoeff = 1
	1058	ENDIF
	1059	!
	1060	IF(western_side) THEN
	1061	ubdy_w(j1:j2) = ubdy_w(j1:j2) + ztcoeff * ptab(i1,j1:j2)
	1062	ENDIF
	1063	IF(eastern_side) THEN
	1064	ubdy_e(j1:j2) = ubdy_e(j1:j2) + ztcoeff * ptab(i1,j1:j2)
	1065	ENDIF
	1066	IF(southern_side) THEN
	1067	ubdy_s(i1:i2) = ubdy_s(i1:i2) + ztcoeff * ptab(i1:i2,j1)
	1068	ENDIF
	1069	IF(northern_side) THEN
	1070	ubdy_n(i1:i2) = ubdy_n(i1:i2) + ztcoeff * ptab(i1:i2,j1)
	1071	ENDIF
	1072	!
	1073	IF( bdy_tinterp == 0 .OR. bdy_tinterp == 2) THEN
	1074	IF(western_side) THEN
[6140]	1075	ubdy_w(j1:j2) = ubdy_w(j1:j2) / (zrhoye2u(i1,j1:j2)) umask(i1,j1:j2,1)
[5656]	1076	ENDIF
	1077	IF(eastern_side) THEN
[6140]	1078	ubdy_e(j1:j2) = ubdy_e(j1:j2) / (zrhoye2u(i1,j1:j2)) umask(i1,j1:j2,1)
[5656]	1079	ENDIF
	1080	IF(southern_side) THEN
[6140]	1081	ubdy_s(i1:i2) = ubdy_s(i1:i2) / (zrhoye2u(i1:i2,j1)) umask(i1:i2,j1,1)
[5656]	1082	ENDIF
	1083	IF(northern_side) THEN
[6140]	1084	ubdy_n(i1:i2) = ubdy_n(i1:i2) / (zrhoye2u(i1:i2,j1)) umask(i1:i2,j1,1)
[5656]	1085	ENDIF
	1086	ENDIF
	1087	ENDIF
	1088	!
	1089	END SUBROUTINE interpunb
[636]	1090
[6140]	1091
	1092	SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before, nb, ndir )
[1605]	1093	!!----------------------------------------------------------------------
[5656]	1094	!! * ROUTINE interpvnb *
[1605]	1095	!!----------------------------------------------------------------------
[6140]	1096	INTEGER , INTENT(in ) :: i1, i2, j1, j2
	1097	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
	1098	LOGICAL , INTENT(in ) :: before
	1099	INTEGER , INTENT(in ) :: nb , ndir
	1100	!
	1101	INTEGER :: ji,jj
	1102	REAL(wp) :: zrhox, zrhot, zt0, zt1, ztcoeff
	1103	LOGICAL :: western_side, eastern_side,northern_side,southern_side
[1605]	1104	!!----------------------------------------------------------------------
[5656]	1105	!
[6140]	1106	IF( before ) THEN
	1107	ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv_n(i1:i2,j1:j2) * vn_b(i1:i2,j1:j2)
[5656]	1108	ELSE
	1109	western_side = (nb == 1).AND.(ndir == 1)
	1110	eastern_side = (nb == 1).AND.(ndir == 2)
	1111	southern_side = (nb == 2).AND.(ndir == 1)
	1112	northern_side = (nb == 2).AND.(ndir == 2)
	1113	zrhox = Agrif_Rhox()
	1114	zrhot = Agrif_rhot()
	1115	! Time indexes bounds for integration
	1116	zt0 = REAL(Agrif_NbStepint() , wp) / zrhot
	1117	zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
	1118	IF( bdy_tinterp == 1 ) THEN
	1119	ztcoeff = zrhot * ( zt1*2._wp ( zt1 - 1._wp) &
[6140]	1120	& - zt0*2._wp ( zt0 - 1._wp) )
[5656]	1121	ELSEIF( bdy_tinterp == 2 ) THEN
	1122	ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp &
[6140]	1123	& - zt0 * ( zt0 - 1._wp)**2._wp )
[5656]	1124	ELSE
	1125	ztcoeff = 1
	1126	ENDIF
	1127	!
	1128	IF(western_side) THEN
	1129	vbdy_w(j1:j2) = vbdy_w(j1:j2) + ztcoeff * ptab(i1,j1:j2)
	1130	ENDIF
	1131	IF(eastern_side) THEN
	1132	vbdy_e(j1:j2) = vbdy_e(j1:j2) + ztcoeff * ptab(i1,j1:j2)
	1133	ENDIF
	1134	IF(southern_side) THEN
	1135	vbdy_s(i1:i2) = vbdy_s(i1:i2) + ztcoeff * ptab(i1:i2,j1)
	1136	ENDIF
	1137	IF(northern_side) THEN
	1138	vbdy_n(i1:i2) = vbdy_n(i1:i2) + ztcoeff * ptab(i1:i2,j1)
	1139	ENDIF
	1140	!
	1141	IF( bdy_tinterp == 0 .OR. bdy_tinterp == 2) THEN
	1142	IF(western_side) THEN
	1143	vbdy_w(j1:j2) = vbdy_w(j1:j2) / (zrhox*e1v(i1,j1:j2)) &
	1144	& * vmask(i1,j1:j2,1)
	1145	ENDIF
	1146	IF(eastern_side) THEN
	1147	vbdy_e(j1:j2) = vbdy_e(j1:j2) / (zrhox*e1v(i1,j1:j2)) &
	1148	& * vmask(i1,j1:j2,1)
	1149	ENDIF
	1150	IF(southern_side) THEN
	1151	vbdy_s(i1:i2) = vbdy_s(i1:i2) / (zrhox*e1v(i1:i2,j1)) &
	1152	& * vmask(i1:i2,j1,1)
	1153	ENDIF
	1154	IF(northern_side) THEN
	1155	vbdy_n(i1:i2) = vbdy_n(i1:i2) / (zrhox*e1v(i1:i2,j1)) &
	1156	& * vmask(i1:i2,j1,1)
	1157	ENDIF
	1158	ENDIF
	1159	ENDIF
	1160	!
	1161	END SUBROUTINE interpvnb
[390]	1162
[6140]	1163
	1164	SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before, nb, ndir )
[1605]	1165	!!----------------------------------------------------------------------
[5656]	1166	!! * ROUTINE interpub2b *
[1605]	1167	!!----------------------------------------------------------------------
[6140]	1168	INTEGER , INTENT(in ) :: i1, i2, j1, j2
	1169	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
	1170	LOGICAL , INTENT(in ) :: before
	1171	INTEGER , INTENT(in ) :: nb , ndir
	1172	!
	1173	INTEGER :: ji,jj
	1174	REAL(wp) :: zrhot, zt0, zt1,zat
	1175	LOGICAL :: western_side, eastern_side,northern_side,southern_side
[1605]	1176	!!----------------------------------------------------------------------
[5656]	1177	IF( before ) THEN
[8993]	1178	IF ( ln_bt_fw ) THEN
	1179	ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2)
	1180	ELSE
	1181	ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
	1182	ENDIF
[5656]	1183	ELSE
	1184	western_side = (nb == 1).AND.(ndir == 1)
	1185	eastern_side = (nb == 1).AND.(ndir == 2)
	1186	southern_side = (nb == 2).AND.(ndir == 1)
	1187	northern_side = (nb == 2).AND.(ndir == 2)
	1188	zrhot = Agrif_rhot()
	1189	! Time indexes bounds for integration
	1190	zt0 = REAL(Agrif_NbStepint() , wp) / zrhot
	1191	zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
	1192	! Polynomial interpolation coefficients:
[6140]	1193	zat = zrhot * ( zt1*2._wp (-2._wp*zt1 + 3._wp) &
	1194	& - zt0*2._wp (-2._wp*zt0 + 3._wp) )
[5656]	1195	!
	1196	IF(western_side ) ubdy_w(j1:j2) = zat * ptab(i1,j1:j2)
	1197	IF(eastern_side ) ubdy_e(j1:j2) = zat * ptab(i1,j1:j2)
	1198	IF(southern_side) ubdy_s(i1:i2) = zat * ptab(i1:i2,j1)
	1199	IF(northern_side) ubdy_n(i1:i2) = zat * ptab(i1:i2,j1)
	1200	ENDIF
	1201	!
	1202	END SUBROUTINE interpub2b
[6140]	1203
[636]	1204
[6140]	1205	SUBROUTINE interpvb2b( ptab, i1, i2, j1, j2, before, nb, ndir )
[4292]	1206	!!----------------------------------------------------------------------
[5656]	1207	!! * ROUTINE interpvb2b *
[4292]	1208	!!----------------------------------------------------------------------
[6140]	1209	INTEGER , INTENT(in ) :: i1, i2, j1, j2
	1210	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
	1211	LOGICAL , INTENT(in ) :: before
	1212	INTEGER , INTENT(in ) :: nb , ndir
	1213	!
	1214	INTEGER :: ji,jj
	1215	REAL(wp) :: zrhot, zt0, zt1,zat
	1216	LOGICAL :: western_side, eastern_side,northern_side,southern_side
[4292]	1217	!!----------------------------------------------------------------------
[5656]	1218	!
	1219	IF( before ) THEN
[8993]	1220	IF ( ln_bt_fw ) THEN
	1221	ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)
	1222	ELSE
	1223	ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
	1224	ENDIF
[5656]	1225	ELSE
	1226	western_side = (nb == 1).AND.(ndir == 1)
	1227	eastern_side = (nb == 1).AND.(ndir == 2)
	1228	southern_side = (nb == 2).AND.(ndir == 1)
	1229	northern_side = (nb == 2).AND.(ndir == 2)
	1230	zrhot = Agrif_rhot()
	1231	! Time indexes bounds for integration
	1232	zt0 = REAL(Agrif_NbStepint() , wp) / zrhot
	1233	zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
	1234	! Polynomial interpolation coefficients:
[6140]	1235	zat = zrhot * ( zt1*2._wp (-2._wp*zt1 + 3._wp) &
	1236	& - zt0*2._wp (-2._wp*zt0 + 3._wp) )
[5656]	1237	!
[6140]	1238	IF(western_side ) vbdy_w(j1:j2) = zat * ptab(i1,j1:j2)
	1239	IF(eastern_side ) vbdy_e(j1:j2) = zat * ptab(i1,j1:j2)
	1240	IF(southern_side) vbdy_s(i1:i2) = zat * ptab(i1:i2,j1)
	1241	IF(northern_side) vbdy_n(i1:i2) = zat * ptab(i1:i2,j1)
[5656]	1242	ENDIF
	1243	!
	1244	END SUBROUTINE interpvb2b
[4292]	1245
[6140]	1246
	1247	SUBROUTINE interpe3t( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )
[5656]	1248	!!----------------------------------------------------------------------
	1249	!! * ROUTINE interpe3t *
	1250	!!----------------------------------------------------------------------
[6140]	1251	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
[5656]	1252	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
[6140]	1253	LOGICAL , INTENT(in ) :: before
	1254	INTEGER , INTENT(in ) :: nb , ndir
[5656]	1255	!
	1256	INTEGER :: ji, jj, jk
	1257	LOGICAL :: western_side, eastern_side, northern_side, southern_side
	1258	REAL(wp) :: ztmpmsk
	1259	!!----------------------------------------------------------------------
	1260	!
[6140]	1261	IF( before ) THEN
	1262	ptab(i1:i2,j1:j2,k1:k2) = tmask(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2)
[5656]	1263	ELSE
	1264	western_side = (nb == 1).AND.(ndir == 1)
	1265	eastern_side = (nb == 1).AND.(ndir == 2)
	1266	southern_side = (nb == 2).AND.(ndir == 1)
	1267	northern_side = (nb == 2).AND.(ndir == 2)
[4292]	1268
[6140]	1269	DO jk = k1, k2
	1270	DO jj = j1, j2
	1271	DO ji = i1, i2
[5656]	1272	! Get velocity mask at boundary edge points:
[6140]	1273	IF( western_side ) ztmpmsk = umask(ji ,jj ,1)
	1274	IF( eastern_side ) ztmpmsk = umask(nlci-2,jj ,1)
	1275	IF( northern_side) ztmpmsk = vmask(ji ,nlcj-2,1)
	1276	IF( southern_side) ztmpmsk = vmask(ji ,2 ,1)
	1277	!
	1278	IF( ABS( ptab(ji,jj,jk) - tmask(ji,jj,jk) * e3t_0(ji,jj,jk) )*ztmpmsk > 1.D-2) THEN
[5656]	1279	IF (western_side) THEN
	1280	WRITE(numout,*) 'ERROR bathymetry merge at the western border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
	1281	ELSEIF (eastern_side) THEN
	1282	WRITE(numout,*) 'ERROR bathymetry merge at the eastern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
	1283	ELSEIF (southern_side) THEN
	1284	WRITE(numout,*) 'ERROR bathymetry merge at the southern border ji,jj,jk', ji+nimpp-1,jj+njmpp-1,jk
	1285	ELSEIF (northern_side) THEN
	1286	WRITE(numout,*) 'ERROR bathymetry merge at the northen border ji,jj,jk', ji+nimpp-1,jj+njmpp-1,jk
	1287	ENDIF
[6140]	1288	WRITE(numout,*) ' ptab(ji,jj,jk), e3t(ji,jj,jk) ', ptab(ji,jj,jk), e3t_0(ji,jj,jk)
[5656]	1289	kindic_agr = kindic_agr + 1
	1290	ENDIF
	1291	END DO
	1292	END DO
	1293	END DO
[6140]	1294	!
[5656]	1295	ENDIF
	1296	!
	1297	END SUBROUTINE interpe3t
	1298
[6140]	1299
	1300	SUBROUTINE interpumsk( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )
[4292]	1301	!!----------------------------------------------------------------------
[5656]	1302	!! * ROUTINE interpumsk *
[4292]	1303	!!----------------------------------------------------------------------
[6140]	1304	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
	1305	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	1306	LOGICAL , INTENT(in ) :: before
	1307	INTEGER , INTENT(in ) :: nb , ndir
[5656]	1308	!
[6140]	1309	INTEGER :: ji, jj, jk
	1310	LOGICAL :: western_side, eastern_side
[4292]	1311	!!----------------------------------------------------------------------
[5656]	1312	!
[6140]	1313	IF( before ) THEN
	1314	ptab(i1:i2,j1:j2,k1:k2) = umask(i1:i2,j1:j2,k1:k2)
[5656]	1315	ELSE
[6140]	1316	western_side = (nb == 1).AND.(ndir == 1)
	1317	eastern_side = (nb == 1).AND.(ndir == 2)
	1318	DO jk = k1, k2
	1319	DO jj = j1, j2
	1320	DO ji = i1, i2
[5656]	1321	! Velocity mask at boundary edge points:
	1322	IF (ABS(ptab(ji,jj,jk) - umask(ji,jj,jk)) > 1.D-2) THEN
	1323	IF (western_side) THEN
	1324	WRITE(numout,*) 'ERROR with umask at the western border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
	1325	WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), umask(ji,jj,jk)
	1326	kindic_agr = kindic_agr + 1
	1327	ELSEIF (eastern_side) THEN
	1328	WRITE(numout,*) 'ERROR with umask at the eastern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
	1329	WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), umask(ji,jj,jk)
	1330	kindic_agr = kindic_agr + 1
	1331	ENDIF
	1332	ENDIF
	1333	END DO
	1334	END DO
[4292]	1335	END DO
[6140]	1336	!
[5656]	1337	ENDIF
	1338	!
	1339	END SUBROUTINE interpumsk
[4292]	1340
[6140]	1341
	1342	SUBROUTINE interpvmsk( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )
[4486]	1343	!!----------------------------------------------------------------------
[5656]	1344	!! * ROUTINE interpvmsk *
[4486]	1345	!!----------------------------------------------------------------------
[6140]	1346	INTEGER , INTENT(in ) :: i1,i2,j1,j2,k1,k2
	1347	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	1348	LOGICAL , INTENT(in ) :: before
	1349	INTEGER , INTENT(in ) :: nb , ndir
[5656]	1350	!
[6140]	1351	INTEGER :: ji, jj, jk
	1352	LOGICAL :: northern_side, southern_side
[4486]	1353	!!----------------------------------------------------------------------
[5656]	1354	!
[6140]	1355	IF( before ) THEN
	1356	ptab(i1:i2,j1:j2,k1:k2) = vmask(i1:i2,j1:j2,k1:k2)
[5656]	1357	ELSE
	1358	southern_side = (nb == 2).AND.(ndir == 1)
	1359	northern_side = (nb == 2).AND.(ndir == 2)
[6140]	1360	DO jk = k1, k2
	1361	DO jj = j1, j2
	1362	DO ji = i1, i2
[5656]	1363	! Velocity mask at boundary edge points:
	1364	IF (ABS(ptab(ji,jj,jk) - vmask(ji,jj,jk)) > 1.D-2) THEN
	1365	IF (southern_side) THEN
	1366	WRITE(numout,*) 'ERROR with vmask at the southern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
	1367	WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), vmask(ji,jj,jk)
	1368	kindic_agr = kindic_agr + 1
	1369	ELSEIF (northern_side) THEN
	1370	WRITE(numout,*) 'ERROR with vmask at the northern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
	1371	WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), vmask(ji,jj,jk)
	1372	kindic_agr = kindic_agr + 1
	1373	ENDIF
	1374	ENDIF
	1375	END DO
	1376	END DO
[4486]	1377	END DO
[6140]	1378	!
[5656]	1379	ENDIF
	1380	!
	1381	END SUBROUTINE interpvmsk
[4486]	1382
[8993]	1383	# if defined key_zdftke \|\| defined key_zdfgls
[5656]	1384
[8999]	1385	SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
[4486]	1386	!!----------------------------------------------------------------------
[5656]	1387	!! * ROUTINE interavm *
[4486]	1388	!!----------------------------------------------------------------------
[8999]	1389	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, m1, m2
	1390	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
[6140]	1391	LOGICAL , INTENT(in ) :: before
[8999]	1392	REAL(wp), DIMENSION(k1:k2) :: tabin
	1393	REAL(wp) :: h_in(k1:k2)
	1394	REAL(wp) :: h_out(1:jpk)
	1395	REAL(wp) :: zrhoxy
	1396	INTEGER :: N_in, N_out, ji, jj, jk
[4486]	1397	!!----------------------------------------------------------------------
[5656]	1398	!
[8999]	1399	zrhoxy = Agrif_rhox()*Agrif_rhoy()
	1400	IF (before) THEN
	1401	DO jk=k1,k2
	1402	DO jj=j1,j2
	1403	DO ji=i1,i2
	1404	ptab(ji,jj,jk,1) = avm_k(ji,jj,jk)
	1405	END DO
	1406	END DO
	1407	END DO
	1408	#ifdef key_vertical
	1409	DO jk=k1,k2
	1410	DO jj=j1,j2
	1411	DO ji=i1,i2
	1412	ptab(ji,jj,jk,2) = wmask(ji,jj,jk) * e1e2t(ji,jj) * e3w_n(ji,jj,jk)
	1413	END DO
	1414	END DO
	1415	END DO
	1416	#else
	1417	ptab(i1:i2,j1:j2,k1:k2,2) = 0._wp
	1418	#endif
	1419	ELSE
	1420	#ifdef key_vertical
	1421	avm_k(i1:i2,j1:j2,1:jpk) = 0.
	1422	DO jj=j1,j2
	1423	DO ji=i1,i2
	1424	N_in = 0
	1425	DO jk=k1,k2 !k2 = jpk of parent grid
	1426	IF (ptab(ji,jj,jk,2) == 0) EXIT
	1427	N_in = N_in + 1
	1428	tabin(jk) = ptab(ji,jj,jk,1)
	1429	h_in(N_in) = ptab(ji,jj,jk,2)/(e1e2t(ji,jj)*zrhoxy)
	1430	END DO
	1431	N_out = 0
	1432	DO jk=1,jpk ! jpk of child grid
	1433	IF (wmask(ji,jj,jk) == 0) EXIT
	1434	N_out = N_out + 1
	1435	h_out(jk) = e3t_n(ji,jj,jk)
	1436	ENDDO
	1437	IF (N_in > 0) THEN
	1438	CALL reconstructandremap(tabin(1:N_in),h_in,avm_k(ji,jj,1:N_out),h_out,N_in,N_out)
	1439	ENDIF
	1440	ENDDO
	1441	ENDDO
	1442	#else
	1443	avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1)
	1444	#endif
[5656]	1445	ENDIF
	1446	!
	1447	END SUBROUTINE interpavm
[4486]	1448
[8993]	1449	# endif /* key_zdftke \|\| key_zdfgls */
[4486]	1450
[390]	1451	#else
[1605]	1452	!!----------------------------------------------------------------------
	1453	!! Empty module no AGRIF zoom
	1454	!!----------------------------------------------------------------------
[636]	1455	CONTAINS
	1456	SUBROUTINE Agrif_OPA_Interp_empty
	1457	WRITE(,) 'agrif_opa_interp : You should not have seen this print! error?'
	1458	END SUBROUTINE Agrif_OPA_Interp_empty
[390]	1459	#endif
[1605]	1460
	1461	!!======================================================================
[636]	1462	END MODULE agrif_opa_interp

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_METO_MERCATOR_2017_agrif/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90 @ 9005

Download in other formats: