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 @ 7977

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

Uncomplete removal of key_vvl
  • Property svn:keywords set to Id
File size: 31.5 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                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / fse3t_n(ji,jj,jk)
328                     END IF
329                  END DO
330               END DO
331            END DO
332         END DO
333      ENDIF
334      !
335   END SUBROUTINE updateTS
336
337   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before )
338      !!---------------------------------------------
339      !!           *** ROUTINE updateu ***
340      !!---------------------------------------------
341#  include "domzgr_substitute.h90"
342      !!
343      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
344      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
345      LOGICAL, INTENT(in) :: before
346      !!
347      INTEGER :: ji, jj, jk
348      REAL(wp) :: zrhoy
349      !!---------------------------------------------
350      !
351      IF (before) THEN
352         zrhoy = Agrif_Rhoy()
353         DO jk=k1,k2
354            DO jj=j1,j2
355               DO ji=i1,i2
356                  tabres(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk)
357                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3u_n(ji,jj,jk)
358               END DO
359            END DO
360         END DO
361         tabres = zrhoy * tabres
362      ELSE
363         DO jk=k1,k2
364            DO jj=j1,j2
365               DO ji=i1,i2
366                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e2u(ji,jj) 
367                  !
368                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
369                     ub(ji,jj,jk) = ub(ji,jj,jk) & 
370                           & + atfp * ( tabres(ji,jj,jk)                     & 
371                           &              - un(ji,jj,jk)*fse3u_a(ji,jj,jk) ) & 
372                           & * umask(ji,jj,jk) / fse3u_b(ji,jj,jk)
373                  ENDIF
374                  !
375                  un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / fse3u_n(ji,jj,jk)
376               END DO
377            END DO
378         END DO
379      ENDIF
380      !
381   END SUBROUTINE updateu
382
383   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before )
384      !!---------------------------------------------
385      !!           *** ROUTINE updatev ***
386      !!---------------------------------------------
387#  include "domzgr_substitute.h90"
388      !!
389      INTEGER :: i1,i2,j1,j2,k1,k2
390      INTEGER :: ji,jj,jk
391      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: tabres
392      LOGICAL :: before
393      !!
394      REAL(wp) :: zrhox
395      !!---------------------------------------------     
396      !
397      IF (before) THEN
398         zrhox = Agrif_Rhox()
399         DO jk=k1,k2
400            DO jj=j1,j2
401               DO ji=i1,i2
402                  tabres(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk)
403                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3v_n(ji,jj,jk)
404               END DO
405            END DO
406         END DO
407         tabres = zrhox * tabres
408      ELSE
409         DO jk=k1,k2
410            DO jj=j1,j2
411               DO ji=i1,i2
412                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e1v(ji,jj)
413                  !
414                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
415                     vb(ji,jj,jk) = vb(ji,jj,jk) & 
416                           & + atfp * ( tabres(ji,jj,jk)                     & 
417                           &              - vn(ji,jj,jk)*fse3v_a(ji,jj,jk) ) & 
418                           & * vmask(ji,jj,jk) / fse3v_b(ji,jj,jk)
419                  ENDIF
420                  !
421                  vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / fse3v_n(ji,jj,jk)
422               END DO
423            END DO
424         END DO
425      ENDIF
426      !
427   END SUBROUTINE updatev
428
429   SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
430      !!---------------------------------------------
431      !!          *** ROUTINE updateu2d ***
432      !!---------------------------------------------
433#  include "domzgr_substitute.h90"
434      !!
435      INTEGER, INTENT(in) :: i1, i2, j1, j2
436      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
437      LOGICAL, INTENT(in) :: before
438      !!
439      INTEGER :: ji, jj, jk
440      REAL(wp) :: zrhoy
441      REAL(wp) :: zcorr
442      !!---------------------------------------------
443      !
444      IF (before) THEN
445         zrhoy = Agrif_Rhoy()
446         DO jj=j1,j2
447            DO ji=i1,i2
448               tabres(ji,jj) = un_b(ji,jj) * hu(ji,jj) * e2u(ji,jj)
449            END DO
450         END DO
451         tabres = zrhoy * tabres
452      ELSE
453         DO jj=j1,j2
454            DO ji=i1,i2
455               tabres(ji,jj) =  tabres(ji,jj) / e2u(ji,jj) 
456               !   
457               ! Update "now" 3d velocities:
458               spgu(ji,jj) = 0.e0
459               DO jk=1,jpkm1
460                  spgu(ji,jj) = spgu(ji,jj) + fse3u_n(ji,jj,jk) * un(ji,jj,jk)
461               END DO
462               !
463               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * hur(ji,jj)
464               DO jk=1,jpkm1             
465                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
466               END DO
467               !
468               ! Update barotropic velocities:
469#if ! defined key_dynspg_ts
470               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
471#else
472               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
473#endif
474                  zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * hur_b(ji,jj)
475                  ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1)
476               END IF             
477               un_b(ji,jj) = tabres(ji,jj) * hur(ji,jj) * umask(ji,jj,1)
478               !       
479               ! Correct "before" velocities to hold correct bt component:
480               spgu(ji,jj) = 0.e0
481               DO jk=1,jpkm1
482                  spgu(ji,jj) = spgu(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk)
483               END DO
484               !
485               zcorr = ub_b(ji,jj) - spgu(ji,jj) * hur_b(ji,jj)
486               DO jk=1,jpkm1             
487                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
488               END DO
489               !
490            END DO
491         END DO
492      ENDIF
493      !
494   END SUBROUTINE updateu2d
495
496   SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
497      !!---------------------------------------------
498      !!          *** ROUTINE updatev2d ***
499      !!---------------------------------------------
500      INTEGER, INTENT(in) :: i1, i2, j1, j2
501      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
502      LOGICAL, INTENT(in) :: before
503      !!
504      INTEGER :: ji, jj, jk
505      REAL(wp) :: zrhox
506      REAL(wp) :: zcorr
507      !!---------------------------------------------
508      !
509      IF (before) THEN
510         zrhox = Agrif_Rhox()
511         DO jj=j1,j2
512            DO ji=i1,i2
513               tabres(ji,jj) = vn_b(ji,jj) * hv(ji,jj) * e1v(ji,jj) 
514            END DO
515         END DO
516         tabres = zrhox * tabres
517      ELSE
518         DO jj=j1,j2
519            DO ji=i1,i2
520               tabres(ji,jj) =  tabres(ji,jj) / e1v(ji,jj) 
521               !   
522               ! Update "now" 3d velocities:
523               spgv(ji,jj) = 0.e0
524               DO jk=1,jpkm1
525                  spgv(ji,jj) = spgv(ji,jj) + fse3v_n(ji,jj,jk) * vn(ji,jj,jk)
526               END DO
527               !
528               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * hvr(ji,jj)
529               DO jk=1,jpkm1             
530                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
531               END DO
532               !
533               ! Update barotropic velocities:
534#if ! defined key_dynspg_ts
535               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
536#else
537               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
538#endif
539                  zcorr = (tabres(ji,jj) - vn_b(ji,jj)*hv_a(ji,jj)) * hvr_b(ji,jj)
540                  vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1)
541               END IF             
542               vn_b(ji,jj) = tabres(ji,jj) * hvr(ji,jj) * vmask(ji,jj,1)
543               !       
544               ! Correct "before" velocities to hold correct bt component:
545               spgv(ji,jj) = 0.e0
546               DO jk=1,jpkm1
547                  spgv(ji,jj) = spgv(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk)
548               END DO
549               !
550               zcorr = vb_b(ji,jj) - spgv(ji,jj) * hvr_b(ji,jj)
551               DO jk=1,jpkm1             
552                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
553               END DO
554               !
555            END DO
556         END DO
557      ENDIF
558      !
559   END SUBROUTINE updatev2d
560
561
562   SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
563      !!---------------------------------------------
564      !!          *** ROUTINE updateSSH ***
565      !!---------------------------------------------
566      INTEGER, INTENT(in) :: i1, i2, j1, j2
567      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
568      LOGICAL, INTENT(in) :: before
569      !!
570      INTEGER :: ji, jj
571      !!---------------------------------------------
572      !
573      IF (before) THEN
574         DO jj=j1,j2
575            DO ji=i1,i2
576               tabres(ji,jj) = sshn(ji,jj)
577            END DO
578         END DO
579      ELSE
580#if ! defined key_dynspg_ts
581         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
582#else
583         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
584#endif
585            DO jj=j1,j2
586               DO ji=i1,i2
587                  sshb(ji,jj) =   sshb(ji,jj) &
588                        & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1)
589               END DO
590            END DO
591         ENDIF
592         !
593         DO jj=j1,j2
594            DO ji=i1,i2
595               sshn(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)
596            END DO
597         END DO
598      ENDIF
599      !
600   END SUBROUTINE updateSSH
601
602   SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
603      !!---------------------------------------------
604      !!          *** ROUTINE updateub2b ***
605      !!---------------------------------------------
606      INTEGER, INTENT(in) :: i1, i2, j1, j2
607      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
608      LOGICAL, INTENT(in) :: before
609      !!
610      INTEGER :: ji, jj
611      REAL(wp) :: zrhoy
612      !!---------------------------------------------
613      !
614      IF (before) THEN
615         zrhoy = Agrif_Rhoy()
616         DO jj=j1,j2
617            DO ji=i1,i2
618               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
619            END DO
620         END DO
621         tabres = zrhoy * tabres
622      ELSE
623         DO jj=j1,j2
624            DO ji=i1,i2
625               ub2_b(ji,jj) = tabres(ji,jj) / e2u(ji,jj)
626            END DO
627         END DO
628      ENDIF
629      !
630   END SUBROUTINE updateub2b
631
632   SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
633      !!---------------------------------------------
634      !!          *** ROUTINE updatevb2b ***
635      !!---------------------------------------------
636      INTEGER, INTENT(in) :: i1, i2, j1, j2
637      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
638      LOGICAL, INTENT(in) :: before
639      !!
640      INTEGER :: ji, jj
641      REAL(wp) :: zrhox
642      !!---------------------------------------------
643      !
644      IF (before) THEN
645         zrhox = Agrif_Rhox()
646         DO jj=j1,j2
647            DO ji=i1,i2
648               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
649            END DO
650         END DO
651         tabres = zrhox * tabres
652      ELSE
653         DO jj=j1,j2
654            DO ji=i1,i2
655               vb2_b(ji,jj) = tabres(ji,jj) / e1v(ji,jj)
656            END DO
657         END DO
658      ENDIF
659      !
660   END SUBROUTINE updatevb2b
661
662
663   SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
664      ! currently not used
665      !!---------------------------------------------
666      !!           *** ROUTINE updateT ***
667      !!---------------------------------------------
668#  include "domzgr_substitute.h90"
669
670      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
671      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
672      LOGICAL, iNTENT(in) :: before
673
674      INTEGER :: ji,jj,jk
675      REAL(wp) :: ztemp
676
677      IF (before) THEN
678         DO jk=k1,k2
679            DO jj=j1,j2
680               DO ji=i1,i2
681                  tabres(ji,jj,jk,1) = e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
682                  tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
683                  tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
684               END DO
685            END DO
686         END DO
687         tabres(:,:,:,1)=tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy()
688         tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
689         tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
690      ELSE
691         DO jk=k1,k2
692            DO jj=j1,j2
693               DO ji=i1,i2
694                  IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
695                     print *,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
696                     print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk)
697                     print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk)
698                     ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
699                     print *,'CORR = ',ztemp-1.
700                     print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, &
701                           tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp
702                     e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
703                     e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
704                  END IF
705               END DO
706            END DO
707         END DO
708      ENDIF
709      !
710   END SUBROUTINE update_scales
711
712# if defined key_zdftke
713   SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
714      !!---------------------------------------------
715      !!           *** ROUTINE updateen ***
716      !!---------------------------------------------
717      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
718      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
719      LOGICAL, INTENT(in) :: before
720      !!---------------------------------------------
721      !
722      IF (before) THEN
723         ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
724      ELSE
725         en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
726      ENDIF
727      !
728   END SUBROUTINE updateEN
729
730
731   SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
732      !!---------------------------------------------
733      !!           *** ROUTINE updateavt ***
734      !!---------------------------------------------
735      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
736      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
737      LOGICAL, INTENT(in) :: before
738      !!---------------------------------------------
739      !
740      IF (before) THEN
741         ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
742      ELSE
743         avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
744      ENDIF
745      !
746   END SUBROUTINE updateAVT
747
748
749   SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
750      !!---------------------------------------------
751      !!           *** ROUTINE updateavm ***
752      !!---------------------------------------------
753      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
754      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
755      LOGICAL, INTENT(in) :: before
756      !!---------------------------------------------
757      !
758      IF (before) THEN
759         ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
760      ELSE
761         avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
762      ENDIF
763      !
764   END SUBROUTINE updateAVM
765
766# endif /* key_zdftke */ 
767
768   SUBROUTINE updatee3t( ptab, i1, i2, j1, j2, k1, k2, before )
769      !!---------------------------------------------
770      !!           *** ROUTINE updatee3t ***
771      !!---------------------------------------------
772      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
773      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
774      LOGICAL, INTENT(in) :: before
775      INTEGER :: ji,jj,jk
776      REAL(wp) :: zcoef
777      !!---------------------------------------------
778      !
779      IF (before) THEN
780!> jc tmp:
781!         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2)
782         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) / e3t_0(i1:i2,j1:j2,k1:k2)
783!< jc tmp:
784      ELSE
785         !
786         ! 1) Updates at before time step:
787         ! -------------------------------
788         !
789!> jc tmp:
790         ptab(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2)
791!< jc tmp:
792
793#if ! defined key_dynspg_ts
794         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
795#else
796         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
797#endif
798            DO jk = 1, jpkm1
799               DO jj=j1,j2
800                  DO ji=i1,i2
801                     fse3t_b(ji,jj,jk) =   fse3t_b(ji,jj,jk) &
802                           & + atfp * ( ptab(ji,jj,jk) - fse3t_n(ji,jj,jk) )
803                  END DO
804               END DO
805            END DO
806            !
807            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)
808            fsdepw_b(i1:i2,j1:j2,1) = 0.0_wp
809            fsdept_b(i1:i2,j1:j2,1) = 0.5_wp * fse3w_b(i1:i2,j1:j2,1)
810            !
811            DO jk = 2, jpkm1
812               DO jj = j1,j2
813                  DO ji = i1,i2           
814                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
815                     fse3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        & 
816                     &                                        ( fse3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  &
817                     &                                    +            0.5_wp * tmask(ji,jj,jk)   *        &
818                     &                                        ( fse3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
819                     fsdepw_b(ji,jj,jk) = fsdepw_b(ji,jj,jk-1) + fse3t_b(ji,jj,jk-1)
820                     fsdept_b(ji,jj,jk) =      zcoef  * ( fsdepw_b(ji,jj,jk  ) + 0.5 * fse3w_b(ji,jj,jk))  &
821                         &                + (1-zcoef) * ( fsdept_b(ji,jj,jk-1) +       fse3w_b(ji,jj,jk)) 
822                  END DO
823               END DO
824            END DO
825         ENDIF       
826         !
827         ! 2) Updates at now time step:
828         ! ----------------------------
829         !
830         ! Update vertical scale factor at T-points:
831         fse3t_n(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
832         !
833         ! Update total depth:
834         ht(i1:i2,j1:j2) = 0._wp
835         DO jk = 1, jpkm1
836            ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + fse3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk)
837         END DO
838         !
839         ! Update vertical scale factor at W-points and depths:
840         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)
841         fsdept_n(i1:i2,j1:j2,1) = 0.5_wp * fse3w_n(i1:i2,j1:j2,1)
842         fsdepw_n(i1:i2,j1:j2,1) = 0.0_wp
843         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
844         !
845         DO jk = 2, jpkm1
846            DO jj = j1,j2
847               DO ji = i1,i2           
848               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
849               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) )   &
850               &                                    +            0.5_wp * tmask(ji,jj,jk)   * ( fse3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
851               fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1)
852               fsdept_n(ji,jj,jk) =      zcoef  * ( fsdepw_n(ji,jj,jk  ) + 0.5 * fse3w_n(ji,jj,jk))  &
853                   &                + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) +       fse3w_n(ji,jj,jk)) 
854               fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh
855               END DO
856            END DO
857         END DO
858         !
859      ENDIF
860      !
861   END SUBROUTINE updatee3t
862
863#else
864CONTAINS
865   SUBROUTINE agrif_opa_update_empty
866      !!---------------------------------------------
867      !!   *** ROUTINE agrif_opa_update_empty ***
868      !!---------------------------------------------
869      WRITE(*,*)  'agrif_opa_update : You should not have seen this print! error?'
870   END SUBROUTINE agrif_opa_update_empty
871#endif
872END MODULE agrif_opa_update
873
Note: See TracBrowser for help on using the repository browser.