Context Navigation

agrif_oce_update.F90 @ 12377

Last change on this file since 12377 was 12377, checked in by acc, 5 years ago

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

Property svn:keywords set to Id

File size: 58.0 KB

Line
1	#undef DECAL_FEEDBACK /* SEPARATION of INTERFACES */
2	#undef DECAL_FEEDBACK_2D /* SEPARATION of INTERFACES (Barotropic mode) */
3	#undef VOL_REFLUX /* VOLUME REFLUXING*/
4
5	MODULE agrif_oce_update
6	!!======================================================================
7	!! * MODULE agrif_oce_interp *
8	!! AGRIF: update package for the ocean dynamics (OPA)
9	!!======================================================================
10	!! History : 2.0 ! 2002-06 (L. Debreu) Original code
11	!! 3.2 ! 2009-04 (R. Benshila)
12	!! 3.6 ! 2014-09 (R. Benshila)
13	!!----------------------------------------------------------------------
14	#if defined key_agrif
15	!!----------------------------------------------------------------------
16	!! 'key_agrif' AGRIF zoom
17	!!----------------------------------------------------------------------
18	USE par_oce
19	USE oce
20	USE dom_oce
21	USE zdf_oce ! vertical physics: ocean variables
22	USE agrif_oce
23	!
24	USE in_out_manager ! I/O manager
25	USE lib_mpp ! MPP library
26	USE domvvl ! Need interpolation routines
27	USE vremap ! Vertical remapping
28
29	IMPLICIT NONE
30	PRIVATE
31
32	PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Agrif_Update_vvl, Agrif_Update_ssh
33	PUBLIC Update_Scales
34
35	!!----------------------------------------------------------------------
36	!! NEMO/NST 4.0 , NEMO Consortium (2018)
37	!! $Id$
38	!! Software governed by the CeCILL license (see ./LICENSE)
39	!!----------------------------------------------------------------------
40	CONTAINS
41
42	SUBROUTINE Agrif_Update_Tra( )
43	!!----------------------------------------------------------------------
44	!! * ROUTINE Agrif_Update_Tra *
45	!!----------------------------------------------------------------------
46	!
47	IF (Agrif_Root()) RETURN
48	!
49	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update tracers from grid Number',Agrif_Fixed()
50
51	#if defined key_vertical
52	! Effect of this has to be carrefully checked
53	! depending on what the nesting tools ensure for
54	! volume conservation:
55	Agrif_UseSpecialValueInUpdate = .FALSE.
56	#else
57	Agrif_UseSpecialValueInUpdate = .TRUE.
58	#endif
59	Agrif_SpecialValueFineGrid = 0._wp
60	!
61	# if ! defined DECAL_FEEDBACK
62	CALL Agrif_Update_Variable(tsn_id, procname=updateTS)
63	! near boundary update:
64	! CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS)
65	# else
66	CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS)
67	! near boundary update:
68	! CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS)
69	# endif
70	!
71	Agrif_UseSpecialValueInUpdate = .FALSE.
72	!
73	!
74	END SUBROUTINE Agrif_Update_Tra
75
76	SUBROUTINE Agrif_Update_Dyn( )
77	!!----------------------------------------------------------------------
78	!! * ROUTINE Agrif_Update_Dyn *
79	!!----------------------------------------------------------------------
80	!
81	IF (Agrif_Root()) RETURN
82	!
83	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update momentum from grid Number',Agrif_Fixed()
84
85	Agrif_UseSpecialValueInUpdate = .FALSE.
86	Agrif_SpecialValueFineGrid = 0.
87	!
88	# if ! defined DECAL_FEEDBACK
89	CALL Agrif_Update_Variable(un_update_id,procname = updateU)
90	CALL Agrif_Update_Variable(vn_update_id,procname = updateV)
91	! near boundary update:
92	! CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU)
93	! CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)
94	# else
95	CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU)
96	CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV)
97	! near boundary update:
98	! CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU)
99	! CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV)
100	# endif
101
102	# if ! defined DECAL_FEEDBACK_2D
103	CALL Agrif_Update_Variable(e1u_id,procname = updateU2d)
104	CALL Agrif_Update_Variable(e2v_id,procname = updateV2d)
105	# else
106	CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d)
107	CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d)
108	# endif
109	!
110	# if ! defined DECAL_FEEDBACK_2D
111	! Account for updated thicknesses at boundary edges
112	IF (.NOT.ln_linssh) THEN
113	! CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_u_bdy)
114	! CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_v_bdy)
115	ENDIF
116	# endif
117	!
118	IF ( ln_dynspg_ts .AND. ln_bt_fw ) THEN
119	! Update time integrated transports
120	# if ! defined DECAL_FEEDBACK_2D
121	CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b)
122	CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b)
123	# else
124	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b)
125	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b)
126	# endif
127	END IF
128	!
129	END SUBROUTINE Agrif_Update_Dyn
130
131	SUBROUTINE Agrif_Update_ssh( )
132	!!---------------------------------------------
133	!! * ROUTINE Agrif_Update_ssh *
134	!!---------------------------------------------
135	!
136	IF (Agrif_Root()) RETURN
137	!
138	Agrif_UseSpecialValueInUpdate = .TRUE.
139	Agrif_SpecialValueFineGrid = 0.
140	# if ! defined DECAL_FEEDBACK_2D
141	CALL Agrif_Update_Variable(sshn_id,procname = updateSSH)
142	# else
143	CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH)
144	# endif
145	!
146	Agrif_UseSpecialValueInUpdate = .FALSE.
147	!
148	# if defined VOL_REFLUX
149	IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN
150	! Refluxing on ssh:
151	# if defined DECAL_FEEDBACK_2D
152	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0, 0/),locupdate2=(/1, 1/),procname = reflux_sshu)
153	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1, 1/),locupdate2=(/0, 0/),procname = reflux_sshv)
154	# else
155	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/-1,-1/),locupdate2=(/ 0, 0/),procname = reflux_sshu)
156	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ 0, 0/),locupdate2=(/-1,-1/),procname = reflux_sshv)
157	# endif
158	END IF
159	# endif
160	!
161	END SUBROUTINE Agrif_Update_ssh
162
163
164	SUBROUTINE Agrif_Update_Tke( )
165	!!---------------------------------------------
166	!! * ROUTINE Agrif_Update_Tke *
167	!!---------------------------------------------
168	!!
169	!
170	IF (Agrif_Root()) RETURN
171	!
172	Agrif_UseSpecialValueInUpdate = .TRUE.
173	Agrif_SpecialValueFineGrid = 0.
174
175	CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN )
176	CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT )
177	CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM )
178
179	Agrif_UseSpecialValueInUpdate = .FALSE.
180
181	END SUBROUTINE Agrif_Update_Tke
182
183
184	SUBROUTINE Agrif_Update_vvl( )
185	!!---------------------------------------------
186	!! * ROUTINE Agrif_Update_vvl *
187	!!---------------------------------------------
188	!
189	IF (Agrif_Root()) RETURN
190	!
191	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step()
192	!
193	Agrif_UseSpecialValueInUpdate = .TRUE.
194	Agrif_SpecialValueFineGrid = 0.
195	!
196	! No interface separation here, update vertical grid at T points
197	! everywhere over the overlapping regions (one account for refluxing in that case):
198	CALL Agrif_Update_Variable(e3t_id, procname=updatee3t)
199	!
200	Agrif_UseSpecialValueInUpdate = .FALSE.
201	!
202	CALL Agrif_ChildGrid_To_ParentGrid()
203	CALL dom_vvl_update_UVF
204	CALL Agrif_ParentGrid_To_ChildGrid()
205	!
206	END SUBROUTINE Agrif_Update_vvl
207
208	SUBROUTINE dom_vvl_update_UVF
209	!!---------------------------------------------
210	!! * ROUTINE dom_vvl_update_UVF *
211	!!---------------------------------------------
212	!!
213	INTEGER :: jk
214	REAL(wp):: zcoef
215	!!---------------------------------------------
216
217	IF (lwp.AND.lk_agrif_debug) Write(,) 'Finalize e3 on grid Number', &
218	& Agrif_Fixed(), 'Step', Agrif_Nb_Step()
219
220	! Save "old" scale factor (prior update) for subsequent asselin correction
221	! of prognostic variables
222	! -----------------------
223	!
224	e3u(:,:,:,Krhs_a) = e3u(:,:,:,Kmm_a)
225	e3v(:,:,:,Krhs_a) = e3v(:,:,:,Kmm_a)
226	! uu(:,:,:,Krhs_a) = e3u(:,:,:,Kbb_a)
227	! vv(:,:,:,Krhs_a) = e3v(:,:,:,Kbb_a)
228	hu(:,:,Krhs_a) = hu(:,:,Kmm_a)
229	hv(:,:,Krhs_a) = hv(:,:,Kmm_a)
230
231	! 1) NOW fields
232	!--------------
233
234	! Vertical scale factor interpolations
235	! ------------------------------------
236	CALL dom_vvl_interpol( e3t(:,:,:,Kmm_a), e3u(:,:,:,Kmm_a) , 'U' )
237	CALL dom_vvl_interpol( e3t(:,:,:,Kmm_a), e3v(:,:,:,Kmm_a) , 'V' )
238	CALL dom_vvl_interpol( e3u(:,:,:,Kmm_a), e3f(:,:,:) , 'F' )
239
240	CALL dom_vvl_interpol( e3u(:,:,:,Kmm_a), e3uw(:,:,:,Kmm_a), 'UW' )
241	CALL dom_vvl_interpol( e3v(:,:,:,Kmm_a), e3vw(:,:,:,Kmm_a), 'VW' )
242
243	! Update total depths:
244	! --------------------
245	hu(:,:,Kmm_a) = 0._wp ! Ocean depth at U-points
246	hv(:,:,Kmm_a) = 0._wp ! Ocean depth at V-points
247	DO jk = 1, jpkm1
248	hu(:,:,Kmm_a) = hu(:,:,Kmm_a) + e3u(:,:,jk,Kmm_a) * umask(:,:,jk)
249	hv(:,:,Kmm_a) = hv(:,:,Kmm_a) + e3v(:,:,jk,Kmm_a) * vmask(:,:,jk)
250	END DO
251	! ! Inverse of the local depth
252	r1_hu(:,:,Kmm_a) = ssumask(:,:) / ( hu(:,:,Kmm_a) + 1._wp - ssumask(:,:) )
253	r1_hv(:,:,Kmm_a) = ssvmask(:,:) / ( hv(:,:,Kmm_a) + 1._wp - ssvmask(:,:) )
254
255
256	! 2) BEFORE fields:
257	!------------------
258	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN
259	!
260	! Vertical scale factor interpolations
261	! ------------------------------------
262	CALL dom_vvl_interpol( e3t(:,:,:,Kbb_a), e3u(:,:,:,Kbb_a), 'U' )
263	CALL dom_vvl_interpol( e3t(:,:,:,Kbb_a), e3v(:,:,:,Kbb_a), 'V' )
264
265	CALL dom_vvl_interpol( e3u(:,:,:,Kbb_a), e3uw(:,:,:,Kbb_a), 'UW' )
266	CALL dom_vvl_interpol( e3v(:,:,:,Kbb_a), e3vw(:,:,:,Kbb_a), 'VW' )
267
268	! Update total depths:
269	! --------------------
270	hu(:,:,Kbb_a) = 0._wp ! Ocean depth at U-points
271	hv(:,:,Kbb_a) = 0._wp ! Ocean depth at V-points
272	DO jk = 1, jpkm1
273	hu(:,:,Kbb_a) = hu(:,:,Kbb_a) + e3u(:,:,jk,Kbb_a) * umask(:,:,jk)
274	hv(:,:,Kbb_a) = hv(:,:,Kbb_a) + e3v(:,:,jk,Kbb_a) * vmask(:,:,jk)
275	END DO
276	! ! Inverse of the local depth
277	r1_hu(:,:,Kbb_a) = ssumask(:,:) / ( hu(:,:,Kbb_a) + 1._wp - ssumask(:,:) )
278	r1_hv(:,:,Kbb_a) = ssvmask(:,:) / ( hv(:,:,Kbb_a) + 1._wp - ssvmask(:,:) )
279	ENDIF
280	!
281	END SUBROUTINE dom_vvl_update_UVF
282
283	#if defined key_vertical
284
285	SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
286	!!----------------------------------------------------------------------
287	!! * ROUTINE updateT *
288	!!---------------------------------------------
289	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
290	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
291	LOGICAL, INTENT(in) :: before
292	!!
293	INTEGER :: ji,jj,jk,jn
294	INTEGER :: N_in, N_out
295	REAL(wp) :: ztb, ztnu, ztno
296	REAL(wp) :: h_in(k1:k2)
297	REAL(wp) :: h_out(1:jpk)
298	REAL(wp) :: tabin(k1:k2,1:jpts)
299	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,1:jpts) :: tabres_child
300	!!---------------------------------------------
301	!
302	IF (before) THEN
303	!jc_alt
304	! AGRIF_SpecialValue = -999._wp
305	DO jn = n1,n2-1
306	DO jk=k1,k2
307	DO jj=j1,j2
308	DO ji=i1,i2
309	!jc_alt
310	! tabres(ji,jj,jk,jn) = (ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) ) &
311	! & * tmask(ji,jj,jk) + (tmask(ji,jj,jk)-1._wp) * 999._wp
312	tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a)
313	END DO
314	END DO
315	END DO
316	END DO
317	DO jk=k1,k2
318	DO jj=j1,j2
319	DO ji=i1,i2
320	!jc_alt
321	! tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) &
322	! & + (tmask(ji,jj,jk) - 1._wp) * 999._wp
323	tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
324	END DO
325	END DO
326	END DO
327	ELSE
328	tabres_child(:,:,:,:) = 0._wp
329	AGRIF_SpecialValue = 0._wp
330	DO jj=j1,j2
331	DO ji=i1,i2
332	N_in = 0
333	DO jk=k1,k2 !k2 = jpk of child grid
334	! jc_alt
335	! IF (tabres(ji,jj,jk,n2) < -900._wp ) EXIT
336	IF (tabres(ji,jj,jk,n2) == 0._wp ) EXIT
337	N_in = N_in + 1
338	tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1)/tabres(ji,jj,jk,n2)
339	h_in(N_in) = tabres(ji,jj,jk,n2)
340	ENDDO
341	N_out = 0
342	DO jk=1,jpk ! jpk of parent grid
343	IF (tmask(ji,jj,jk) == 0 ) EXIT ! TODO: Will not work with ISF
344	N_out = N_out + 1
345	h_out(N_out) = e3t(ji,jj,jk,Kmm_a)
346	ENDDO
347	IF (N_in*N_out > 0) THEN !Remove this?
348	CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),tabres_child(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts)
349	ENDIF
350	ENDDO
351	ENDDO
352
353	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
354	! Add asselin part
355	DO jn = 1,jpts
356	DO jk = 1, jpkm1
357	DO jj = j1, j2
358	DO ji = i1, i2
359	IF( tabres_child(ji,jj,jk,jn) /= 0._wp ) THEN
360	ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used
361	ztnu = tabres_child(ji,jj,jk,jn) * e3t(ji,jj,jk,Kmm_a)
362	ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a)
363	ts(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) &
364	& * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a)
365	ENDIF
366	END DO
367	END DO
368	END DO
369	END DO
370	ENDIF
371	DO jn = 1,jpts
372	DO jk = 1, jpkm1
373	DO jj = j1, j2
374	DO ji = i1, i2
375	IF( tabres_child(ji,jj,jk,jn) /= 0._wp ) THEN
376	ts(ji,jj,jk,jn,Kmm_a) = tabres_child(ji,jj,jk,jn)
377	END IF
378	END DO
379	END DO
380	END DO
381	END DO
382	!
383	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
384	ts(i1:i2,j1:j2,1:jpkm1,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,1:jpkm1,1:jpts,Kmm_a)
385	ENDIF
386	ENDIF
387	!
388	END SUBROUTINE updateTS
389
390	# else
391
392	SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
393	!!---------------------------------------------
394	!! * ROUTINE updateT *
395	!!---------------------------------------------
396	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
397	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
398	LOGICAL, INTENT(in) :: before
399	!!
400	INTEGER :: ji,jj,jk,jn
401	REAL(wp) :: ztb, ztnu, ztno
402	!!---------------------------------------------
403	!
404	IF (before) THEN
405	DO jn = 1,jpts
406	DO jk=k1,k2
407	DO jj=j1,j2
408	DO ji=i1,i2
409	!> jc tmp
410	tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) / e3t_0(ji,jj,jk)
411	! tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a)
412	!< jc tmp
413	END DO
414	END DO
415	END DO
416	END DO
417	ELSE
418	!> jc tmp
419	DO jn = 1,jpts
420	tabres(i1:i2,j1:j2,k1:k2,jn) = tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) &
421	& * tmask(i1:i2,j1:j2,k1:k2)
422	ENDDO
423	!< jc tmp
424	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
425	! Add asselin part
426	DO jn = 1,jpts
427	DO jk = k1, k2
428	DO jj = j1, j2
429	DO ji = i1, i2
430	IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN
431	ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used
432	ztnu = tabres(ji,jj,jk,jn)
433	ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a)
434	ts(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) &
435	& * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a)
436	ENDIF
437	END DO
438	END DO
439	END DO
440	END DO
441	ENDIF
442	DO jn = 1,jpts
443	DO jk=k1,k2
444	DO jj=j1,j2
445	DO ji=i1,i2
446	IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN
447	ts(ji,jj,jk,jn,Kmm_a) = tabres(ji,jj,jk,jn) / e3t(ji,jj,jk,Kmm_a)
448	END IF
449	END DO
450	END DO
451	END DO
452	END DO
453	!
454	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
455	ts(i1:i2,j1:j2,k1:k2,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,k1:k2,1:jpts,Kmm_a)
456	ENDIF
457	!
458	ENDIF
459	!
460	END SUBROUTINE updateTS
461
462	# endif
463
464	# if defined key_vertical
465
466	SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
467	!!---------------------------------------------
468	!! * ROUTINE updateu *
469	!!---------------------------------------------
470	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
471	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
472	LOGICAL , INTENT(in ) :: before
473	!
474	INTEGER :: ji, jj, jk
475	REAL(wp):: zrhoy, zub, zunu, zuno
476	! VERTICAL REFINEMENT BEGIN
477	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: tabres_child
478	REAL(wp) :: h_in(k1:k2)
479	REAL(wp) :: h_out(1:jpk)
480	INTEGER :: N_in, N_out
481	REAL(wp) :: h_diff, excess, thick
482	REAL(wp) :: tabin(k1:k2)
483	! VERTICAL REFINEMENT END
484	!!---------------------------------------------
485	!
486	IF( before ) THEN
487	zrhoy = Agrif_Rhoy()
488	!jc_alt
489	! AGRIF_SpecialValue = -999._wp
490	DO jk=k1,k2
491	DO jj=j1,j2
492	DO ji=i1,i2
493	!jc_alt
494	! tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * umask(ji,jj,jk) * uu(ji,jj,jk,Kmm_a) &
495	! & + (umask(ji,jj,jk)-1._wp)*999._wp
496	tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * umask(ji,jj,jk) * uu(ji,jj,jk,Kmm_a)
497	!jc_alt
498	! tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) &
499	! & + (umask(ji,jj,jk)-1._wp)*999._wp
500	tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a)
501	END DO
502	END DO
503	END DO
504	ELSE
505	tabres_child(:,:,:) = 0.
506	AGRIF_SpecialValue = 0._wp
507	DO jj=j1,j2
508	DO ji=i1,i2
509	N_in = 0
510	h_in(:) = 0._wp
511	tabin(:) = 0._wp
512	DO jk=k1,k2 !k2=jpk of child grid
513	!jc_alt
514	! IF( tabres(ji,jj,jk,2) < -900._wp) EXIT
515	IF( tabres(ji,jj,jk,2) == 0.) EXIT
516	N_in = N_in + 1
517	tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)
518	h_in(N_in) = tabres(ji,jj,jk,2)/e2u(ji,jj)
519	ENDDO
520	N_out = 0
521	DO jk=1,jpk
522	IF (umask(ji,jj,jk) == 0) EXIT
523	N_out = N_out + 1
524	h_out(N_out) = e3u(ji,jj,jk,Kmm_a)
525	ENDDO
526	IF (N_in * N_out > 0) THEN
527	h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
528	excess = 0._wp
529	IF (h_diff < -1.e-4) THEN
530	!Even if bathy at T points match it's possible for the U points to be deeper in the child grid.
531	!In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.
532	DO jk=N_in,1,-1
533	thick = MIN(-1*h_diff, h_in(jk))
534	excess = excess + tabin(jk)thicke2u(ji,jj)
535	tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))
536	h_diff = h_diff + thick
537	IF ( h_diff == 0) THEN
538	N_in = jk
539	h_in(jk) = h_in(jk) - thick
540	EXIT
541	ENDIF
542	ENDDO
543	ENDIF
544	CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1)
545	tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e2u(ji,jj)*h_out(N_out))
546	ENDIF
547	ENDDO
548	ENDDO
549	!
550	DO jk=1,jpk
551	DO jj=j1,j2
552	DO ji=i1,i2
553	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
554	zub = uu(ji,jj,jk,Kbb_a) * e3u(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used
555	zuno = uu(ji,jj,jk,Kmm_a) * e3u(ji,jj,jk,Krhs_a)
556	zunu = tabres_child(ji,jj,jk) * e3u(ji,jj,jk,Kmm_a)
557	uu(ji,jj,jk,Kbb_a) = ( zub + atfp * ( zunu - zuno) ) &
558	& * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a)
559	ENDIF
560	!
561	uu(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)
562	END DO
563	END DO
564	END DO
565	!
566	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
567	uu(i1:i2,j1:j2,1:jpkm1,Kbb_a) = uu(i1:i2,j1:j2,1:jpkm1,Kmm_a)
568	ENDIF
569	!
570	ENDIF
571	!
572	END SUBROUTINE updateu
573
574	#else
575
576	SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
577	!!---------------------------------------------
578	!! * ROUTINE updateu *
579	!!---------------------------------------------
580	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
581	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
582	LOGICAL , INTENT(in ) :: before
583	!
584	INTEGER :: ji, jj, jk
585	REAL(wp) :: zrhoy, zub, zunu, zuno
586	!!---------------------------------------------
587	!
588	IF( before ) THEN
589	zrhoy = Agrif_Rhoy()
590	DO jk = k1, k2
591	tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u(i1:i2,j1:j2,jk,Kmm_a) * uu(i1:i2,j1:j2,jk,Kmm_a)
592	END DO
593	ELSE
594	DO jk=k1,k2
595	DO jj=j1,j2
596	DO ji=i1,i2
597	tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e2u(ji,jj)
598	!
599	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
600	zub = uu(ji,jj,jk,Kbb_a) * e3u(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used
601	zuno = uu(ji,jj,jk,Kmm_a) * e3u(ji,jj,jk,Krhs_a)
602	zunu = tabres(ji,jj,jk,1)
603	uu(ji,jj,jk,Kbb_a) = ( zub + atfp * ( zunu - zuno) ) &
604	& * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a)
605	ENDIF
606	!
607	uu(ji,jj,jk,Kmm_a) = tabres(ji,jj,jk,1) * umask(ji,jj,jk) / e3u(ji,jj,jk,Kmm_a)
608	END DO
609	END DO
610	END DO
611	!
612	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
613	uu(i1:i2,j1:j2,k1:k2,Kbb_a) = uu(i1:i2,j1:j2,k1:k2,Kmm_a)
614	ENDIF
615	!
616	ENDIF
617	!
618	END SUBROUTINE updateu
619
620	# endif
621
622	SUBROUTINE correct_u_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
623	!!---------------------------------------------
624	!! * ROUTINE correct_u_bdy *
625	!!---------------------------------------------
626	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
627	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
628	LOGICAL , INTENT(in ) :: before
629	INTEGER , INTENT(in) :: nb, ndir
630	!!
631	LOGICAL :: western_side, eastern_side
632	!
633	INTEGER :: jj, jk
634	REAL(wp) :: zcor
635	!!---------------------------------------------
636	!
637	IF( .NOT.before ) THEN
638	!
639	western_side = (nb == 1).AND.(ndir == 1)
640	eastern_side = (nb == 1).AND.(ndir == 2)
641	!
642	IF (western_side) THEN
643	DO jj=j1,j2
644	zcor = uu_b(i1-1,jj,Kmm_a) * hu(i1-1,jj,Krhs_a) * r1_hu(i1-1,jj,Kmm_a) - uu_b(i1-1,jj,Kmm_a)
645	uu_b(i1-1,jj,Kmm_a) = uu_b(i1-1,jj,Kmm_a) + zcor
646	DO jk=1,jpkm1
647	uu(i1-1,jj,jk,Kmm_a) = uu(i1-1,jj,jk,Kmm_a) + zcor * umask(i1-1,jj,jk)
648	END DO
649	END DO
650	ENDIF
651	!
652	IF (eastern_side) THEN
653	DO jj=j1,j2
654	zcor = uu_b(i2+1,jj,Kmm_a) * hu(i2+1,jj,Krhs_a) * r1_hu(i2+1,jj,Kmm_a) - uu_b(i2+1,jj,Kmm_a)
655	uu_b(i2+1,jj,Kmm_a) = uu_b(i2+1,jj,Kmm_a) + zcor
656	DO jk=1,jpkm1
657	uu(i2+1,jj,jk,Kmm_a) = uu(i2+1,jj,jk,Kmm_a) + zcor * umask(i2+1,jj,jk)
658	END DO
659	END DO
660	ENDIF
661	!
662	ENDIF
663	!
664	END SUBROUTINE correct_u_bdy
665
666	# if defined key_vertical
667
668	SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
669	!!---------------------------------------------
670	!! * ROUTINE updatev *
671	!!---------------------------------------------
672	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
673	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
674	LOGICAL , INTENT(in ) :: before
675	!
676	INTEGER :: ji, jj, jk
677	REAL(wp) :: zrhox, zvb, zvnu, zvno
678	! VERTICAL REFINEMENT BEGIN
679	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: tabres_child
680	REAL(wp) :: h_in(k1:k2)
681	REAL(wp) :: h_out(1:jpk)
682	INTEGER :: N_in, N_out
683	REAL(wp) :: h_diff, excess, thick
684	REAL(wp) :: tabin(k1:k2)
685	! VERTICAL REFINEMENT END
686	!!---------------------------------------------
687	!
688	IF( before ) THEN
689	zrhox = Agrif_Rhox()
690	!jc_alt
691	! AGRIF_SpecialValue = -999._wp
692	DO jk=k1,k2
693	DO jj=j1,j2
694	DO ji=i1,i2
695	!jc_alt
696	! tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) * vv(ji,jj,jk,Kmm_a) &
697	! & + (vmask(ji,jj,jk)-1._wp) * 999._wp
698	tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) * vv(ji,jj,jk,Kmm_a)
699	!jc_alt
700	! tabres(ji,jj,jk,2) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) &
701	! & + (vmask(ji,jj,jk)-1._wp) * 999._wp
702	tabres(ji,jj,jk,2) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk)
703	END DO
704	END DO
705	END DO
706	ELSE
707	tabres_child(:,:,:) = 0.
708	AGRIF_SpecialValue = 0._wp
709	DO jj=j1,j2
710	DO ji=i1,i2
711	N_in = 0
712	DO jk=k1,k2
713	!jc_alt
714	! IF (tabres(ji,jj,jk,2) < -900._wp) EXIT
715	IF (tabres(ji,jj,jk,2) == 0) EXIT
716	N_in = N_in + 1
717	tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)
718	h_in(N_in) = tabres(ji,jj,jk,2)/e1v(ji,jj)
719	ENDDO
720	N_out = 0
721	DO jk=1,jpk
722	IF (vmask(ji,jj,jk) == 0) EXIT
723	N_out = N_out + 1
724	h_out(N_out) = e3v(ji,jj,jk,Kmm_a)
725	ENDDO
726	IF (N_in * N_out > 0) THEN
727	h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
728	excess = 0._wp
729	IF (h_diff < -1.e-4) then
730	!Even if bathy at T points match it's possible for the V points to be deeper in the child grid.
731	!In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.
732	DO jk=N_in,1,-1
733	thick = MIN(-1*h_diff, h_in(jk))
734	excess = excess + tabin(jk)thicke2u(ji,jj)
735	tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))
736	h_diff = h_diff + thick
737	IF ( h_diff == 0) THEN
738	N_in = jk
739	h_in(jk) = h_in(jk) - thick
740	EXIT
741	ENDIF
742	ENDDO
743	ENDIF
744	CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1)
745	tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e1v(ji,jj)*h_out(N_out))
746	ENDIF
747	ENDDO
748	ENDDO
749	!
750	DO jk=1,jpkm1
751	DO jj=j1,j2
752	DO ji=i1,i2
753	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
754	zvb = vv(ji,jj,jk,Kbb_a) * e3v(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used
755	zvno = vv(ji,jj,jk,Kmm_a) * e3v(ji,jj,jk,Krhs_a)
756	zvnu = tabres_child(ji,jj,jk) * e3v(ji,jj,jk,Kmm_a)
757	vv(ji,jj,jk,Kbb_a) = ( zvb + atfp * ( zvnu - zvno) ) &
758	& * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a)
759	ENDIF
760	!
761	vv(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)
762	END DO
763	END DO
764	END DO
765	!
766	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
767	vv(i1:i2,j1:j2,1:jpkm1,Kbb_a) = vv(i1:i2,j1:j2,1:jpkm1,Kmm_a)
768	ENDIF
769	!
770	ENDIF
771	!
772	END SUBROUTINE updatev
773
774	# else
775
776	SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before)
777	!!---------------------------------------------
778	!! * ROUTINE updatev *
779	!!---------------------------------------------
780	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
781	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
782	LOGICAL , INTENT(in ) :: before
783	!
784	INTEGER :: ji, jj, jk
785	REAL(wp) :: zrhox, zvb, zvnu, zvno
786	!!---------------------------------------------
787	!
788	IF (before) THEN
789	zrhox = Agrif_Rhox()
790	DO jk=k1,k2
791	DO jj=j1,j2
792	DO ji=i1,i2
793	tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)
794	END DO
795	END DO
796	END DO
797	ELSE
798	DO jk=k1,k2
799	DO jj=j1,j2
800	DO ji=i1,i2
801	tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e1v(ji,jj)
802	!
803	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
804	zvb = vv(ji,jj,jk,Kbb_a) * e3v(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used
805	zvno = vv(ji,jj,jk,Kmm_a) * e3v(ji,jj,jk,Krhs_a)
806	zvnu = tabres(ji,jj,jk,1)
807	vv(ji,jj,jk,Kbb_a) = ( zvb + atfp * ( zvnu - zvno) ) &
808	& * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a)
809	ENDIF
810	!
811	vv(ji,jj,jk,Kmm_a) = tabres(ji,jj,jk,1) * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kmm_a)
812	END DO
813	END DO
814	END DO
815	!
816	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
817	vv(i1:i2,j1:j2,k1:k2,Kbb_a) = vv(i1:i2,j1:j2,k1:k2,Kmm_a)
818	ENDIF
819	!
820	ENDIF
821	!
822	END SUBROUTINE updatev
823
824	# endif
825
826	SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
827	!!---------------------------------------------
828	!! * ROUTINE correct_u_bdy *
829	!!---------------------------------------------
830	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
831	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
832	LOGICAL , INTENT(in ) :: before
833	INTEGER , INTENT(in) :: nb, ndir
834	!!
835	LOGICAL :: southern_side, northern_side
836	!
837	INTEGER :: ji, jk
838	REAL(wp) :: zcor
839	!!---------------------------------------------
840	!
841	IF( .NOT.before ) THEN
842	!
843	southern_side = (nb == 2).AND.(ndir == 1)
844	northern_side = (nb == 2).AND.(ndir == 2)
845	!
846	IF (southern_side) THEN
847	DO ji=i1,i2
848	zcor = vv_b(ji,j1-1,Kmm_a) * hv(ji,j1-1,Krhs_a) * r1_hv(ji,j1-1,Kmm_a) - vv_b(ji,j1-1,Kmm_a)
849	vv_b(ji,j1-1,Kmm_a) = vv_b(ji,j1-1,Kmm_a) + zcor
850	DO jk=1,jpkm1
851	vv(ji,j1-1,jk,Kmm_a) = vv(ji,j1-1,jk,Kmm_a) + zcor * vmask(ji,j1-1,jk)
852	END DO
853	END DO
854	ENDIF
855	!
856	IF (northern_side) THEN
857	DO ji=i1,i2
858	zcor = vv_b(ji,j2+1,Kmm_a) * hv(ji,j2+1,Krhs_a) * r1_hv(ji,j2+1,Kmm_a) - vv_b(ji,j2+1,Kmm_a)
859	vv_b(ji,j2+1,Kmm_a) = vv_b(ji,j2+1,Kmm_a) + zcor
860	DO jk=1,jpkm1
861	vv(ji,j2+1,jk,Kmm_a) = vv(ji,j2+1,jk,Kmm_a) + zcor * vmask(ji,j2+1,jk)
862	END DO
863	END DO
864	ENDIF
865	!
866	ENDIF
867	!
868	END SUBROUTINE correct_v_bdy
869
870
871	SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
872	!!----------------------------------------------------------------------
873	!! * ROUTINE updateu2d *
874	!!----------------------------------------------------------------------
875	INTEGER , INTENT(in ) :: i1, i2, j1, j2
876	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
877	LOGICAL , INTENT(in ) :: before
878	!!
879	INTEGER :: ji, jj, jk
880	REAL(wp) :: zrhoy
881	REAL(wp) :: zcorr
882	!!---------------------------------------------
883	!
884	IF( before ) THEN
885	zrhoy = Agrif_Rhoy()
886	DO jj=j1,j2
887	DO ji=i1,i2
888	tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm_a) * hu(ji,jj,Kmm_a) * e2u(ji,jj)
889	END DO
890	END DO
891	ELSE
892	DO jj=j1,j2
893	DO ji=i1,i2
894	tabres(ji,jj) = tabres(ji,jj) * r1_e2u(ji,jj)
895	!
896	! Update "now" 3d velocities:
897	spgu(ji,jj) = 0._wp
898	DO jk=1,jpkm1
899	spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a)
900	END DO
901	!
902	zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu(ji,jj,Kmm_a)
903	DO jk=1,jpkm1
904	uu(ji,jj,jk,Kmm_a) = uu(ji,jj,jk,Kmm_a) + zcorr * umask(ji,jj,jk)
905	END DO
906	!
907	! Update barotropic velocities:
908	IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
909	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
910	zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a)
911	uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + atfp * zcorr * umask(ji,jj,1)
912	END IF
913	ENDIF
914	uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1)
915	!
916	! Correct "before" velocities to hold correct bt component:
917	spgu(ji,jj) = 0.e0
918	DO jk=1,jpkm1
919	spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kbb_a) * uu(ji,jj,jk,Kbb_a)
920	END DO
921	!
922	zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu(ji,jj,Kbb_a)
923	DO jk=1,jpkm1
924	uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) + zcorr * umask(ji,jj,jk)
925	END DO
926	!
927	END DO
928	END DO
929	!
930	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
931	uu_b(i1:i2,j1:j2,Kbb_a) = uu_b(i1:i2,j1:j2,Kmm_a)
932	ENDIF
933	ENDIF
934	!
935	END SUBROUTINE updateu2d
936
937
938	SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
939	!!----------------------------------------------------------------------
940	!! * ROUTINE updatev2d *
941	!!----------------------------------------------------------------------
942	INTEGER , INTENT(in ) :: i1, i2, j1, j2
943	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
944	LOGICAL , INTENT(in ) :: before
945	!
946	INTEGER :: ji, jj, jk
947	REAL(wp) :: zrhox, zcorr
948	!!----------------------------------------------------------------------
949	!
950	IF( before ) THEN
951	zrhox = Agrif_Rhox()
952	DO jj=j1,j2
953	DO ji=i1,i2
954	tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm_a) * hv(ji,jj,Kmm_a) * e1v(ji,jj)
955	END DO
956	END DO
957	ELSE
958	DO jj=j1,j2
959	DO ji=i1,i2
960	tabres(ji,jj) = tabres(ji,jj) * r1_e1v(ji,jj)
961	!
962	! Update "now" 3d velocities:
963	spgv(ji,jj) = 0.e0
964	DO jk=1,jpkm1
965	spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)
966	END DO
967	!
968	zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv(ji,jj,Kmm_a)
969	DO jk=1,jpkm1
970	vv(ji,jj,jk,Kmm_a) = vv(ji,jj,jk,Kmm_a) + zcorr * vmask(ji,jj,jk)
971	END DO
972	!
973	! Update barotropic velocities:
974	IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
975	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
976	zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a)
977	vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + atfp * zcorr * vmask(ji,jj,1)
978	END IF
979	ENDIF
980	vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1)
981	!
982	! Correct "before" velocities to hold correct bt component:
983	spgv(ji,jj) = 0.e0
984	DO jk=1,jpkm1
985	spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kbb_a) * vv(ji,jj,jk,Kbb_a)
986	END DO
987	!
988	zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv(ji,jj,Kbb_a)
989	DO jk=1,jpkm1
990	vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) + zcorr * vmask(ji,jj,jk)
991	END DO
992	!
993	END DO
994	END DO
995	!
996	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
997	vv_b(i1:i2,j1:j2,Kbb_a) = vv_b(i1:i2,j1:j2,Kmm_a)
998	ENDIF
999	!
1000	ENDIF
1001	!
1002	END SUBROUTINE updatev2d
1003
1004
1005	SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
1006	!!----------------------------------------------------------------------
1007	!! * ROUTINE updateSSH *
1008	!!----------------------------------------------------------------------
1009	INTEGER , INTENT(in ) :: i1, i2, j1, j2
1010	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
1011	LOGICAL , INTENT(in ) :: before
1012	!!
1013	INTEGER :: ji, jj
1014	!!----------------------------------------------------------------------
1015	!
1016	IF( before ) THEN
1017	DO jj=j1,j2
1018	DO ji=i1,i2
1019	tabres(ji,jj) = ssh(ji,jj,Kmm_a)
1020	END DO
1021	END DO
1022	ELSE
1023	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
1024	DO jj=j1,j2
1025	DO ji=i1,i2
1026	ssh(ji,jj,Kbb_a) = ssh(ji,jj,Kbb_a) &
1027	& + atfp * ( tabres(ji,jj) - ssh(ji,jj,Kmm_a) ) * tmask(ji,jj,1)
1028	END DO
1029	END DO
1030	ENDIF
1031	!
1032	DO jj=j1,j2
1033	DO ji=i1,i2
1034	ssh(ji,jj,Kmm_a) = tabres(ji,jj) * tmask(ji,jj,1)
1035	END DO
1036	END DO
1037	!
1038	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
1039	ssh(i1:i2,j1:j2,Kbb_a) = ssh(i1:i2,j1:j2,Kmm_a)
1040	ENDIF
1041	!
1042
1043	ENDIF
1044	!
1045	END SUBROUTINE updateSSH
1046
1047
1048	SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
1049	!!----------------------------------------------------------------------
1050	!! * ROUTINE updateub2b *
1051	!!----------------------------------------------------------------------
1052	INTEGER , INTENT(in) :: i1, i2, j1, j2
1053	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
1054	LOGICAL , INTENT(in) :: before
1055	!!
1056	INTEGER :: ji, jj
1057	REAL(wp) :: zrhoy, za1, zcor
1058	!!---------------------------------------------
1059	!
1060	IF (before) THEN
1061	zrhoy = Agrif_Rhoy()
1062	DO jj=j1,j2
1063	DO ji=i1,i2
1064	tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
1065	END DO
1066	END DO
1067	tabres = zrhoy * tabres
1068	ELSE
1069	!
1070	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
1071	!
1072	za1 = 1._wp / REAL(Agrif_rhot(), wp)
1073	DO jj=j1,j2
1074	DO ji=i1,i2
1075	zcor=tabres(ji,jj) - ub2_b(ji,jj)
1076	! Update time integrated fluxes also in case of multiply nested grids:
1077	ub2_i_b(ji,jj) = ub2_i_b(ji,jj) + za1 * zcor
1078	! Update corrective fluxes:
1079	un_bf(ji,jj) = un_bf(ji,jj) + zcor
1080	! Update half step back fluxes:
1081	ub2_b(ji,jj) = tabres(ji,jj)
1082	END DO
1083	END DO
1084	ENDIF
1085	!
1086	END SUBROUTINE updateub2b
1087
1088	SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
1089	!!---------------------------------------------
1090	!! * ROUTINE reflux_sshu *
1091	!!---------------------------------------------
1092	INTEGER, INTENT(in) :: i1, i2, j1, j2
1093	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
1094	LOGICAL, INTENT(in) :: before
1095	INTEGER, INTENT(in) :: nb, ndir
1096	!!
1097	LOGICAL :: western_side, eastern_side
1098	INTEGER :: ji, jj
1099	REAL(wp) :: zrhoy, za1, zcor
1100	!!---------------------------------------------
1101	!
1102	IF (before) THEN
1103	zrhoy = Agrif_Rhoy()
1104	DO jj=j1,j2
1105	DO ji=i1,i2
1106	tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
1107	END DO
1108	END DO
1109	tabres = zrhoy * tabres
1110	ELSE
1111	!
1112	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
1113	!
1114	western_side = (nb == 1).AND.(ndir == 1)
1115	eastern_side = (nb == 1).AND.(ndir == 2)
1116	!
1117	IF (western_side) THEN
1118	DO jj=j1,j2
1119	zcor = rdt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj))
1120	ssh(i1 ,jj,Kmm_a) = ssh(i1 ,jj,Kmm_a) + zcor
1121	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(i1 ,jj,Kbb_a) = ssh(i1 ,jj,Kbb_a) + atfp * zcor
1122	END DO
1123	ENDIF
1124	IF (eastern_side) THEN
1125	DO jj=j1,j2
1126	zcor = - rdt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj))
1127	ssh(i2+1,jj,Kmm_a) = ssh(i2+1,jj,Kmm_a) + zcor
1128	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(i2+1,jj,Kbb_a) = ssh(i2+1,jj,Kbb_a) + atfp * zcor
1129	END DO
1130	ENDIF
1131	!
1132	ENDIF
1133	!
1134	END SUBROUTINE reflux_sshu
1135
1136	SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
1137	!!----------------------------------------------------------------------
1138	!! * ROUTINE updatevb2b *
1139	!!----------------------------------------------------------------------
1140	INTEGER , INTENT(in ) :: i1, i2, j1, j2
1141	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
1142	LOGICAL , INTENT(in ) :: before
1143	!!
1144	INTEGER :: ji, jj
1145	REAL(wp) :: zrhox, za1, zcor
1146	!!---------------------------------------------
1147	!
1148	IF( before ) THEN
1149	zrhox = Agrif_Rhox()
1150	DO jj=j1,j2
1151	DO ji=i1,i2
1152	tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj)
1153	END DO
1154	END DO
1155	tabres = zrhox * tabres
1156	ELSE
1157	!
1158	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
1159	!
1160	za1 = 1._wp / REAL(Agrif_rhot(), wp)
1161	DO jj=j1,j2
1162	DO ji=i1,i2
1163	zcor=tabres(ji,jj) - vb2_b(ji,jj)
1164	! Update time integrated fluxes also in case of multiply nested grids:
1165	vb2_i_b(ji,jj) = vb2_i_b(ji,jj) + za1 * zcor
1166	! Update corrective fluxes:
1167	vn_bf(ji,jj) = vn_bf(ji,jj) + zcor
1168	! Update half step back fluxes:
1169	vb2_b(ji,jj) = tabres(ji,jj)
1170	END DO
1171	END DO
1172	ENDIF
1173	!
1174	END SUBROUTINE updatevb2b
1175
1176	SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
1177	!!---------------------------------------------
1178	!! * ROUTINE reflux_sshv *
1179	!!---------------------------------------------
1180	INTEGER, INTENT(in) :: i1, i2, j1, j2
1181	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
1182	LOGICAL, INTENT(in) :: before
1183	INTEGER, INTENT(in) :: nb, ndir
1184	!!
1185	LOGICAL :: southern_side, northern_side
1186	INTEGER :: ji, jj
1187	REAL(wp) :: zrhox, za1, zcor
1188	!!---------------------------------------------
1189	!
1190	IF (before) THEN
1191	zrhox = Agrif_Rhox()
1192	DO jj=j1,j2
1193	DO ji=i1,i2
1194	tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj)
1195	END DO
1196	END DO
1197	tabres = zrhox * tabres
1198	ELSE
1199	!
1200	tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
1201	!
1202	southern_side = (nb == 2).AND.(ndir == 1)
1203	northern_side = (nb == 2).AND.(ndir == 2)
1204	!
1205	IF (southern_side) THEN
1206	DO ji=i1,i2
1207	zcor = rdt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1))
1208	ssh(ji,j1 ,Kmm_a) = ssh(ji,j1 ,Kmm_a) + zcor
1209	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(ji,j1 ,Kbb_a) = ssh(ji,j1,Kbb_a) + atfp * zcor
1210	END DO
1211	ENDIF
1212	IF (northern_side) THEN
1213	DO ji=i1,i2
1214	zcor = - rdt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2))
1215	ssh(ji,j2+1,Kmm_a) = ssh(ji,j2+1,Kmm_a) + zcor
1216	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(ji,j2+1,Kbb_a) = ssh(ji,j2+1,Kbb_a) + atfp * zcor
1217	END DO
1218	ENDIF
1219	!
1220	ENDIF
1221	!
1222	END SUBROUTINE reflux_sshv
1223
1224	SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
1225	!
1226	! ====>>>>>>>>>> currently not used
1227	!
1228	!!----------------------------------------------------------------------
1229	!! * ROUTINE updateT *
1230	!!----------------------------------------------------------------------
1231	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2
1232	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
1233	LOGICAL , INTENT(in ) :: before
1234	!!
1235	INTEGER :: ji,jj,jk
1236	REAL(wp) :: ztemp
1237	!!----------------------------------------------------------------------
1238
1239	IF (before) THEN
1240	DO jk=k1,k2
1241	DO jj=j1,j2
1242	DO ji=i1,i2
1243	tabres(ji,jj,jk,1) = e1t(ji,jj)e2t(ji,jj)tmask(ji,jj,jk)
1244	tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
1245	tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
1246	END DO
1247	END DO
1248	END DO
1249	tabres(:,:,:,1)=tabres(:,:,:,1)Agrif_Rhox()Agrif_Rhoy()
1250	tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
1251	tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
1252	ELSE
1253	DO jk=k1,k2
1254	DO jj=j1,j2
1255	DO ji=i1,i2
1256	IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
1257	print ,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)e2t(ji,jj)*tmask(ji,jj,jk)
1258	print ,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)tmask(ji,jj,jk)
1259	print ,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)tmask(ji,jj,jk)
1260	ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
1261	print *,'CORR = ',ztemp-1.
1262	print ,'NEW VALS = ',tabres(ji,jj,jk,2)ztemp,tabres(ji,jj,jk,3)*ztemp, &
1263	tabres(ji,jj,jk,2)ztemptabres(ji,jj,jk,3)*ztemp
1264	e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
1265	e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
1266	END IF
1267	END DO
1268	END DO
1269	END DO
1270	ENDIF
1271	!
1272	END SUBROUTINE update_scales
1273
1274
1275	SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
1276	!!----------------------------------------------------------------------
1277	!! * ROUTINE updateen *
1278	!!----------------------------------------------------------------------
1279	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
1280	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
1281	LOGICAL , INTENT(in ) :: before
1282	!!----------------------------------------------------------------------
1283	!
1284	IF( before ) THEN
1285	ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
1286	ELSE
1287	en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
1288	ENDIF
1289	!
1290	END SUBROUTINE updateEN
1291
1292
1293	SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
1294	!!----------------------------------------------------------------------
1295	!! * ROUTINE updateavt *
1296	!!----------------------------------------------------------------------
1297	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
1298	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
1299	LOGICAL , INTENT(in ) :: before
1300	!!----------------------------------------------------------------------
1301	!
1302	IF( before ) THEN ; ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
1303	ELSE ; avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
1304	ENDIF
1305	!
1306	END SUBROUTINE updateAVT
1307
1308
1309	SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
1310	!!---------------------------------------------
1311	!! * ROUTINE updateavm *
1312	!!----------------------------------------------------------------------
1313	INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2
1314	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
1315	LOGICAL , INTENT(in ) :: before
1316	!!----------------------------------------------------------------------
1317	!
1318	IF( before ) THEN ; ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
1319	ELSE ; avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
1320	ENDIF
1321	!
1322	END SUBROUTINE updateAVM
1323
1324	SUBROUTINE updatee3t(ptab_dum, i1, i2, j1, j2, k1, k2, before )
1325	!!---------------------------------------------
1326	!! * ROUTINE updatee3t *
1327	!!---------------------------------------------
1328	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ptab_dum
1329	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
1330	LOGICAL, INTENT(in) :: before
1331	!
1332	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptab
1333	INTEGER :: ji,jj,jk
1334	REAL(wp) :: zcoef
1335	!!---------------------------------------------
1336	!
1337	IF (.NOT.before) THEN
1338	!
1339	ALLOCATE(ptab(i1:i2,j1:j2,1:jpk))
1340	!
1341	! Update e3t from ssh (z* case only)
1342	DO jk = 1, jpkm1
1343	DO jj=j1,j2
1344	DO ji=i1,i2
1345	ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + ssh(ji,jj,Kmm_a) &
1346	& *ssmask(ji,jj)/(ht_0(ji,jj)-1._wp + ssmask(ji,jj)))
1347	END DO
1348	END DO
1349	END DO
1350	!
1351	! 1) Updates at BEFORE time step:
1352	! -------------------------------
1353	!
1354	! Save "old" scale factor (prior update) for subsequent asselin correction
1355	! of prognostic variables
1356	e3t(i1:i2,j1:j2,1:jpkm1,Krhs_a) = e3t(i1:i2,j1:j2,1:jpkm1,Kmm_a)
1357
1358	! One should also save e3t(:,:,:,Kbb_a), but lacking of workspace...
1359	! hdiv(i1:i2,j1:j2,1:jpkm1) = e3t(i1:i2,j1:j2,1:jpkm1,Kbb_a)
1360
1361	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN
1362	DO jk = 1, jpkm1
1363	DO jj=j1,j2
1364	DO ji=i1,i2
1365	e3t(ji,jj,jk,Kbb_a) = e3t(ji,jj,jk,Kbb_a) &
1366	& + atfp * ( ptab(ji,jj,jk) - e3t(ji,jj,jk,Kmm_a) )
1367	END DO
1368	END DO
1369	END DO
1370	!
1371	e3w (i1:i2,j1:j2,1,Kbb_a) = e3w_0(i1:i2,j1:j2,1) + e3t(i1:i2,j1:j2,1,Kbb_a) - e3t_0(i1:i2,j1:j2,1)
1372	gdepw(i1:i2,j1:j2,1,Kbb_a) = 0.0_wp
1373	gdept(i1:i2,j1:j2,1,Kbb_a) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kbb_a)
1374	!
1375	DO jk = 2, jpk
1376	DO jj = j1,j2
1377	DO ji = i1,i2
1378	zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
1379	e3w(ji,jj,jk,Kbb_a) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * &
1380	& ( e3t(ji,jj,jk-1,Kbb_a) - e3t_0(ji,jj,jk-1) ) &
1381	& + 0.5_wp * tmask(ji,jj,jk) * &
1382	& ( e3t(ji,jj,jk ,Kbb_a) - e3t_0(ji,jj,jk ) )
1383	gdepw(ji,jj,jk,Kbb_a) = gdepw(ji,jj,jk-1,Kbb_a) + e3t(ji,jj,jk-1,Kbb_a)
1384	gdept(ji,jj,jk,Kbb_a) = zcoef * ( gdepw(ji,jj,jk ,Kbb_a) + 0.5 * e3w(ji,jj,jk,Kbb_a)) &
1385	& + (1-zcoef) * ( gdept(ji,jj,jk-1,Kbb_a) + e3w(ji,jj,jk,Kbb_a))
1386	END DO
1387	END DO
1388	END DO
1389	!
1390	ENDIF
1391	!
1392	! 2) Updates at NOW time step:
1393	! ----------------------------
1394	!
1395	! Update vertical scale factor at T-points:
1396	e3t(i1:i2,j1:j2,1:jpkm1,Kmm_a) = ptab(i1:i2,j1:j2,1:jpkm1)
1397	!
1398	! Update total depth:
1399	ht(i1:i2,j1:j2) = 0._wp
1400	DO jk = 1, jpkm1
1401	ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + e3t(i1:i2,j1:j2,jk,Kmm_a) * tmask(i1:i2,j1:j2,jk)
1402	END DO
1403	!
1404	! Update vertical scale factor at W-points and depths:
1405	e3w (i1:i2,j1:j2,1,Kmm_a) = e3w_0(i1:i2,j1:j2,1) + e3t(i1:i2,j1:j2,1,Kmm_a) - e3t_0(i1:i2,j1:j2,1)
1406	gdept(i1:i2,j1:j2,1,Kmm_a) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kmm_a)
1407	gdepw(i1:i2,j1:j2,1,Kmm_a) = 0.0_wp
1408	gde3w(i1:i2,j1:j2,1) = gdept(i1:i2,j1:j2,1,Kmm_a) - (ht(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh
1409	!
1410	DO jk = 2, jpk
1411	DO jj = j1,j2
1412	DO ji = i1,i2
1413	zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
1414	e3w(ji,jj,jk,Kmm_a) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( e3t(ji,jj,jk-1,Kmm_a) - e3t_0(ji,jj,jk-1) ) &
1415	& + 0.5_wp * tmask(ji,jj,jk) * ( e3t(ji,jj,jk ,Kmm_a) - e3t_0(ji,jj,jk ) )
1416	gdepw(ji,jj,jk,Kmm_a) = gdepw(ji,jj,jk-1,Kmm_a) + e3t(ji,jj,jk-1,Kmm_a)
1417	gdept(ji,jj,jk,Kmm_a) = zcoef * ( gdepw(ji,jj,jk ,Kmm_a) + 0.5 * e3w(ji,jj,jk,Kmm_a)) &
1418	& + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm_a) + e3w(ji,jj,jk,Kmm_a))
1419	gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm_a) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh
1420	END DO
1421	END DO
1422	END DO
1423	!
1424	IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
1425	e3t (i1:i2,j1:j2,1:jpk,Kbb_a) = e3t (i1:i2,j1:j2,1:jpk,Kmm_a)
1426	e3w (i1:i2,j1:j2,1:jpk,Kbb_a) = e3w (i1:i2,j1:j2,1:jpk,Kmm_a)
1427	gdepw(i1:i2,j1:j2,1:jpk,Kbb_a) = gdepw(i1:i2,j1:j2,1:jpk,Kmm_a)
1428	gdept(i1:i2,j1:j2,1:jpk,Kbb_a) = gdept(i1:i2,j1:j2,1:jpk,Kmm_a)
1429	ENDIF
1430	!
1431	DEALLOCATE(ptab)
1432	ENDIF
1433	!
1434	END SUBROUTINE updatee3t
1435
1436	#else
1437	!!----------------------------------------------------------------------
1438	!! Empty module no AGRIF zoom
1439	!!----------------------------------------------------------------------
1440	CONTAINS
1441	SUBROUTINE agrif_oce_update_empty
1442	WRITE(,) 'agrif_oce_update : You should not have seen this print! error?'
1443	END SUBROUTINE agrif_oce_update_empty
1444	#endif
1445
1446	!!======================================================================
1447	END MODULE agrif_oce_update
1448

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

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

Context Navigation

source: NEMO/trunk/src/NST/agrif_oce_update.F90 @ 12377

Download in other formats: