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.

agrif_opa_update.F90 in branches/2017/dev_r8624_AGRIF3_VVL/NEMOGCM/NEMO/NST_SRC – NEMO

Context Navigation

source: branches/2017/dev_r8624_AGRIF3_VVL/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90 @ 8901

Last change on this file since 8901 was 8901, checked in by jchanut, 6 years ago
AGRIF: Remove update frequency parameter from namelist - #1965
Property svn:keywords set to `Id`
File size: 41.2 KB

Rev	Line
[5656]	1	#define TWO_WAY /* TWO WAY NESTING */
[8803]	2	#undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/
	3	#undef VOL_REFLUX /* VOLUME REFLUXING*/
[5656]	4
[636]	5	MODULE agrif_opa_update
[7646]	6	#if defined key_agrif
[636]	7	USE par_oce
	8	USE oce
	9	USE dom_oce
[782]	10	USE agrif_oce
[2715]	11	USE in_out_manager ! I/O manager
	12	USE lib_mpp
[3294]	13	USE wrk_nemo
[5656]	14	USE zdf_oce ! vertical physics: ocean variables
[8741]	15	USE domvvl ! Need interpolation routines
[3294]	16
[636]	17	IMPLICIT NONE
	18	PRIVATE
[390]	19
[8803]	20	PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Update_Scales, Agrif_Update_vvl, Agrif_Update_ssh
[8741]	21
[5656]	22	# if defined key_zdftke
	23	PUBLIC Agrif_Update_Tke
	24	# endif
[1156]	25	!!----------------------------------------------------------------------
[5656]	26	!! NEMO/NST 3.6 , NEMO Consortium (2010)
[1156]	27	!! $Id$
[2528]	28	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1156]	29	!!----------------------------------------------------------------------
[636]	30	CONTAINS
	31
[8741]	32	SUBROUTINE Agrif_Update_Tra( )
[636]	33	!!---------------------------------------------
	34	!! * ROUTINE Agrif_Update_Tra *
	35	!!---------------------------------------------
[5656]	36	!
	37	IF (Agrif_Root()) RETURN
	38	!
	39	#if defined TWO_WAY
[8901]	40	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update tracers from grid Number',Agrif_Fixed()
[636]	41
[390]	42	Agrif_UseSpecialValueInUpdate = .TRUE.
	43	Agrif_SpecialValueFineGrid = 0.
[5656]	44	!
	45	# if ! defined DECAL_FEEDBACK
[8901]	46	CALL Agrif_Update_Variable(tsn_id, procname=updateTS)
	47	! near boundary update:
	48	! CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS)
[5656]	49	# else
[8901]	50	CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS)
	51	! near boundary update:
	52	! CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS)
[5656]	53	# endif
	54	!
[390]	55	Agrif_UseSpecialValueInUpdate = .FALSE.
[5656]	56	!
[390]	57	#endif
[5656]	58	!
[636]	59	END SUBROUTINE Agrif_Update_Tra
[390]	60
[8741]	61	SUBROUTINE Agrif_Update_Dyn( )
[636]	62	!!---------------------------------------------
	63	!! * ROUTINE Agrif_Update_Dyn *
	64	!!---------------------------------------------
[5656]	65	!
	66	IF (Agrif_Root()) RETURN
	67	!
[390]	68	#if defined TWO_WAY
[8901]	69	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update momentum from grid Number',Agrif_Fixed()
[390]	70
[5656]	71	Agrif_UseSpecialValueInUpdate = .FALSE.
	72	Agrif_SpecialValueFineGrid = 0.
	73	!
	74	# if ! defined DECAL_FEEDBACK
[8901]	75	CALL Agrif_Update_Variable(un_update_id,procname = updateU)
	76	CALL Agrif_Update_Variable(vn_update_id,procname = updateV)
	77	! near boundary update:
	78	! CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU)
	79	! CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)
[5656]	80	# else
[8901]	81	CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU)
	82	CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV)
	83	! near boundary update:
	84	! CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU)
	85	! CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV)
[5656]	86	# endif
[390]	87
[5656]	88	# if ! defined DECAL_FEEDBACK
	89	CALL Agrif_Update_Variable(e1u_id,procname = updateU2d)
	90	CALL Agrif_Update_Variable(e2v_id,procname = updateV2d)
	91	# else
	92	CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d)
	93	CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d)
	94	# endif
[8803]	95	!
[8898]	96	# if ! defined DECAL_FEEDBACK
	97	! Account for updated thicknesses at boundary edges
	98	IF (.NOT.ln_linssh) THEN
	99	CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_u_bdy)
	100	CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_v_bdy)
	101	ENDIF
	102	# endif
[8803]	103	!
[5930]	104	IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN
[4486]	105	! Update time integrated transports
[5656]	106	# if ! defined DECAL_FEEDBACK
[8803]	107	CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b)
	108	CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b)
[5656]	109	# else
[8803]	110	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b)
	111	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b)
[5656]	112	# endif
	113	END IF
[8803]	114	#endif
[5656]	115	!
[8803]	116	END SUBROUTINE Agrif_Update_Dyn
	117
	118	SUBROUTINE Agrif_Update_ssh( )
	119	!!---------------------------------------------
	120	!! * ROUTINE Agrif_Update_ssh *
	121	!!---------------------------------------------
	122	!
	123	IF (Agrif_Root()) RETURN
[5656]	124	!
[8803]	125	#if defined TWO_WAY
	126	!
[5656]	127	Agrif_UseSpecialValueInUpdate = .TRUE.
[636]	128	Agrif_SpecialValueFineGrid = 0.
[5656]	129	# if ! defined DECAL_FEEDBACK
	130	CALL Agrif_Update_Variable(sshn_id,procname = updateSSH)
	131	# else
	132	CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH)
	133	# endif
