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.
dynnxt.F90 in branches/dev_004_VVL/NEMO/OPA_SRC/DYN – NEMO

source: branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynnxt.F90 @ 1361

Last change on this file since 1361 was 1361, checked in by rblod, 15 years ago

dev004_VVL: new organisation see ticket #389

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 8.9 KB
Line 
1MODULE dynnxt
2   !!======================================================================
3   !!                       ***  MODULE  dynnxt  ***
4   !! Ocean dynamics: time stepping
5   !!======================================================================
6   !!======================================================================
7   !! History :  OPA  !  1987-02  (P. Andrich, D. L Hostis)  Original code
8   !!                 !  1990-10  (C. Levy, G. Madec)
9   !!            7.0  !  1993-03  (M. Guyon)  symetrical conditions
10   !!            8.0  !  1997-02  (G. Madec & M. Imbard)  opa, release 8.0
11   !!            8.2  !  1997-04  (A. Weaver)  Euler forward step
12   !!             -   !  1997-06  (G. Madec)  lateral boudary cond., lbc routine
13   !!    NEMO    1.0  !  2002-08  (G. Madec)  F90: Free form and module
14   !!             -   !  2002-10  (C. Talandier, A-M. Treguier) Open boundary cond.
15   !!            2.0  !  2005-11  (V. Garnier) Surface pressure gradient organization
16   !!            2.3  !  2007-07  (D. Storkey) Calls to BDY routines.
17   !!            3.1  !  2009-02  (G. Madec)  re-introduce the vvl option
18   !!----------------------------------------------------------------------
19 
20   !!----------------------------------------------------------------------
21   !!   dyn_nxt      : update the horizontal velocity from the momentum trend
22   !!----------------------------------------------------------------------
23   !! * Modules used
24   USE oce             ! ocean dynamics and tracers
25   USE dom_oce         ! ocean space and time domain
26   USE in_out_manager  ! I/O manager
27   USE obc_oce         ! ocean open boundary conditions
28   USE obcdyn          ! open boundary condition for momentum (obc_dyn routine)
29   USE obcdyn_bt       ! 2D open boundary condition for momentum (obc_dyn_bt routine)
30   USE obcvol          ! ocean open boundary condition (obc_vol routines)
31   USE bdy_oce         ! unstructured open boundary conditions
32   USE bdydta          ! unstructured open boundary conditions
33   USE bdydyn          ! unstructured open boundary conditions
34   USE dynspg_oce      ! type of surface pressure gradient
35   USE lbclnk          ! lateral boundary condition (or mpp link)
36   USE prtctl          ! Print control
37   USE agrif_opa_update
38   USE agrif_opa_interp
39   USE domvvl          ! variable volume
40
41   IMPLICIT NONE
42   PRIVATE
43
44   !! * Accessibility
45   PUBLIC dyn_nxt                ! routine called by step.F90
46   !! * Substitutions
47#  include "domzgr_substitute.h90"
48   !!----------------------------------------------------------------------
49   !!  OPA 9.0 , LOCEAN-IPSL (2005)
50   !! $Id$
51   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
52   !!----------------------------------------------------------------------
53
54CONTAINS
55
56   SUBROUTINE dyn_nxt ( kt )
57      !!----------------------------------------------------------------------
58      !!                  ***  ROUTINE dyn_nxt  ***
59      !!                   
60      !! ** Purpose :   Compute the after horizontal velocity from the
61      !!      momentum trend.
62      !!
63      !! ** Method  :   Apply lateral boundary conditions on the trends (ua,va)
64      !!      through calls to routine lbc_lnk.
65      !!      After velocity is compute using a leap-frog scheme environment:
66      !!         (ua,va) = (ub,vb) + 2 rdt (ua,va)
67      !!      Note that if lk_dynspg_flt=T, the time stepping has already been
68      !!      performed in dynspg module
69      !!      Time filter applied on now horizontal velocity to avoid the
70      !!      divergence of two consecutive time-steps and swap of dynamics
71      !!      arrays to start the next time step:
72      !!         (ub,vb) = (un,vn) + atfp [ (ub,vb) + (ua,va) - 2 (un,vn) ]
73      !!         (un,vn) = (ua,va)
74      !!
75      !! ** Action : - Update ub,vb arrays, the before horizontal velocity
76      !!             - Update un,vn arrays, the now horizontal velocity
77      !!
78      !!----------------------------------------------------------------------
79      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
80      !!
81      INTEGER  ::   jk   ! dummy loop indices
82      REAL(wp) ::   z2dt         ! temporary scalar
83      !!----------------------------------------------------------------------
84
85      IF( kt == nit000 ) THEN
86         IF(lwp) WRITE(numout,*)
87         IF(lwp) WRITE(numout,*) 'dyn_nxt : time stepping'
88         IF(lwp) WRITE(numout,*) '~~~~~~~'
89      ENDIF
90
91      ! Local constant initialization
92      z2dt = 2. * rdt
93      IF( neuler == 0 .AND. kt == nit000 )  z2dt = rdt
94
95      !! Explicit physics with thickness weighted updates
96
97      ! Lateral boundary conditions on ( ua, va )
98      CALL lbc_lnk( ua, 'U', -1. )
99      CALL lbc_lnk( va, 'V', -1. )
100
101      ! Next velocity
102      ! -------------
103#if defined key_dynspg_flt
104      ! Leap-frog time stepping already done in dynspg_flt.F routine
105#else
106      IF( lk_vvl ) THEN                                  ! Varying levels
107         !RB_vvl scale factors are wrong at this point
108         DO jk = 1, jpkm1
109            ua(ji,jj,jk) = (        ub(:,:,jk) * fse3u(:,:,jk)     &
110               &           + z2dt * ua(:,:,jk) * fse3u(:,:,jk) )   &
111               &         / fse3u(:,:,jk) * umask(:,:,jk)
112            va(ji,jj,jk) = (        vb(:,:,jk) * fse3v(:,:,jk)     &
113               &           + z2dt * va(:,:,jk) * fse3v(:,:,jk) )   &
114               &         / fse3v(:,:,jk) * vmask(:,:,jk)
115         END DO
116      ELSE
117         DO jk = 1, jpkm1
118                  ! Leap-frog time stepping
119            ua(:,:,jk) = ( ub(:,:,jk) + z2dt * ua(:,:,jk) ) * umask(:,:,jk)
120            va(:,:,jk) = ( vb(:,:,jk) + z2dt * va(:,:,jk) ) * vmask(:,:,jk)
121         END DO
122      ENDIF
123
124# if defined key_obc
125      ! Update (ua,va) along open boundaries (only in the rigid-lid case)
126      CALL obc_dyn( kt )
127
128      IF ( lk_dynspg_exp .OR. lk_dynspg_ts ) THEN
129         !Flather boundary condition :
130         !        - Update sea surface height on each open boundary
131         !                 sshn (= after ssh) for explicit case
132         !                 sshn_b (= after ssha_b) for time-splitting case
133         !        - Correct the barotropic velocities
134         CALL obc_dyn_bt( kt )
135
136         !Boundary conditions on sshn ( after ssh)
137         CALL lbc_lnk( sshn, 'T', 1. )
138
139         IF(ln_ctl) THEN         ! print sum trends (used for debugging)
140            CALL prt_ctl(tab2d_1=sshn, clinfo1=' ssh      : ', mask1=tmask)
141         ENDIF
142
143         IF ( ln_vol_cst ) CALL obc_vol( kt )
144
145      ENDIF
146
147# elif defined key_bdy 
148      ! Update (ua,va) along open boundaries (for exp or ts options).
149      IF ( lk_dynspg_exp .or. lk_dynspg_ts ) THEN
150 
151         CALL bdy_dyn_frs( kt )
152
153         IF ( ln_bdy_fla ) THEN
154
155            ua_e(:,:)=0.0
156            va_e(:,:)=0.0
157
158            ! Set these variables for use in bdy_dyn_fla
159            hu_e(:,:) = hu(:,:)
160            hv_e(:,:) = hv(:,:)
161
162            DO jk = 1, jpkm1
163               !! Vertically integrated momentum trends
164               ua_e(:,:) = ua_e(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ua(:,:,jk)
165               va_e(:,:) = va_e(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * va(:,:,jk)
166            END DO
167
168            DO jk = 1 , jpkm1
169               ua(:,:,jk) = ua(:,:,jk) - ua_e(:,:) * hur(:,:)
170               va(:,:,jk) = va(:,:,jk) - va_e(:,:) * hvr(:,:)
171            END DO
172
173            CALL bdy_dta_bt( kt+1, 0)
174            CALL bdy_dyn_fla
175
176         ENDIF
177
178         DO jk = 1 , jpkm1
179            ua(:,:,jk) = ua(:,:,jk) + ua_e(:,:) * hur(:,:)
180            va(:,:,jk) = va(:,:,jk) + va_e(:,:) * hvr(:,:)
181         END DO
182
183      ENDIF
184
185# endif
186# if defined key_agrif
187      CALL Agrif_dyn( kt )
188# endif
189#endif
190
191      ! Time filter and swap of dynamics arrays
192      ! ------------------------------------------
193      IF( neuler == 0 .AND. kt == nit000 ) THEN
194         DO jk = 1, jpkm1                                 
195            ub(:,:,jk) = un(:,:,jk) 
196            vb(:,:,jk) = vn(:,:,jk)
197            un(:,:,jk) = ua(:,:,jk)
198            vn(:,:,jk) = va(:,:,jk)
199         END DO
200      ELSE
201         IF( lk_vvl ) THEN                                ! Varying levels
202            ! Not done
203         ELSE                                             ! Fixed levels
204!RB_vvl : should be done as in tranxt ?
205            DO jk = 1, jpkm1                              ! filter applied on velocities
206               ub(:,:,jk) = atfp * ( ub(:,:,jk) + ua(:,:,jk) ) + atfp1 * un(:,:,jk)
207               vb(:,:,jk) = atfp * ( vb(:,:,jk) + va(:,:,jk) ) + atfp1 * vn(:,:,jk)
208               un(:,:,jk) = ua(:,:,jk)
209               vn(:,:,jk) = va(:,:,jk)
210            END DO
211         ENDIF
212      ENDIF
213
214      IF(ln_ctl)   THEN
215         CALL prt_ctl(tab3d_1=un, clinfo1=' nxt  - Un: ', mask1=umask, &
216            &         tab3d_2=vn, clinfo2=' Vn: ', mask2=vmask)
217      ENDIF
218
219#if defined key_agrif
220      IF (.NOT.Agrif_Root())    CALL Agrif_Update_Dyn( kt )
221#endif     
222
223   END SUBROUTINE dyn_nxt
224
225   !!======================================================================
226END MODULE dynnxt
Note: See TracBrowser for help on using the repository browser.