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.
wzvmod.F90 in trunk/NEMO/OPA_SRC/DYN – NEMO

source: trunk/NEMO/OPA_SRC/DYN/wzvmod.F90 @ 911

Last change on this file since 911 was 911, checked in by ctlod, 16 years ago

Implementation of the BDY package, see ticket: #126
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.3 KB
Line 
1MODULE wzvmod
2   !!==============================================================================
3   !!                       ***  MODULE  wzvmod  ***
4   !! Ocean diagnostic variable : vertical velocity
5   !!==============================================================================
6   !! History :   5.0  !  90-10  (C. Levy, G. Madec)  Original code
7   !!             7.0  !  96-01  (G. Madec)  Statement function for e3
8   !!             8.5  !  02-07  (G. Madec)  Free form, F90
9   !!              "   !  07-07  (D. Storkey) Zero zhdiv at open boundary (BDY)
10   !!----------------------------------------------------------------------
11   !!   wzv        : Compute the vertical velocity
12   !!----------------------------------------------------------------------
13   !! * Modules used
14   USE oce             ! ocean dynamics and tracers variables
15   USE dom_oce         ! ocean space and time domain variables
16   USE sbc_oce         ! surface boundary condition: ocean
17   USE domvvl          ! Variable volume
18   USE in_out_manager  ! I/O manager
19   USE prtctl          ! Print control
20   USE phycst
21   USE bdy_oce         ! unstructured open boundaries
22   USE lbclnk          ! ocean lateral boundary condition (or mpp link)
23
24   IMPLICIT NONE
25   PRIVATE
26
27   !! * Routine accessibility
28   PUBLIC wzv          ! routine called by step.F90 and inidtr.F90
29
30   !! * Substitutions
31#  include "domzgr_substitute.h90"
32   !!----------------------------------------------------------------------
33   !!  OPA 9.0 , LOCEAN-IPSL (2005)
34   !! $Id$
35   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
36   !!----------------------------------------------------------------------
37
38CONTAINS
39
40   SUBROUTINE wzv( kt )
41      !!----------------------------------------------------------------------
42      !!                    ***  ROUTINE wzv  ***
43      !!
44      !! ** Purpose :   Compute the now vertical velocity after the array swap
45      !!
46      !! ** Method  :   Using the incompressibility hypothesis, the vertical
47      !!      velocity is computed by integrating the horizontal divergence
48      !!      from the bottom to the surface.
49      !!        The boundary conditions are w=0 at the bottom (no flux) and,
50      !!      in regid-lid case, w=0 at the sea surface.
51      !!
52      !! ** action  :   wn array : the now vertical velocity
53      !!----------------------------------------------------------------------
54      !! * Arguments
55      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
56
57      !! * Local declarations
58      INTEGER  ::     jgrd, jb           ! temporary scalars
59      INTEGER  ::           jk           ! dummy loop indices
60      !! Variable volume
61      INTEGER  ::   ji, jj               ! dummy loop indices
62      REAL(wp) ::   z2dt, zraur          ! temporary scalar
63      REAL(wp), DIMENSION (jpi,jpj) ::   zssha, zun, zvn, zhdiv
64      !!----------------------------------------------------------------------
65
66      IF( kt == nit000 ) THEN
67         IF(lwp) WRITE(numout,*)
68         IF(lwp) WRITE(numout,*) 'wzv     : vertical velocity from continuity eq.'
69         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
70
71         ! bottom boundary condition: w=0 (set once for all)
72         wn(:,:,jpk) = 0.e0
73      ENDIF
74
75      IF( lk_vvl ) THEN                ! Variable volume
76         !
77         z2dt = 2. * rdt                                       ! time step: leap-frog
78         IF( neuler == 0 .AND. kt == nit000 ) z2dt = rdt       ! time step: Euler if restart from rest
79         zraur  = 1. / rauw
80
81         ! Vertically integrated quantities
82         ! --------------------------------
83         zun(:,:) = 0.e0
84         zvn(:,:) = 0.e0
85         !
86         DO jk = 1, jpkm1             ! Vertically integrated transports (now)
87            zun(:,:) = zun(:,:) + fse3u(:,:,jk) * un(:,:,jk)
88            zvn(:,:) = zvn(:,:) + fse3v(:,:,jk) * vn(:,:,jk)
89         END DO
90
91         ! Horizontal divergence of barotropic transports
92         !--------------------------------------------------
93         zhdiv(:,:) = 0.e0
94         DO jj = 2, jpjm1
95            DO ji = 2, jpim1   ! vector opt.
96               zhdiv(ji,jj) = (  e2u(ji  ,jj  ) * zun(ji  ,jj  )     &
97                  &            - e2u(ji-1,jj  ) * zun(ji-1,jj  )     &
98                  &            + e1v(ji  ,jj  ) * zvn(ji  ,jj  )     &
99                  &            - e1v(ji  ,jj-1) * zvn(ji  ,jj-1) )   &
100                  &           / ( e1t(ji,jj) * e2t(ji,jj) )
101            END DO
102         END DO
103
104#if defined key_obc && ( key_dynspg_exp || key_dynspg_ts )
105         ! open boundaries (div must be zero behind the open boundary)
106         !  mpp remark: The zeroing of hdiv can probably be extended to 1->jpi/jpj for the correct row/column
107         IF( lp_obc_east  )   zhdiv(nie0p1:nie1p1,nje0  :nje1)   = 0.e0    ! east
108         IF( lp_obc_west  )   zhdiv(niw0  :niw1  ,njw0  :njw1)   = 0.e0    ! west
109         IF( lp_obc_north )   zhdiv(nin0  :nin1  ,njn0p1:njn1p1) = 0.e0    ! north
110         IF( lp_obc_south )   zhdiv(nis0  :nis1  ,njs0  :njs1)   = 0.e0    ! south
111#endif
112
113#if defined key_bdy || defined key_bdy_tides
114         jgrd=1 !: tracer grid.
115         DO jb = 1, nblenrim(jgrd)
116           ji = nbi(jb,jgrd)
117           jj = nbj(jb,jgrd)
118           zhdiv(ji, jj) = 0.e0
119           zhdiv(ji, jj) = 0.e0
120           zhdiv(ji, jj) = 0.e0
121         END DO
122#endif
123
124         CALL lbc_lnk( zhdiv, 'T', 1. )
125
126         ! Sea surface elevation time stepping
127         ! -----------------------------------
128         zssha(:,:) = sshb(:,:) - z2dt * ( zraur * emp(:,:) + zhdiv(:,:) ) * tmask(:,:,1)
129
130         ! Vertical velocity computed from bottom
131         ! --------------------------------------
132         DO jk = jpkm1, 1, -1
133            wn(:,:,jk) = wn(:,:,jk+1) - fse3t(:,:,jk) * hdivn(:,:,jk) &
134              &        - ( zssha(:,:) - sshb(:,:) ) * fsve3t(:,:,jk) * mut(:,:,jk) / z2dt
135         END DO
136
137      ELSE                             ! Fixed volume
138
139         ! Vertical velocity computed from bottom
140         ! --------------------------------------
141         DO jk = jpkm1, 1, -1
142            wn(:,:,jk) = wn(:,:,jk+1) - fse3t(:,:,jk) * hdivn(:,:,jk)
143         END DO
144     
145      ENDIF
146
147      IF(ln_ctl)   CALL prt_ctl(tab3d_1=wn, clinfo1=' w**2 -   : ', mask1=wn)
148
149   END SUBROUTINE wzv
150
151   !!======================================================================
152END MODULE wzvmod
Note: See TracBrowser for help on using the repository browser.