[8803]	134	!
[636]	135	Agrif_UseSpecialValueInUpdate = .FALSE.
[8803]	136	!
	137	# if defined DECAL_FEEDBACK && defined VOL_REFLUX
	138	IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN
	139	! Refluxing on ssh:
	140	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0, 0/),locupdate2=(/1, 1/),procname = reflux_sshu)
	141	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1, 1/),locupdate2=(/0, 0/),procname = reflux_sshv)
	142	END IF
	143	# endif
	144	!
[390]	145	#endif
[5656]	146	!
[8803]	147	END SUBROUTINE Agrif_Update_ssh
[636]	148
[5656]	149	# if defined key_zdftke
[6140]	150
[8762]	151	SUBROUTINE Agrif_Update_Tke( kt )
[636]	152	!!---------------------------------------------
[5656]	153	!! * ROUTINE Agrif_Update_Tke *
[636]	154	!!---------------------------------------------
[5656]	155	!!
[8762]	156	INTEGER, INTENT(in) :: kt
[8741]	157	!
	158	IF (Agrif_Root()) RETURN
[5656]	159	!
	160	IF( ( Agrif_NbStepint() .NE. 0 ) .AND. (kt /= 0) ) RETURN
	161	# if defined TWO_WAY
[636]	162
[5656]	163	Agrif_UseSpecialValueInUpdate = .TRUE.
	164	Agrif_SpecialValueFineGrid = 0.
[636]	165
[5656]	166	CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN )
	167	CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT )
	168	CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM )
	169
	170	Agrif_UseSpecialValueInUpdate = .FALSE.
	171
	172	# endif
	173
	174	END SUBROUTINE Agrif_Update_Tke
[6140]	175
[5656]	176	# endif /* key_zdftke */
	177
[8741]	178	SUBROUTINE Agrif_Update_vvl( )
[636]	179	!!---------------------------------------------
[8741]	180	!! * ROUTINE Agrif_Update_vvl *
	181	!!---------------------------------------------
	182	!
	183	IF (Agrif_Root()) RETURN
	184	!
	185	#if defined TWO_WAY
	186	!
	187	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step()
	188	!
	189	Agrif_UseSpecialValueInUpdate = .TRUE.
	190	Agrif_SpecialValueFineGrid = 0.
	191	!
[8803]	192	! No interface separation here, update vertical grid at T points
	193	! everywhere over the overlapping regions (one account for refluxing in that case):
	194	CALL Agrif_Update_Variable(e3t_id, procname=updatee3t)
[8741]	195	!
	196	Agrif_UseSpecialValueInUpdate = .FALSE.
	197	!
	198	CALL Agrif_ChildGrid_To_ParentGrid()
	199	CALL dom_vvl_update_UVF
	200	CALL Agrif_ParentGrid_To_ChildGrid()
	201	!
	202	#endif
	203	!
	204	END SUBROUTINE Agrif_Update_vvl
	205
	206	SUBROUTINE dom_vvl_update_UVF
	207	!!---------------------------------------------
	208	!! * ROUTINE dom_vvl_update_UVF *
	209	!!---------------------------------------------
	210	!!
	211	INTEGER :: jk
	212	REAL(wp):: zcoef
	213	!!---------------------------------------------
	214
	215	IF (lwp.AND.lk_agrif_debug) Write(,) 'Finalize e3 on grid Number', &
	216	& Agrif_Fixed(), 'Step', Agrif_Nb_Step()
	217
	218	! Save "old" scale factor (prior update) for subsequent asselin correction
[8803]	219	! of prognostic variables
	220	! -----------------------
[8741]	221	!
	222	e3u_a(:,:,:) = e3u_n(:,:,:)
	223	e3v_a(:,:,:) = e3v_n(:,:,:)
	224	! ua(:,:,:) = e3u_b(:,:,:)
	225	! va(:,:,:) = e3v_b(:,:,:)
	226	hu_a(:,:) = hu_n(:,:)
	227	hv_a(:,:) = hv_n(:,:)
	228
	229	! 1) NOW fields
	230	!--------------
	231
	232	! Vertical scale factor interpolations
	233	! ------------------------------------
	234	CALL dom_vvl_interpol( e3t_n(:,:,:), e3u_n(:,:,:) , 'U' )
	235	CALL dom_vvl_interpol( e3t_n(:,:,:), e3v_n(:,:,:) , 'V' )
	236	CALL dom_vvl_interpol( e3u_n(:,:,:), e3f_n(:,:,:) , 'F' )
	237
	238	CALL dom_vvl_interpol( e3u_n(:,:,:), e3uw_n(:,:,:), 'UW' )
	239	CALL dom_vvl_interpol( e3v_n(:,:,:), e3vw_n(:,:,:), 'VW' )
	240
	241	! Update total depths:
	242	! --------------------
	243	hu_n(:,:) = 0._wp ! Ocean depth at U-points
	244	hv_n(:,:) = 0._wp ! Ocean depth at V-points
	245	DO jk = 1, jpkm1
	246	hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk)
	247	hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk)
	248	END DO
	249	! ! Inverse of the local depth
	250	r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) )
	251	r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) )
	252
	253
	254	! 2) BEFORE fields:
	255	!------------------
[8803]	256	! IF ( (.NOT.(lk_agrif_fstep.AND.(neuler==0)).AND.(ln_dynspg_exp)) &
	257	! & .OR.(.NOT.(lk_agrif_fstep.AND.(neuler==0)).AND.(ln_dynspg_ts &
	258	! & .AND.(.NOT.ln_bt_fw)))) THEN
	259	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN
