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_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/NST_SRC – NEMO

source: branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90 @ 7971

Last change on this file since 7971 was 7971, checked in by jchanut, 7 years ago

Add zstar coordinate with AGRIF
  • Property svn:keywords set to Id
File size: 31.6 KB
Line 
1#define TWO_WAY        /* TWO WAY NESTING */
2#undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/
3 
4MODULE agrif_opa_update
5#if defined key_agrif  && ! defined key_offline
6   USE par_oce
7   USE oce
8   USE dom_oce
9   USE agrif_oce
10   USE in_out_manager  ! I/O manager
11   USE lib_mpp
12   USE wrk_nemo 
13   USE dynspg_oce
14   USE zdf_oce        ! vertical physics: ocean variables
15   USE domvvl         ! Need interpolation routines
16
17   IMPLICIT NONE
18   PRIVATE
19
20   PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Update_Scales, Agrif_Update_vvl
21
22# if defined key_zdftke
23   PUBLIC Agrif_Update_Tke
24# endif
25   !!----------------------------------------------------------------------
26   !! NEMO/NST 3.6 , NEMO Consortium (2010)
27   !! $Id$
28   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
29   !!----------------------------------------------------------------------
30
31CONTAINS
32
33   RECURSIVE SUBROUTINE Agrif_Update_Tra( )
34      !!---------------------------------------------
35      !!   *** ROUTINE Agrif_Update_Tra ***
36      !!---------------------------------------------
37      !
38      IF (Agrif_Root()) RETURN
39      !
40#if defined TWO_WAY 
41      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers  from grid Number',Agrif_Fixed(), 'nbcline', nbcline
42
43      Agrif_UseSpecialValueInUpdate = .TRUE.
44      Agrif_SpecialValueFineGrid = 0.
45      !
46      IF (MOD(nbcline,nbclineupdate) == 0) THEN
47# if ! defined DECAL_FEEDBACK
48         CALL Agrif_Update_Variable(tsn_id, procname=updateTS)
49# else
50         CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS)
51# endif
52      ELSE
53# if ! defined DECAL_FEEDBACK
54         CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS)
55# else
56         CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS)
57# endif
58      ENDIF
59      !
60      Agrif_UseSpecialValueInUpdate = .FALSE.
61      !
62      IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update:
63         CALL Agrif_ChildGrid_To_ParentGrid()
64         CALL Agrif_Update_Tra()
65         CALL Agrif_ParentGrid_To_ChildGrid()
66      ENDIF
67      !
68#endif
69      !
70   END SUBROUTINE Agrif_Update_Tra
71
72   RECURSIVE SUBROUTINE Agrif_Update_Dyn( )
73      !!---------------------------------------------
74      !!   *** ROUTINE Agrif_Update_Dyn ***
75      !!---------------------------------------------
76      !
77      IF (Agrif_Root()) RETURN
78      !
79#if defined TWO_WAY
80      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed(), 'nbcline', nbcline
81
82      Agrif_UseSpecialValueInUpdate = .FALSE.
83      Agrif_SpecialValueFineGrid = 0.
84      !     
85      IF (mod(nbcline,nbclineupdate) == 0) THEN
86# if ! defined DECAL_FEEDBACK
87         CALL Agrif_Update_Variable(un_update_id,procname = updateU)
88         CALL Agrif_Update_Variable(vn_update_id,procname = updateV)
89# else
90         CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU)
91         CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV)
92# endif
93      ELSE
94# if ! defined DECAL_FEEDBACK
95         CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU)
96         CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)         
97# else
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      ENDIF
102
103# if ! defined DECAL_FEEDBACK
104      CALL Agrif_Update_Variable(e1u_id,procname = updateU2d)
105      CALL Agrif_Update_Variable(e2v_id,procname = updateV2d) 
106# else
107      CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d)
108      CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d) 
109# endif
110
111# if defined key_dynspg_ts
112      IF (ln_bt_fw) THEN
113         ! Update time integrated transports
114         IF (mod(nbcline,nbclineupdate) == 0) THEN
115#  if ! defined DECAL_FEEDBACK
116            CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b)
117            CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b)
118#  else
119            CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b)
120            CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b)
121#  endif
122         ELSE
123#  if ! defined DECAL_FEEDBACK
124            CALL Agrif_Update_Variable(ub2b_update_id,locupdate=(/0,1/),procname = updateub2b)
125            CALL Agrif_Update_Variable(vb2b_update_id,locupdate=(/0,1/),procname = updatevb2b)
126#  else
127            CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateub2b)
128            CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updatevb2b)
129#  endif
130         ENDIF
131      END IF
132# endif
133      !
134      nbcline = nbcline + 1
135      !
136      Agrif_UseSpecialValueInUpdate = .TRUE.
137      Agrif_SpecialValueFineGrid = 0.
138# if ! defined DECAL_FEEDBACK
139      CALL Agrif_Update_Variable(sshn_id,procname = updateSSH)
140# else
141      CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH)
142# endif
143      Agrif_UseSpecialValueInUpdate = .FALSE.
144      !
145#endif
146      !
147      ! Do recursive update:
148      IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update:
149         CALL Agrif_ChildGrid_To_ParentGrid()
150         CALL Agrif_Update_Dyn()
151         CALL Agrif_ParentGrid_To_ChildGrid()
152      ENDIF
153      !
154   END SUBROUTINE Agrif_Update_Dyn
155
156# if defined key_zdftke
157   SUBROUTINE Agrif_Update_Tke( kt )
158      !!---------------------------------------------
159      !!   *** ROUTINE Agrif_Update_Tke ***
160      !!---------------------------------------------
161      !!
162      INTEGER, INTENT(in) :: kt
163      !       
164      IF( ( Agrif_NbStepint() .NE. 0 ) .AND. (kt /= 0) ) RETURN
165#  if defined TWO_WAY
166
167      Agrif_UseSpecialValueInUpdate = .TRUE.
168      Agrif_SpecialValueFineGrid = 0.
169
170      CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN  )
171      CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT )
172      CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM )
173
174      Agrif_UseSpecialValueInUpdate = .FALSE.
175
176#  endif
177     
178   END SUBROUTINE Agrif_Update_Tke
179# endif /* key_zdftke */
180
181   RECURSIVE SUBROUTINE Agrif_Update_vvl( )
182      !!---------------------------------------------
183      !!   *** ROUTINE Agrif_Update_vvl ***
184      !!---------------------------------------------
185      !
186      IF (Agrif_Root()) RETURN
187      !
188#if defined TWO_WAY 
189      !
190      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step()
191      !
192      Agrif_UseSpecialValueInUpdate = .FALSE.
193      Agrif_SpecialValueFineGrid = 0.
194      !
195# if ! defined DECAL_FEEDBACK
196      CALL Agrif_Update_Variable(e3t_id, procname=updatee3t)
197# else
198      CALL Agrif_Update_Variable(e3t_id, locupdate=(/1,0/), procname=updatee3t)
199# endif 
200      !
201      Agrif_UseSpecialValueInUpdate = .FALSE.
202      !
203      CALL Agrif_ChildGrid_To_ParentGrid()
204      CALL dom_vvl_update_UVF
205      CALL Agrif_ParentGrid_To_ChildGrid()
206      !
207      IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update:
208         CALL Agrif_ChildGrid_To_ParentGrid()
209         CALL Agrif_Update_vvl()
210         CALL Agrif_ParentGrid_To_ChildGrid()
211      ENDIF
212      !
213#endif
214      !
215   END SUBROUTINE Agrif_Update_vvl
216
217   SUBROUTINE dom_vvl_update_UVF
218      !!---------------------------------------------
219      !!       *** ROUTINE dom_vvl_update_UVF ***
220      !!---------------------------------------------
221#  include "domzgr_substitute.h90"
222      !!
223      INTEGER :: jk
224      REAL(wp):: zcoef
225      !!---------------------------------------------
226
227      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Finalize e3 on grid Number', &
228                  & Agrif_Fixed(), 'Step', Agrif_Nb_Step()
229
230      CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3u_n(:,:,:),  'U'  )
231      CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3v_n(:,:,:),  'V'  )
232      CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3f_n(:,:,:),  'F'  )
233
234      CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3uw_n(:,:,:), 'UW' )
235      CALL dom_vvl_interpol( fse3v_n(:,:,:), fse3vw_n(:,:,:), 'VW' )
236
237      ! Update total depths:
238      hu(:,:) = 0._wp                        ! Ocean depth at U-points
239      hv(:,:) = 0._wp                        ! Ocean depth at V-points
240      DO jk = 1, jpkm1
241         hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk)
242         hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk)
243      END DO
244      !                                      ! Inverse of the local depth
245      hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:)
246      hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:)
247
248#if ! defined key_dynspg_ts
249      IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
250#else
251      IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
252#endif
253         CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3u_b(:,:,:),  'U'  )
254         CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3v_b(:,:,:),  'V'  )
255
256         CALL dom_vvl_interpol( fse3u_b(:,:,:), fse3uw_b(:,:,:), 'UW' )
257         CALL dom_vvl_interpol( fse3v_b(:,:,:), fse3vw_b(:,:,:), 'VW' )
258
259         hu_b(:,:) = 0._wp                        ! Ocean depth at U-points
260         hv_b(:,:) = 0._wp                        ! Ocean depth at V-points
261         DO jk = 1, jpkm1
262            hu_b(:,:) = hu_b(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk)
263            hv_b(:,:) = hv_b(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk)
264         END DO
265         !                                      ! Inverse of the local depth
266         hur_b(:,:) = 1._wp / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:)
267         hvr_b(:,:) = 1._wp / ( hv_b(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:)
268      ENDIF
269
270      !
271   END SUBROUTINE dom_vvl_update_UVF
272
273   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
274      !!---------------------------------------------
275      !!           *** ROUTINE updateT ***
276      !!---------------------------------------------
277      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
278      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
279      LOGICAL, INTENT(in) :: before
280      !!
281      INTEGER :: ji,jj,jk,jn
282      !!---------------------------------------------
283      !
284      !
285      IF (before) THEN
286         DO jn = n1,n2
287            DO jk=k1,k2
288               DO jj=j1,j2
289                  DO ji=i1,i2
290!> jc tmp
291                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * fse3t_n(ji,jj,jk) / e3t_0(ji,jj,jk)
292!                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * fse3t_n(ji,jj,jk)
293!< jc tmp
294                  END DO
295               END DO
296            END DO
297         END DO
298      ELSE
299!> jc tmp
300         DO jn = n1,n2
301            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2)
302         ENDDO
303!< jc tmp
304
305         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
306            ! Add asselin part
307            DO jn = n1,n2
308               DO jk=k1,k2
309                  DO jj=j1,j2
310                     DO ji=i1,i2
311                        IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
312                           tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn) & 
313                                 & + atfp * ( tabres(ji,jj,jk,jn) &
314                                 &             - tsn(ji,jj,jk,jn)*fse3t_a(ji,jj,jk) ) & 
315                                 & * tmask(ji,jj,jk) / fse3t_b(ji,jj,jk)
316                        ENDIF
317                     ENDDO
318                  ENDDO
319               ENDDO
320            ENDDO
321         ENDIF
322         DO jn = n1,n2
323            DO jk=k1,k2
324               DO jj=j1,j2
325                  DO ji=i1,i2
326                     IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN 
327#if defined key_vvl
328                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / fse3t_n(ji,jj,jk)
329#else
330                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn)
331#endif
332                     END IF
333                  END DO
334               END DO
335            END DO
336         END DO
337      ENDIF
338      !
339   END SUBROUTINE updateTS
340
341   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before )
342      !!---------------------------------------------
343      !!           *** ROUTINE updateu ***
344      !!---------------------------------------------
345#  include "domzgr_substitute.h90"
346      !!
347      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
348      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
349      LOGICAL, INTENT(in) :: before
350      !!
351      INTEGER :: ji, jj, jk
352      REAL(wp) :: zrhoy
353      !!---------------------------------------------
354      !
355      IF (before) THEN
356         zrhoy = Agrif_Rhoy()
357         DO jk=k1,k2
358            DO jj=j1,j2
359               DO ji=i1,i2
360                  tabres(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk)
361                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3u_n(ji,jj,jk)
362               END DO
363            END DO
364         END DO
365         tabres = zrhoy * tabres
366      ELSE
367         DO jk=k1,k2
368            DO jj=j1,j2
369               DO ji=i1,i2
370                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e2u(ji,jj) 
371                  !
372                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
373                     ub(ji,jj,jk) = ub(ji,jj,jk) & 
374                           & + atfp * ( tabres(ji,jj,jk)                     & 
375                           &              - un(ji,jj,jk)*fse3u_a(ji,jj,jk) ) & 
376                           & * umask(ji,jj,jk) / fse3u_b(ji,jj,jk)
377                  ENDIF
378                  !
379                  un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / fse3u_n(ji,jj,jk)
380               END DO
381            END DO
382         END DO
383      ENDIF
384      !
385   END SUBROUTINE updateu
386
387   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before )
388      !!---------------------------------------------
389      !!           *** ROUTINE updatev ***
390      !!---------------------------------------------
391#  include "domzgr_substitute.h90"
392      !!
393      INTEGER :: i1,i2,j1,j2,k1,k2
394      INTEGER :: ji,jj,jk
395      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: tabres
396      LOGICAL :: before
397      !!
398      REAL(wp) :: zrhox
399      !!---------------------------------------------     
400      !
401      IF (before) THEN
402         zrhox = Agrif_Rhox()
403         DO jk=k1,k2
404            DO jj=j1,j2
405               DO ji=i1,i2
406                  tabres(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk)
407                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3v_n(ji,jj,jk)
408               END DO
409            END DO
410         END DO
411         tabres = zrhox * tabres
412      ELSE
413         DO jk=k1,k2
414            DO jj=j1,j2
415               DO ji=i1,i2
416                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e1v(ji,jj)
417                  !
418                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
419                     vb(ji,jj,jk) = vb(ji,jj,jk) & 
420                           & + atfp * ( tabres(ji,jj,jk)                     & 
421                           &              - vn(ji,jj,jk)*fse3v_a(ji,jj,jk) ) & 
422                           & * vmask(ji,jj,jk) / fse3v_b(ji,jj,jk)
423                  ENDIF
424                  !
425                  vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / fse3v_n(ji,jj,jk)
426               END DO
427            END DO
428         END DO
429      ENDIF
430      !
431   END SUBROUTINE updatev
432
433   SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
434      !!---------------------------------------------
435      !!          *** ROUTINE updateu2d ***
436      !!---------------------------------------------
437#  include "domzgr_substitute.h90"
438      !!
439      INTEGER, INTENT(in) :: i1, i2, j1, j2
440      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
441      LOGICAL, INTENT(in) :: before
442      !!
443      INTEGER :: ji, jj, jk
444      REAL(wp) :: zrhoy
445      REAL(wp) :: zcorr
446      !!---------------------------------------------
447      !
448      IF (before) THEN
449         zrhoy = Agrif_Rhoy()
450         DO jj=j1,j2
451            DO ji=i1,i2
452               tabres(ji,jj) = un_b(ji,jj) * hu(ji,jj) * e2u(ji,jj)
453            END DO
454         END DO
455         tabres = zrhoy * tabres
456      ELSE
457         DO jj=j1,j2
458            DO ji=i1,i2
459               tabres(ji,jj) =  tabres(ji,jj) / e2u(ji,jj) 
460               !   
461               ! Update "now" 3d velocities:
462               spgu(ji,jj) = 0.e0
463               DO jk=1,jpkm1
464                  spgu(ji,jj) = spgu(ji,jj) + fse3u_n(ji,jj,jk) * un(ji,jj,jk)
465               END DO
466               !
467               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * hur(ji,jj)
468               DO jk=1,jpkm1             
469                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
470               END DO
471               !
472               ! Update barotropic velocities:
473#if ! defined key_dynspg_ts
474               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
475#else
476               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
477#endif
478                  zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * hur_b(ji,jj)
479                  ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1)
480               END IF             
481               un_b(ji,jj) = tabres(ji,jj) * hur(ji,jj) * umask(ji,jj,1)
482               !       
483               ! Correct "before" velocities to hold correct bt component:
484               spgu(ji,jj) = 0.e0
485               DO jk=1,jpkm1
486                  spgu(ji,jj) = spgu(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk)
487               END DO
488               !
489               zcorr = ub_b(ji,jj) - spgu(ji,jj) * hur_b(ji,jj)
490               DO jk=1,jpkm1             
491                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
492               END DO
493               !
494            END DO
495         END DO
496      ENDIF
497      !
498   END SUBROUTINE updateu2d
499
500   SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
501      !!---------------------------------------------
502      !!          *** ROUTINE updatev2d ***
503      !!---------------------------------------------
504      INTEGER, INTENT(in) :: i1, i2, j1, j2
505      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
506      LOGICAL, INTENT(in) :: before
507      !!
508      INTEGER :: ji, jj, jk
509      REAL(wp) :: zrhox
510      REAL(wp) :: zcorr
511      !!---------------------------------------------
512      !
513      IF (before) THEN
514         zrhox = Agrif_Rhox()
515         DO jj=j1,j2
516            DO ji=i1,i2
517               tabres(ji,jj) = vn_b(ji,jj) * hv(ji,jj) * e1v(ji,jj) 
518            END DO
519         END DO
520         tabres = zrhox * tabres
521      ELSE
522         DO jj=j1,j2
523            DO ji=i1,i2
524               tabres(ji,jj) =  tabres(ji,jj) / e1v(ji,jj) 
525               !   
526               ! Update "now" 3d velocities:
527               spgv(ji,jj) = 0.e0
528               DO jk=1,jpkm1
529                  spgv(ji,jj) = spgv(ji,jj) + fse3v_n(ji,jj,jk) * vn(ji,jj,jk)
530               END DO
531               !
532               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * hvr(ji,jj)
533               DO jk=1,jpkm1             
534                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
535               END DO
536               !
537               ! Update barotropic velocities:
538#if ! defined key_dynspg_ts
539               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
540#else
541               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
542#endif
543                  zcorr = (tabres(ji,jj) - vn_b(ji,jj)*hv_a(ji,jj)) * hvr_b(ji,jj)
544                  vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1)
545               END IF             
546               vn_b(ji,jj) = tabres(ji,jj) * hvr(ji,jj) * vmask(ji,jj,1)
547               !       
548               ! Correct "before" velocities to hold correct bt component:
549               spgv(ji,jj) = 0.e0
550               DO jk=1,jpkm1
551                  spgv(ji,jj) = spgv(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk)
552               END DO
553               !
554               zcorr = vb_b(ji,jj) - spgv(ji,jj) * hvr_b(ji,jj)
555               DO jk=1,jpkm1             
556                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
557               END DO
558               !
559            END DO
560         END DO
561      ENDIF
562      !
563   END SUBROUTINE updatev2d
564
565
566   SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
567      !!---------------------------------------------
568      !!          *** ROUTINE updateSSH ***
569      !!---------------------------------------------
570      INTEGER, INTENT(in) :: i1, i2, j1, j2
571      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
572      LOGICAL, INTENT(in) :: before
573      !!
574      INTEGER :: ji, jj
575      !!---------------------------------------------
576      !
577      IF (before) THEN
578         DO jj=j1,j2
579            DO ji=i1,i2
580               tabres(ji,jj) = sshn(ji,jj)
581            END DO
582         END DO
583      ELSE
584#if ! defined key_dynspg_ts
585         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
586#else
587         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
588#endif
589            DO jj=j1,j2
590               DO ji=i1,i2
591                  sshb(ji,jj) =   sshb(ji,jj) &
592                        & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1)
593               END DO
594            END DO
595         ENDIF
596         !
597         DO jj=j1,j2
598            DO ji=i1,i2
599               sshn(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)
600            END DO
601         END DO
602      ENDIF
603      !
604   END SUBROUTINE updateSSH
605
606   SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
607      !!---------------------------------------------
608      !!          *** ROUTINE updateub2b ***
609      !!---------------------------------------------
610      INTEGER, INTENT(in) :: i1, i2, j1, j2
611      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
612      LOGICAL, INTENT(in) :: before
613      !!
614      INTEGER :: ji, jj
615      REAL(wp) :: zrhoy
616      !!---------------------------------------------
617      !
618      IF (before) THEN
619         zrhoy = Agrif_Rhoy()
620         DO jj=j1,j2
621            DO ji=i1,i2
622               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
623            END DO
624         END DO
625         tabres = zrhoy * tabres
626      ELSE
627         DO jj=j1,j2
628            DO ji=i1,i2
629               ub2_b(ji,jj) = tabres(ji,jj) / e2u(ji,jj)
630            END DO
631         END DO
632      ENDIF
633      !
634   END SUBROUTINE updateub2b
635
636   SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
637      !!---------------------------------------------
638      !!          *** ROUTINE updatevb2b ***
639      !!---------------------------------------------
640      INTEGER, INTENT(in) :: i1, i2, j1, j2
641      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
642      LOGICAL, INTENT(in) :: before
643      !!
644      INTEGER :: ji, jj
645      REAL(wp) :: zrhox
646      !!---------------------------------------------
647      !
648      IF (before) THEN
649         zrhox = Agrif_Rhox()
650         DO jj=j1,j2
651            DO ji=i1,i2
652               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
653            END DO
654         END DO
655         tabres = zrhox * tabres
656      ELSE
657         DO jj=j1,j2
658            DO ji=i1,i2
659               vb2_b(ji,jj) = tabres(ji,jj) / e1v(ji,jj)
660            END DO
661         END DO
662      ENDIF
663      !
664   END SUBROUTINE updatevb2b
665
666
667   SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
668      ! currently not used
669      !!---------------------------------------------
670      !!           *** ROUTINE updateT ***
671      !!---------------------------------------------
672#  include "domzgr_substitute.h90"
673
674      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
675      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
676      LOGICAL, iNTENT(in) :: before
677
678      INTEGER :: ji,jj,jk
679      REAL(wp) :: ztemp
680
681      IF (before) THEN
682         DO jk=k1,k2
683            DO jj=j1,j2
684               DO ji=i1,i2
685                  tabres(ji,jj,jk,1) = e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
686                  tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
687                  tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
688               END DO
689            END DO
690         END DO
691         tabres(:,:,:,1)=tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy()
692         tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
693         tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
694      ELSE
695         DO jk=k1,k2
696            DO jj=j1,j2
697               DO ji=i1,i2
698                  IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
699                     print *,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
700                     print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk)
701                     print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk)
702                     ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
703                     print *,'CORR = ',ztemp-1.
704                     print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, &
705                           tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp
706                     e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
707                     e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
708                  END IF
709               END DO
710            END DO
711         END DO
712      ENDIF
713      !
714   END SUBROUTINE update_scales
715
716# if defined key_zdftke
717   SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
718      !!---------------------------------------------
719      !!           *** ROUTINE updateen ***
720      !!---------------------------------------------
721      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
722      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
723      LOGICAL, INTENT(in) :: before
724      !!---------------------------------------------
725      !
726      IF (before) THEN
727         ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
728      ELSE
729         en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
730      ENDIF
731      !
732   END SUBROUTINE updateEN
733
734
735   SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
736      !!---------------------------------------------
737      !!           *** ROUTINE updateavt ***
738      !!---------------------------------------------
739      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
740      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
741      LOGICAL, INTENT(in) :: before
742      !!---------------------------------------------
743      !
744      IF (before) THEN
745         ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
746      ELSE
747         avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
748      ENDIF
749      !
750   END SUBROUTINE updateAVT
751
752
753   SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
754      !!---------------------------------------------
755      !!           *** ROUTINE updateavm ***
756      !!---------------------------------------------
757      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
758      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
759      LOGICAL, INTENT(in) :: before
760      !!---------------------------------------------
761      !
762      IF (before) THEN
763         ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
764      ELSE
765         avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
766      ENDIF
767      !
768   END SUBROUTINE updateAVM
769
770# endif /* key_zdftke */ 
771
772   SUBROUTINE updatee3t( ptab, i1, i2, j1, j2, k1, k2, before )
773      !!---------------------------------------------
774      !!           *** ROUTINE updatee3t ***
775      !!---------------------------------------------
776      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
777      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
778      LOGICAL, INTENT(in) :: before
779      INTEGER :: ji,jj,jk
780      REAL(wp) :: zcoef
781      !!---------------------------------------------
782      !
783      IF (before) THEN
784!> jc tmp:
785!         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2)
786         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) / e3t_0(i1:i2,j1:j2,k1:k2)
787!< jc tmp:
788      ELSE
789         !
790         ! 1) Updates at before time step:
791         ! -------------------------------
792         !
793!> jc tmp:
794         ptab(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2)
795!< jc tmp:
796
797#if ! defined key_dynspg_ts
798         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
799#else
800         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
801#endif
802            DO jk = 1, jpkm1
803               DO jj=j1,j2
804                  DO ji=i1,i2
805                     fse3t_b(ji,jj,jk) =   fse3t_b(ji,jj,jk) &
806                           & + atfp * ( ptab(ji,jj,jk) - fse3t_n(ji,jj,jk) )
807                  END DO
808               END DO
809            END DO
810            !
811            fse3w_b (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
812            fsdepw_b(i1:i2,j1:j2,1) = 0.0_wp
813            fsdept_b(i1:i2,j1:j2,1) = 0.5_wp * fse3w_b(i1:i2,j1:j2,1)
814            !
815            DO jk = 2, jpkm1
816               DO jj = j1,j2
817                  DO ji = i1,i2           
818                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
819                     fse3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        & 
820                     &                                        ( fse3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  &
821                     &                                    +            0.5_wp * tmask(ji,jj,jk)   *        &
822                     &                                        ( fse3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
823                     fsdepw_b(ji,jj,jk) = fsdepw_b(ji,jj,jk-1) + fse3t_b(ji,jj,jk-1)
824                     fsdept_b(ji,jj,jk) =      zcoef  * ( fsdepw_b(ji,jj,jk  ) + 0.5 * fse3w_b(ji,jj,jk))  &
825                         &                + (1-zcoef) * ( fsdept_b(ji,jj,jk-1) +       fse3w_b(ji,jj,jk)) 
826                  END DO
827               END DO
828            END DO
829         ENDIF       
830         !
831         ! 2) Updates at now time step:
832         ! ----------------------------
833         !
834         ! Update vertical scale factor at T-points:
835         fse3t_n(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
836         !
837         ! Update total depth:
838         ht(i1:i2,j1:j2) = 0._wp
839         DO jk = 1, jpkm1
840            ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + fse3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk)
841         END DO
842         !
843         ! Update vertical scale factor at W-points and depths:
844         fse3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
845         fsdept_n(i1:i2,j1:j2,1) = 0.5_wp * fse3w_n(i1:i2,j1:j2,1)
846         fsdepw_n(i1:i2,j1:j2,1) = 0.0_wp
847         fsde3w_n(i1:i2,j1:j2,1) = fsdept_n(i1:i2,j1:j2,1) - (ht(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh
848         !
849         DO jk = 2, jpkm1
850            DO jj = j1,j2
851               DO ji = i1,i2           
852               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
853               fse3w_n(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( fse3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )   &
854               &                                    +            0.5_wp * tmask(ji,jj,jk)   * ( fse3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
855               fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1)
856               fsdept_n(ji,jj,jk) =      zcoef  * ( fsdepw_n(ji,jj,jk  ) + 0.5 * fse3w_n(ji,jj,jk))  &
857                   &                + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) +       fse3w_n(ji,jj,jk)) 
858               fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh
859               END DO
860            END DO
861         END DO
862         !
863      ENDIF
864      !
865   END SUBROUTINE updatee3t
866
867#else
868CONTAINS
869   SUBROUTINE agrif_opa_update_empty
870      !!---------------------------------------------
871      !!   *** ROUTINE agrif_opa_update_empty ***
872      !!---------------------------------------------
873      WRITE(*,*)  'agrif_opa_update : You should not have seen this print! error?'
874   END SUBROUTINE agrif_opa_update_empty
875#endif
876END MODULE agrif_opa_update
877
Note: See TracBrowser for help on using the repository browser.