[8741]	260	!
	261	! Vertical scale factor interpolations
	262	! ------------------------------------
	263	CALL dom_vvl_interpol( e3t_b(:,:,:), e3u_b(:,:,:), 'U' )
	264	CALL dom_vvl_interpol( e3t_b(:,:,:), e3v_b(:,:,:), 'V' )
	265
	266	CALL dom_vvl_interpol( e3u_b(:,:,:), e3uw_b(:,:,:), 'UW' )
	267	CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' )
	268
	269	! Update total depths:
	270	! --------------------
	271	hu_b(:,:) = 0._wp ! Ocean depth at U-points
	272	hv_b(:,:) = 0._wp ! Ocean depth at V-points
	273	DO jk = 1, jpkm1
	274	hu_b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk)
	275	hv_b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk)
	276	END DO
	277	! ! Inverse of the local depth
	278	r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )
	279	r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) )
	280	ENDIF
	281	!
	282	END SUBROUTINE dom_vvl_update_UVF
	283
[8762]	284	SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
[8741]	285	!!---------------------------------------------
[636]	286	!! * ROUTINE updateT *
	287	!!---------------------------------------------
[3294]	288	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
	289	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
[5656]	290	LOGICAL, INTENT(in) :: before
	291	!!
[3294]	292	INTEGER :: ji,jj,jk,jn
[8741]	293	REAL(wp) :: ztb, ztnu, ztno
[5656]	294	!!---------------------------------------------
	295	!
[636]	296	IF (before) THEN
[3294]	297	DO jn = n1,n2
	298	DO jk=k1,k2
	299	DO jj=j1,j2
	300	DO ji=i1,i2
[8741]	301	!> jc tmp
	302	tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) * e3t_n(ji,jj,jk) / e3t_0(ji,jj,jk)
	303	! tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) * e3t_n(ji,jj,jk)
	304	!< jc tmp
[3294]	305	END DO
[636]	306	END DO
	307	END DO
	308	END DO
	309	ELSE
[8741]	310	!> jc tmp
	311	DO jn = n1,n2
	312	tabres(i1:i2,j1:j2,k1:k2,jn) = tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) &
	313	& * tmask(i1:i2,j1:j2,k1:k2)
	314	ENDDO
	315	!< jc tmp
[4491]	316	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
[5656]	317	! Add asselin part
[4491]	318	DO jn = n1,n2
	319	DO jk=k1,k2
	320	DO jj=j1,j2
	321	DO ji=i1,i2
	322	IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
[8741]	323	ztb = tsb(ji,jj,jk,jn) * e3t_b(ji,jj,jk) ! fse3t_b prior update should be used
	324	ztnu = tabres(ji,jj,jk,jn)
	325	ztno = tsn(ji,jj,jk,jn) * e3t_a(ji,jj,jk)
	326	tsb(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) ) &
	327	& * tmask(ji,jj,jk) / e3t_b(ji,jj,jk)
[4491]	328	ENDIF
	329	ENDDO
	330	ENDDO
	331	ENDDO
	332	ENDDO
	333	ENDIF
[3294]	334	DO jn = n1,n2
	335	DO jk=k1,k2
	336	DO jj=j1,j2
	337	DO ji=i1,i2
	338	IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
[8741]	339	tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / e3t_n(ji,jj,jk)
[3294]	340	END IF
	341	END DO
[636]	342	END DO
	343	END DO
	344	END DO
[8741]	345	!
	346	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	347	tsb(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2)
	348	ENDIF
	349	!
[636]	350	ENDIF
[5656]	351	!
[3294]	352	END SUBROUTINE updateTS
[390]	353
[6140]	354
[636]	355	SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before )
	356	!!---------------------------------------------
	357	!! * ROUTINE updateu *
	358	!!---------------------------------------------
[6140]	359	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
[636]	360	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
[6140]	361	LOGICAL , INTENT(in ) :: before
	362	!
[8741]	363	INTEGER :: ji, jj, jk
	364	REAL(wp) :: zrhoy, zub, zunu, zuno
[5656]	365	!!---------------------------------------------
	366	!
[6140]	367	IF( before ) THEN
[636]	368	zrhoy = Agrif_Rhoy()
[6140]	369	DO jk = k1, k2
	370	tabres(i1:i2,j1:j2,jk) = zrhoy * e2u(i1:i2,j1:j2) * e3u_n(i1:i2,j1:j2,jk) * un(i1:i2,j1:j2,jk)
[636]	371	END DO
	372	ELSE
[390]	373	DO jk=k1,k2
[636]	374	DO jj=j1,j2
	375	DO ji=i1,i2
[8741]	376	tabres(ji,jj,jk) = tabres(ji,jj,jk) * r1_e2u(ji,jj)
[4491]	377	!
	378	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
[8898]	379	zub = ub(ji,jj,jk) * e3u_b(ji,jj,jk) ! fse3t_b prior update should be used
[8741]	380	zuno = un(ji,jj,jk) * e3u_a(ji,jj,jk)
	381	zunu = tabres(ji,jj,jk)
	382	ub(ji,jj,jk) = ( zub + atfp * ( zunu - zuno) ) &
	383	& * umask(ji,jj,jk) / e3u_b(ji,jj,jk)
[4491]	384	ENDIF
	385	!
[8741]	386	un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / e3u_n(ji,jj,jk)
[636]	387	END DO
	388	END DO
	389	END DO
[8741]	390	!
	391	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	392	ub(i1:i2,j1:j2,k1:k2) = un(i1:i2,j1:j2,k1:k2)
	393	ENDIF
	394	!
[636]	395	ENDIF
[5656]	396	!
[636]	397	END SUBROUTINE updateu
[390]	398
[8898]	399	SUBROUTINE correct_u_bdy( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir )
	400	!!---------------------------------------------
	401	!! * ROUTINE correct_u_bdy *
	402	!!---------------------------------------------
	403	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
	404	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
	405	LOGICAL , INTENT(in ) :: before
	406	INTEGER , INTENT(in) :: nb, ndir
	407	!!
	408	LOGICAL :: western_side, eastern_side
	409	!
	410	INTEGER :: jj, jk
	411	REAL(wp) :: zcor
	412	!!---------------------------------------------
	413	!
	414	IF( .NOT.before ) THEN
	415	!
	416	western_side = (nb == 1).AND.(ndir == 1)
	417	eastern_side = (nb == 1).AND.(ndir == 2)
	418	!
	419	IF (western_side) THEN
	420	DO jj=j1,j2
	421	zcor = un_b(i1-1,jj) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - un_b(i1-1,jj)
	422	un_b(i1-1,jj) = un_b(i1-1,jj) + zcor
	423	DO jk=1,jpkm1
	424	un(i1-1,jj,jk) = un(i1-1,jj,jk) + zcor * umask(i1-1,jj,jk)
	425	END DO
	426	END DO
	427	ENDIF
	428	!
	429	IF (eastern_side) THEN
	430	DO jj=j1,j2
	431	zcor = un_b(i2+1,jj) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - un_b(i2+1,jj)
	432	un_b(i2+1,jj) = un_b(i2+1,jj) + zcor
	433	DO jk=1,jpkm1
	434	un(i2+1,jj,jk) = un(i2+1,jj,jk) + zcor * umask(i2+1,jj,jk)
	435	END DO
	436	END DO
	437	ENDIF
	438	!
	439	ENDIF
	440	!
	441	END SUBROUTINE correct_u_bdy
[6140]	442
[8898]	443
[8741]	444	SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before)
[636]	445	!!---------------------------------------------
	446	!! * ROUTINE updatev *
	447	!!---------------------------------------------
[8741]	448	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
	449	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
	450	LOGICAL , INTENT(in ) :: before
	451	!
	452	INTEGER :: ji, jj, jk
	453	REAL(wp) :: zrhox, zvb, zvnu, zvno
[5656]	454	!!---------------------------------------------
	455	!
[636]	456	IF (before) THEN
	457	zrhox = Agrif_Rhox()
[390]	458	DO jk=k1,k2
[636]	459	DO jj=j1,j2
	460	DO ji=i1,i2
[6140]	461	tabres(ji,jj,jk) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)
[636]	462	END DO
	463	END DO
	464	END DO
	465	ELSE
[390]	466	DO jk=k1,k2
[636]	467	DO jj=j1,j2
	468	DO ji=i1,i2
[8741]	469	tabres(ji,jj,jk) = tabres(ji,jj,jk) * r1_e1v(ji,jj)
[4491]	470	!
	471	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
[8898]	472	zvb = vb(ji,jj,jk) * e3v_b(ji,jj,jk) ! fse3t_b prior update should be used
[8741]	473	zvno = vn(ji,jj,jk) * e3v_a(ji,jj,jk)
	474	zvnu = tabres(ji,jj,jk)
	475	vb(ji,jj,jk) = ( zvb + atfp * ( zvnu - zvno) ) &
	476	& * vmask(ji,jj,jk) / e3v_b(ji,jj,jk)
[4491]	477	ENDIF
	478	!
[8741]	479	vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / e3v_n(ji,jj,jk)
[636]	480	END DO
	481	END DO
	482	END DO
[8741]	483	!
	484	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	485	vb(i1:i2,j1:j2,k1:k2) = vn(i1:i2,j1:j2,k1:k2)
	486	ENDIF
	487	!
[636]	488	ENDIF
[5656]	489	!
[636]	490	END SUBROUTINE updatev
[390]	491
[8898]	492	SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir )
	493	!!---------------------------------------------
	494	!! * ROUTINE correct_u_bdy *
	495	!!---------------------------------------------
	496	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
	497	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
	498	LOGICAL , INTENT(in ) :: before
	499	INTEGER , INTENT(in) :: nb, ndir
	500	!!
	501	LOGICAL :: southern_side, northern_side
	502	!
	503	INTEGER :: ji, jk
	504	REAL(wp) :: zcor
	505	!!---------------------------------------------
	506	!
	507	IF( .NOT.before ) THEN
	508	!
	509	southern_side = (nb == 2).AND.(ndir == 1)
	510	northern_side = (nb == 2).AND.(ndir == 2)
	511	!
	512	IF (southern_side) THEN
	513	DO ji=i1,i2
	514	zcor = vn_b(ji,j1-1) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vn_b(ji,j1-1)
	515	vn_b(ji,j1-1) = vn_b(ji,j1-1) + zcor
	516	DO jk=1,jpkm1
	517	vn(ji,j1-1,jk) = vn(ji,j1-1,jk) + zcor * vmask(ji,j1-1,jk)
	518	END DO
	519	END DO
	520	ENDIF
	521	!
	522	IF (northern_side) THEN
	523	DO ji=i1,i2
	524	zcor = vn_b(ji,j2+1) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vn_b(ji,j2+1)
	525	vn_b(ji,j2+1) = vn_b(ji,j2+1) + zcor
	526	DO jk=1,jpkm1
	527	vn(ji,j2+1,jk) = vn(ji,j2+1,jk) + zcor * vmask(ji,j2+1,jk)
	528	END DO
	529	END DO
	530	ENDIF
	531	!
	532	ENDIF
	533	!
	534	END SUBROUTINE correct_v_bdy
[6140]	535
[8898]	536
[636]	537	SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
	538	!!---------------------------------------------
	539	!! * ROUTINE updateu2d *
	540	!!---------------------------------------------
	541	INTEGER, INTENT(in) :: i1, i2, j1, j2
	542	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	543	LOGICAL, INTENT(in) :: before
[5656]	544	!!
[8741]	545	INTEGER :: ji, jj, jk
[636]	546	REAL(wp) :: zrhoy
[4491]	547	REAL(wp) :: zcorr
[5656]	548	!!---------------------------------------------
	549	!
[636]	550	IF (before) THEN
	551	zrhoy = Agrif_Rhoy()
	552	DO jj=j1,j2
	553	DO ji=i1,i2
[6140]	554	tabres(ji,jj) = zrhoy * un_b(ji,jj) * hu_n(ji,jj) * e2u(ji,jj)
[636]	555	END DO
	556	END DO
	557	ELSE
	558	DO jj=j1,j2
	559	DO ji=i1,i2
[8741]	560	tabres(ji,jj) = tabres(ji,jj) * r1_e2u(ji,jj)
[4491]	561	!
	562	! Update "now" 3d velocities:
[6140]	563	spgu(ji,jj) = 0._wp
[4491]	564	DO jk=1,jpkm1
[6140]	565	spgu(ji,jj) = spgu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk)
[4491]	566	END DO
	567	!
[8741]	568	zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj)
[4491]	569	DO jk=1,jpkm1
	570	un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)
	571	END DO
	572	!
[4486]	573	! Update barotropic velocities:
[5930]	574	IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
	575	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
[8741]	576	zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * r1_hu_b(ji,jj)
[5930]	577	ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1)
	578	END IF
[8741]	579	ENDIF
	580	un_b(ji,jj) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1)
[4491]	581	!
	582	! Correct "before" velocities to hold correct bt component:
	583	spgu(ji,jj) = 0.e0
	584	DO jk=1,jpkm1
[6140]	585	spgu(ji,jj) = spgu(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk)
[4491]	586	END DO
	587	!
[8741]	588	zcorr = ub_b(ji,jj) - spgu(ji,jj) * r1_hu_b(ji,jj)
[4491]	589	DO jk=1,jpkm1
	590	ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)
	591	END DO
	592	!
[636]	593	END DO
	594	END DO
[8741]	595	!
	596	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	597	ub_b(i1:i2,j1:j2) = un_b(i1:i2,j1:j2)
	598	ENDIF
[636]	599	ENDIF
[5656]	600	!
[636]	601	END SUBROUTINE updateu2d
[390]	602
[6140]	603
[636]	604	SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
	605	!!---------------------------------------------
	606	!! * ROUTINE updatev2d *
	607	!!---------------------------------------------
	608	INTEGER, INTENT(in) :: i1, i2, j1, j2
	609	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	610	LOGICAL, INTENT(in) :: before
[5656]	611	!!
[8741]	612	INTEGER :: ji, jj, jk
[636]	613	REAL(wp) :: zrhox
[4491]	614	REAL(wp) :: zcorr
[5656]	615	!!---------------------------------------------
	616	!
[636]	617	IF (before) THEN
	618	zrhox = Agrif_Rhox()
	619	DO jj=j1,j2
	620	DO ji=i1,i2
[6140]	621	tabres(ji,jj) = zrhox * vn_b(ji,jj) * hv_n(ji,jj) * e1v(ji,jj)
[636]	622	END DO
	623	END DO
	624	ELSE
	625	DO jj=j1,j2
	626	DO ji=i1,i2
[8741]	627	tabres(ji,jj) = tabres(ji,jj) * r1_e1v(ji,jj)
[4491]	628	!
	629	! Update "now" 3d velocities:
	630	spgv(ji,jj) = 0.e0
	631	DO jk=1,jpkm1
[6140]	632	spgv(ji,jj) = spgv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk)
[4491]	633	END DO
	634	!
[8741]	635	zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj)
[4491]	636	DO jk=1,jpkm1
	637	vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)
	638	END DO
	639	!
[4486]	640	! Update barotropic velocities:
[5930]	641	IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
	642	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
[8741]	643	zcorr = (tabres(ji,jj) - vn_b(ji,jj) * hv_a(ji,jj)) * r1_hv_b(ji,jj)
[5930]	644	vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1)
	645	END IF
	646	ENDIF
[8741]	647	vn_b(ji,jj) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1)
[4491]	648	!
	649	! Correct "before" velocities to hold correct bt component:
	650	spgv(ji,jj) = 0.e0
	651	DO jk=1,jpkm1
[6140]	652	spgv(ji,jj) = spgv(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk)
[4491]	653	END DO
	654	!
[8741]	655	zcorr = vb_b(ji,jj) - spgv(ji,jj) * r1_hv_b(ji,jj)
[4491]	656	DO jk=1,jpkm1
	657	vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)
	658	END DO
	659	!
[636]	660	END DO
	661	END DO
[8741]	662	!
	663	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	664	vb_b(i1:i2,j1:j2) = vn_b(i1:i2,j1:j2)
	665	ENDIF
	666	!
[636]	667	ENDIF
[5656]	668	!
[636]	669	END SUBROUTINE updatev2d
[390]	670
[5656]	671
[8762]	672	SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
[636]	673	!!---------------------------------------------
	674	!! * ROUTINE updateSSH *
	675	!!---------------------------------------------
	676	INTEGER, INTENT(in) :: i1, i2, j1, j2
	677	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	678	LOGICAL, INTENT(in) :: before
[5656]	679	!!
[636]	680	INTEGER :: ji, jj
[5656]	681	!!---------------------------------------------
	682	!
[636]	683	IF (before) THEN
	684	DO jj=j1,j2
	685	DO ji=i1,i2
[4486]	686	tabres(ji,jj) = sshn(ji,jj)
[636]	687	END DO
	688	END DO
	689	ELSE
[8898]	690	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
[8803]	691	DO jj=j1,j2
	692	DO ji=i1,i2
	693	sshb(ji,jj) = sshb(ji,jj) &
	694	& + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1)
[4491]	695	END DO
[8803]	696	END DO
[4491]	697	ENDIF
[6140]	698	!
[636]	699	DO jj=j1,j2
	700	DO ji=i1,i2
[4486]	701	sshn(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)
[636]	702	END DO
	703	END DO
[8741]	704	!
	705	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	706	sshb(i1:i2,j1:j2) = sshn(i1:i2,j1:j2)
	707	ENDIF
	708	!
[8762]	709
[636]	710	ENDIF
[5656]	711	!
[636]	712	END SUBROUTINE updateSSH
	713
[6140]	714
[8803]	715	SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
[4486]	716	!!---------------------------------------------
	717	!! * ROUTINE updateub2b *
	718	!!---------------------------------------------
	719	INTEGER, INTENT(in) :: i1, i2, j1, j2
	720	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	721	LOGICAL, INTENT(in) :: before
[8803]	722	!!
	723	INTEGER :: ji, jj
	724	REAL(wp) :: zrhoy, za1, zcor
	725	!!---------------------------------------------
	726	!
	727	IF (before) THEN
	728	zrhoy = Agrif_Rhoy()
	729	DO jj=j1,j2
	730	DO ji=i1,i2
	731	tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
	732	END DO
	733	END DO
	734	tabres = zrhoy * tabres
	735	ELSE
	736	!
	737	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
	738	!
	739	za1 = 1._wp / REAL(Agrif_rhot(), wp)
	740	DO jj=j1,j2
	741	DO ji=i1,i2
	742	zcor=tabres(ji,jj) - ub2_b(ji,jj)
	743	! Update time integrated fluxes also in case of multiply nested grids:
	744	ub2_i_b(ji,jj) = ub2_i_b(ji,jj) + za1 * zcor
	745	! Update corrective fluxes:
[8898]	746	un_bf(ji,jj) = un_bf(ji,jj) + zcor
[8803]	747	! Update half step back fluxes:
	748	ub2_b(ji,jj) = tabres(ji,jj)
	749	END DO
	750	END DO
	751	ENDIF
	752	!
	753	END SUBROUTINE updateub2b
	754
	755	SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
	756	!!---------------------------------------------
	757	!! * ROUTINE reflux_sshu *
	758	!!---------------------------------------------
	759	INTEGER, INTENT(in) :: i1, i2, j1, j2
	760	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	761	LOGICAL, INTENT(in) :: before
[8762]	762	INTEGER, INTENT(in) :: nb, ndir
[5656]	763	!!
[8762]	764	LOGICAL :: western_side, eastern_side
[4486]	765	INTEGER :: ji, jj
[8803]	766	REAL(wp) :: zrhoy, za1, zcor
[5656]	767	!!---------------------------------------------
	768	!
[4486]	769	IF (before) THEN
	770	zrhoy = Agrif_Rhoy()
	771	DO jj=j1,j2
	772	DO ji=i1,i2
	773	tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
	774	END DO
	775	END DO
	776	tabres = zrhoy * tabres
	777	ELSE
[8762]	778	!
	779	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
	780	!
	781	western_side = (nb == 1).AND.(ndir == 1)
	782	eastern_side = (nb == 1).AND.(ndir == 2)
	783	!
	784	IF (western_side) THEN
	785	DO jj=j1,j2
[8803]	786	zcor = rdt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj))
	787	sshn(i1 ,jj) = sshn(i1 ,jj) + zcor
	788	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i1 ,jj) = sshb(i1 ,jj) + atfp * zcor
[8762]	789	END DO
	790	ENDIF
	791	IF (eastern_side) THEN
	792	DO jj=j1,j2
[8803]	793	zcor = - rdt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj))
	794	sshn(i2+1,jj) = sshn(i2+1,jj) + zcor
	795	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i2+1,jj) = sshb(i2+1,jj) + atfp * zcor
[8762]	796	END DO
	797	ENDIF
[8803]	798	!
	799	ENDIF
	800	!
	801	END SUBROUTINE reflux_sshu
	802
	803	SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
	804	!!---------------------------------------------
	805	!! * ROUTINE updatevb2b *
	806	!!---------------------------------------------
	807	INTEGER, INTENT(in) :: i1, i2, j1, j2
	808	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	809	LOGICAL, INTENT(in) :: before
	810	!!
	811	INTEGER :: ji, jj
	812	REAL(wp) :: zrhox, za1, zcor
	813	!!---------------------------------------------
	814	!
	815	IF (before) THEN
	816	zrhox = Agrif_Rhox()
	817	DO jj=j1,j2
	818	DO ji=i1,i2
	819	tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj)
	820	END DO
	821	END DO
	822	tabres = zrhox * tabres
	823	ELSE
[8762]	824	!
[8803]	825	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
	826	!
[8741]	827	za1 = 1._wp / REAL(Agrif_rhot(), wp)
[4486]	828	DO jj=j1,j2
	829	DO ji=i1,i2
[8803]	830	zcor=tabres(ji,jj) - vb2_b(ji,jj)
[8762]	831	! Update time integrated fluxes also in case of multiply nested grids:
[8803]	832	vb2_i_b(ji,jj) = vb2_i_b(ji,jj) + za1 * zcor
	833	! Update corrective fluxes:
[8898]	834	vn_bf(ji,jj) = vn_bf(ji,jj) + zcor
[8762]	835	! Update half step back fluxes:
[8803]	836	vb2_b(ji,jj) = tabres(ji,jj)
[4486]	837	END DO
	838	END DO
	839	ENDIF
[5656]	840	!
[8803]	841	END SUBROUTINE updatevb2b
[4486]	842
[8803]	843	SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
[4486]	844	!!---------------------------------------------
[8803]	845	!! * ROUTINE reflux_sshv *
[4486]	846	!!---------------------------------------------
	847	INTEGER, INTENT(in) :: i1, i2, j1, j2
	848	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	849	LOGICAL, INTENT(in) :: before
[8762]	850	INTEGER, INTENT(in) :: nb, ndir
[5656]	851	!!
[8762]	852	LOGICAL :: southern_side, northern_side
[4486]	853	INTEGER :: ji, jj
[8803]	854	REAL(wp) :: zrhox, za1, zcor
[5656]	855	!!---------------------------------------------
	856	!
[4486]	857	IF (before) THEN
	858	zrhox = Agrif_Rhox()
	859	DO jj=j1,j2
	860	DO ji=i1,i2
	861	tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj)
	862	END DO
	863	END DO
	864	tabres = zrhox * tabres
	865	ELSE
[8762]	866	!
	867	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
	868	!
	869	southern_side = (nb == 2).AND.(ndir == 1)
	870	northern_side = (nb == 2).AND.(ndir == 2)
	871	!
	872	IF (southern_side) THEN
	873	DO ji=i1,i2
[8803]	874	zcor = rdt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1))
	875	sshn(ji,j1 ) = sshn(ji,j1 ) + zcor
	876	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j1 ) = sshb(ji,j1) + atfp * zcor
[8762]	877	END DO
	878	ENDIF
	879	IF (northern_side) THEN
	880	DO ji=i1,i2
[8803]	881	zcor = - rdt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2))
	882	sshn(ji,j2+1) = sshn(ji,j2+1) + zcor
	883	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j2+1) = sshb(ji,j2+1) + atfp * zcor
[8762]	884	END DO
	885	ENDIF
	886	!
[4486]	887	ENDIF
[5656]	888	!
[8803]	889	END SUBROUTINE reflux_sshv
[4486]	890
[5656]	891	SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
	892	! currently not used
	893	!!---------------------------------------------
	894	!! * ROUTINE updateT *
	895	!!---------------------------------------------
	896	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
	897	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
	898	LOGICAL, iNTENT(in) :: before
[6140]	899	!
[5656]	900	INTEGER :: ji,jj,jk
	901	REAL(wp) :: ztemp
[6140]	902	!!---------------------------------------------
[5656]	903
	904	IF (before) THEN
	905	DO jk=k1,k2
	906	DO jj=j1,j2
	907	DO ji=i1,i2
	908	tabres(ji,jj,jk,1) = e1t(ji,jj)e2t(ji,jj)tmask(ji,jj,jk)
	909	tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
	910	tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
	911	END DO
	912	END DO
	913	END DO
	914	tabres(:,:,:,1)=tabres(:,:,:,1)Agrif_Rhox()Agrif_Rhoy()
	915	tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
	916	tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
	917	ELSE
	918	DO jk=k1,k2
	919	DO jj=j1,j2
	920	DO ji=i1,i2
	921	IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
	922	print ,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)e2t(ji,jj)*tmask(ji,jj,jk)
	923	print ,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)tmask(ji,jj,jk)
	924	print ,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)tmask(ji,jj,jk)
	925	ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
	926	print *,'CORR = ',ztemp-1.
	927	print ,'NEW VALS = ',tabres(ji,jj,jk,2)ztemp,tabres(ji,jj,jk,3)*ztemp, &
	928	tabres(ji,jj,jk,2)ztemptabres(ji,jj,jk,3)*ztemp
	929	e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
	930	e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
	931	END IF
	932	END DO
	933	END DO
	934	END DO
	935	ENDIF
	936	!
	937	END SUBROUTINE update_scales
	938
	939	# if defined key_zdftke
[6140]	940
[5656]	941	SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
	942	!!---------------------------------------------
	943	!! * ROUTINE updateen *
	944	!!---------------------------------------------
	945	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	946	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	947	LOGICAL, INTENT(in) :: before
	948	!!---------------------------------------------
	949	!
	950	IF (before) THEN
	951	ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
	952	ELSE
	953	en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
	954	ENDIF
	955	!
	956	END SUBROUTINE updateEN
	957
	958
	959	SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
	960	!!---------------------------------------------
	961	!! * ROUTINE updateavt *
	962	!!---------------------------------------------
	963	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	964	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	965	LOGICAL, INTENT(in) :: before
	966	!!---------------------------------------------
	967	!
	968	IF (before) THEN
	969	ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
	970	ELSE
	971	avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
	972	ENDIF
	973	!
	974	END SUBROUTINE updateAVT
	975
	976
	977	SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
	978	!!---------------------------------------------
	979	!! * ROUTINE updateavm *
	980	!!---------------------------------------------
	981	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	982	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	983	LOGICAL, INTENT(in) :: before
	984	!!---------------------------------------------
	985	!
	986	IF (before) THEN
	987	ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
	988	ELSE
	989	avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
	990	ENDIF
	991	!
	992	END SUBROUTINE updateAVM
	993
	994	# endif /* key_zdftke */
	995
[8803]	996	SUBROUTINE updatee3t(ptab_dum, i1, i2, j1, j2, k1, k2, before )
[8741]	997	!!---------------------------------------------
	998	!! * ROUTINE updatee3t *
	999	!!---------------------------------------------
[8803]	1000	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ptab_dum
[8741]	1001	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	1002	LOGICAL, INTENT(in) :: before
[8803]	1003	!
	1004	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptab
[8741]	1005	INTEGER :: ji,jj,jk
	1006	REAL(wp) :: zcoef
	1007	!!---------------------------------------------
	1008	!
[8803]	1009	IF (.NOT.before) THEN
[8741]	1010	!
[8803]	1011	ALLOCATE(ptab(i1:i2,j1:j2,1:jpk))
[8741]	1012	!
[8803]	1013	! Update e3t from ssh (z* case only)
[8762]	1014	DO jk = 1, jpkm1
	1015	DO jj=j1,j2
	1016	DO ji=i1,i2
[8803]	1017	ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + sshn(ji,jj) &
	1018	& *ssmask(ji,jj)/(ht_0(ji,jj)-1._wp + ssmask(ji,jj)))
[8762]	1019	END DO
	1020	END DO
	1021	END DO
[8803]	1022	!
	1023	! 1) Updates at BEFORE time step:
	1024	! -------------------------------
	1025	!
[8741]	1026	! Save "old" scale factor (prior update) for subsequent asselin correction
[8803]	1027	! of prognostic variables
	1028	e3t_a(i1:i2,j1:j2,1:jpkm1) = e3t_n(i1:i2,j1:j2,1:jpkm1)
[8741]	1029
[8803]	1030	! One should also save e3t_b, but lacking of workspace...
	1031	! hdivn(i1:i2,j1:j2,1:jpkm1) = e3t_b(i1:i2,j1:j2,1:jpkm1)
	1032
[8898]	1033	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN
[8741]	1034	DO jk = 1, jpkm1
	1035	DO jj=j1,j2
	1036	DO ji=i1,i2
	1037	e3t_b(ji,jj,jk) = e3t_b(ji,jj,jk) &
	1038	& + atfp * ( ptab(ji,jj,jk) - e3t_n(ji,jj,jk) )
	1039	END DO
	1040	END DO
	1041	END DO
	1042	!
	1043	e3w_b (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + e3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
	1044	gdepw_b(i1:i2,j1:j2,1) = 0.0_wp
	1045	gdept_b(i1:i2,j1:j2,1) = 0.5_wp * e3w_b(i1:i2,j1:j2,1)
	1046	!
	1047	DO jk = 2, jpk
	1048	DO jj = j1,j2
	1049	DO ji = i1,i2
	1050	zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
	1051	e3w_b(ji,jj,jk) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * &
	1052	& ( e3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) ) &
	1053	& + 0.5_wp * tmask(ji,jj,jk) * &
	1054	& ( e3t_b(ji,jj,jk ) - e3t_0(ji,jj,jk ) )
	1055	gdepw_b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1)
	1056	gdept_b(ji,jj,jk) = zcoef * ( gdepw_b(ji,jj,jk ) + 0.5 * e3w_b(ji,jj,jk)) &
	1057	& + (1-zcoef) * ( gdept_b(ji,jj,jk-1) + e3w_b(ji,jj,jk))
	1058	END DO
	1059	END DO
	1060	END DO
	1061	!
	1062	ENDIF
	1063	!
	1064	! 2) Updates at NOW time step:
	1065	! ----------------------------
	1066	!
	1067	! Update vertical scale factor at T-points:
[8803]	1068	e3t_n(i1:i2,j1:j2,1:jpkm1) = ptab(i1:i2,j1:j2,1:jpkm1)
[8741]	1069	!
	1070	! Update total depth:
	1071	ht_n(i1:i2,j1:j2) = 0._wp
	1072	DO jk = 1, jpkm1
	1073	ht_n(i1:i2,j1:j2) = ht_n(i1:i2,j1:j2) + e3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk)
	1074	END DO
	1075	!
	1076	! Update vertical scale factor at W-points and depths:
	1077	e3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + e3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
	1078	gdept_n(i1:i2,j1:j2,1) = 0.5_wp * e3w_n(i1:i2,j1:j2,1)
	1079	gdepw_n(i1:i2,j1:j2,1) = 0.0_wp
	1080	gde3w_n(i1:i2,j1:j2,1) = gdept_n(i1:i2,j1:j2,1) - (ht_n(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh
	1081	!
	1082	DO jk = 2, jpk
	1083	DO jj = j1,j2
	1084	DO ji = i1,i2
	1085	zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
	1086	e3w_n(ji,jj,jk) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( e3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) ) &
	1087	& + 0.5_wp * tmask(ji,jj,jk) * ( e3t_n(ji,jj,jk ) - e3t_0(ji,jj,jk ) )
	1088	gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1)
	1089	gdept_n(ji,jj,jk) = zcoef * ( gdepw_n(ji,jj,jk ) + 0.5 * e3w_n(ji,jj,jk)) &
	1090	& + (1-zcoef) * ( gdept_n(ji,jj,jk-1) + e3w_n(ji,jj,jk))
	1091	gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - (ht_n(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh
	1092	END DO
	1093	END DO
	1094	END DO
	1095	!
	1096	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
	1097	e3t_b (i1:i2,j1:j2,1:jpk) = e3t_n (i1:i2,j1:j2,1:jpk)
	1098	e3w_b (i1:i2,j1:j2,1:jpk) = e3w_n (i1:i2,j1:j2,1:jpk)
	1099	gdepw_b(i1:i2,j1:j2,1:jpk) = gdepw_n(i1:i2,j1:j2,1:jpk)
	1100	gdept_b(i1:i2,j1:j2,1:jpk) = gdept_n(i1:i2,j1:j2,1:jpk)
	1101	ENDIF
	1102	!
[8803]	1103	DEALLOCATE(ptab)
[8741]	1104	ENDIF
	1105	!
	1106	END SUBROUTINE updatee3t
	1107
[390]	1108	#else
[636]	1109	CONTAINS
	1110	SUBROUTINE agrif_opa_update_empty
	1111	!!---------------------------------------------
	1112	!! * ROUTINE agrif_opa_update_empty *
	1113	!!---------------------------------------------
	1114	WRITE(,) 'agrif_opa_update : You should not have seen this print! error?'
	1115	END SUBROUTINE agrif_opa_update_empty
[390]	1116	#endif
[636]	1117	END MODULE agrif_opa_update
[4486]	1118

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

Download in other formats: