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.
bdyvol.F90 in branches/2015/dev_merge_2015/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

source: branches/2015/dev_merge_2015/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90 @ 6060

Last change on this file since 6060 was 6060, checked in by timgraham, 8 years ago

Merged dev_r5836_noc2_VVL_BY_DEFAULT into branch

  • Property svn:keywords set to Id
File size: 8.7 KB
Line 
1MODULE bdyvol
2   !!======================================================================
3   !!                       ***  MODULE  bdyvol  ***
4   !! Ocean dynamic :  Volume constraint when unstructured boundary
5   !!                  and filtered free surface are used
6   !!======================================================================
7   !! History :  1.0  !  2005-01  (J. Chanut, A. Sellar)  Original code
8   !!             -   !  2006-01  (J. Chanut) Bug correction
9   !!            3.0  !  2008-04  (NEMO team)  add in the reference version
10   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge
11   !!----------------------------------------------------------------------
12#if defined key_bdy
13   !!----------------------------------------------------------------------
14   !!   'key_bdy'                     unstructured open boundary conditions
15   !!----------------------------------------------------------------------
16   USE oce            ! ocean dynamics and tracers
17   USE bdy_oce        ! ocean open boundary conditions
18   USE sbc_oce        ! ocean surface boundary conditions
19   USE dom_oce        ! ocean space and time domain
20   USE phycst         ! physical constants
21   USE sbcisf         ! ice shelf
22   !
23   USE in_out_manager ! I/O manager
24   USE lib_mpp        ! for mppsum
25   USE timing         ! Timing
26   USE lib_fortran    ! Fortran routines library
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC   bdy_vol    ! called by ???
32
33   !!----------------------------------------------------------------------
34   !! NEMO/OPA 3.7 , NEMO Consortium (2015)
35   !! $Id$
36   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
37   !!----------------------------------------------------------------------
38CONTAINS
39
40   SUBROUTINE bdy_vol( kt )
41      !!----------------------------------------------------------------------
42      !!                      ***  ROUTINE bdyvol  ***
43      !!
44      !! ** Purpose :   This routine controls the volume of the system.
45      !!      A correction velocity is calculated to correct the total transport
46      !!      through the unstructured OBC.
47      !!      The total depth used is constant (H0) to be consistent with the
48      !!      linear free surface coded in OPA 8.2    <<<=== !!gm  ???? true ????
49      !!
50      !! ** Method  :   The correction velocity (zubtpecor here) is defined calculating
51      !!      the total transport through all open boundaries (trans_bdy) minus
52      !!      the cumulate E-P flux (z_cflxemp) divided by the total lateral
53      !!      surface (bdysurftot) of the unstructured boundary.
54      !!         zubtpecor = [trans_bdy - z_cflxemp ]*(1./bdysurftot)
55      !!      with z_cflxemp => sum of (Evaporation minus Precipitation)
56      !!                       over all the domain in m3/s at each time step.
57      !!      z_cflxemp < 0 when precipitation dominate
58      !!      z_cflxemp > 0 when evaporation dominate
59      !!
60      !!      There are 2 options (user's desiderata):
61      !!         1/ The volume changes according to E-P, this is the default
62      !!            option. In this case the cumulate E-P flux are setting to
63      !!            zero (z_cflxemp=0) to calculate the correction velocity. So
64      !!            it will only balance the flux through open boundaries.
65      !!            (set nn_volctl to 0 in tne namelist for this option)
66      !!         2/ The volume is constant even with E-P flux. In this case
67      !!            the correction velocity must balance both the flux
68      !!            through open boundaries and the ones through the free
69      !!            surface.
70      !!            (set nn_volctl to 1 in tne namelist for this option)
71      !!----------------------------------------------------------------------
72      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
73      !
74      INTEGER  ::   ji, jj, jk, jb, jgrd
75      INTEGER  ::   ib_bdy, ii, ij
76      REAL(wp) ::   zubtpecor, z_cflxemp, ztranst
77      TYPE(OBC_INDEX), POINTER :: idx
78      !!-----------------------------------------------------------------------------
79      !
80      IF( nn_timing == 1 )   CALL timing_start('bdy_vol')
81      !
82      IF( ln_vol ) THEN
83      !
84      IF( kt == nit000 ) THEN
85         IF(lwp) WRITE(numout,*)
86         IF(lwp) WRITE(numout,*)'bdy_vol : Correction of velocities along unstructured OBC'
87         IF(lwp) WRITE(numout,*)'~~~~~~~'
88      END IF 
89
90      ! Calculate the cumulate surface Flux z_cflxemp (m3/s) over all the domain
91      ! -----------------------------------------------------------------------
92!!gm replace these lines :
93      z_cflxemp = SUM ( ( emp(:,:)-rnf(:,:)+fwfisf(:,:) ) * bdytmask(:,:) * e1e2t(:,:) ) / rau0
94      IF( lk_mpp )   CALL mpp_sum( z_cflxemp )     ! sum over the global domain
95!!gm   by :
96!!gm      z_cflxemp = glob_sum(  ( emp(:,:)-rnf(:,:)+fwfisf(:,:) ) * bdytmask(:,:) * e1e2t(:,:)  ) / rau0
97!!gm
98
99      ! Transport through the unstructured open boundary
100      ! ------------------------------------------------
101      zubtpecor = 0._wp
102      DO ib_bdy = 1, nb_bdy
103         idx => idx_bdy(ib_bdy)
104         !
105         jgrd = 2                               ! cumulate u component contribution first
106         DO jb = 1, idx%nblenrim(jgrd)
107            DO jk = 1, jpkm1
108               ii = idx%nbi(jb,jgrd)
109               ij = idx%nbj(jb,jgrd)
110               zubtpecor = zubtpecor + idx%flagu(jb,jgrd) * ua(ii,ij, jk) * e2u(ii,ij) * e3u_n(ii,ij,jk)
111            END DO
112         END DO
113         jgrd = 3                               ! then add v component contribution
114         DO jb = 1, idx%nblenrim(jgrd)
115            DO jk = 1, jpkm1
116               ii = idx%nbi(jb,jgrd)
117               ij = idx%nbj(jb,jgrd)
118               zubtpecor = zubtpecor + idx%flagv(jb,jgrd) * va(ii,ij, jk) * e1v(ii,ij) * e3v_n(ii,ij,jk) 
119            END DO
120         END DO
121         !
122      END DO
123      IF( lk_mpp )   CALL mpp_sum( zubtpecor )   ! sum over the global domain
124
125      ! The normal velocity correction
126      ! ------------------------------
127      IF( nn_volctl==1 ) THEN   ;   zubtpecor = ( zubtpecor - z_cflxemp ) / bdysurftot 
128      ELSE                      ;   zubtpecor =   zubtpecor               / bdysurftot
129      END IF
130
131      ! Correction of the total velocity on the unstructured boundary to respect the mass flux conservation
132      ! -------------------------------------------------------------
133      ztranst = 0._wp
134      DO ib_bdy = 1, nb_bdy
135         idx => idx_bdy(ib_bdy)
136         !
137         jgrd = 2                               ! correct u component
138         DO jb = 1, idx%nblenrim(jgrd)
139            DO jk = 1, jpkm1
140               ii = idx%nbi(jb,jgrd)
141               ij = idx%nbj(jb,jgrd)
142               ua(ii,ij,jk) = ua(ii,ij,jk) - idx%flagu(jb,jgrd) * zubtpecor * umask(ii,ij,jk)
143               ztranst = ztranst + idx%flagu(jb,jgrd) * ua(ii,ij,jk) * e2u(ii,ij) * e3u_n(ii,ij,jk)
144            END DO
145         END DO
146         jgrd = 3                              ! correct v component
147         DO jb = 1, idx%nblenrim(jgrd)
148            DO jk = 1, jpkm1
149               ii = idx%nbi(jb,jgrd)
150               ij = idx%nbj(jb,jgrd)
151               va(ii,ij,jk) = va(ii,ij,jk) -idx%flagv(jb,jgrd) * zubtpecor * vmask(ii,ij,jk)
152               ztranst = ztranst + idx%flagv(jb,jgrd) * va(ii,ij,jk) * e1v(ii,ij) * e3v_n(ii,ij,jk)
153            END DO
154         END DO
155         !
156      END DO
157      IF( lk_mpp )   CALL mpp_sum( ztranst )   ! sum over the global domain
158 
159      ! Check the cumulated transport through unstructured OBC once barotropic velocities corrected
160      ! ------------------------------------------------------
161      IF( lwp .AND. MOD( kt, nwrite ) == 0 ) THEN
162         IF(lwp) WRITE(numout,*)
163         IF(lwp) WRITE(numout,*)'bdy_vol : time step :', kt
164         IF(lwp) WRITE(numout,*)'~~~~~~~ '
165         IF(lwp) WRITE(numout,*)'          cumulate flux EMP             =', z_cflxemp  , ' (m3/s)'
166         IF(lwp) WRITE(numout,*)'          total lateral surface of OBC  =', bdysurftot, '(m2)'
167         IF(lwp) WRITE(numout,*)'          correction velocity zubtpecor =', zubtpecor , '(m/s)'
168         IF(lwp) WRITE(numout,*)'          cumulated transport ztranst   =', ztranst   , '(m3/s)'
169      END IF 
170      !
171      IF( nn_timing == 1 )   CALL timing_stop('bdy_vol')
172      !
173      END IF ! ln_vol
174      !
175   END SUBROUTINE bdy_vol
176
177#else
178   !!----------------------------------------------------------------------
179   !!   Dummy module                   NO Unstruct Open Boundary Conditions
180   !!----------------------------------------------------------------------
181CONTAINS
182   SUBROUTINE bdy_vol( kt )        ! Empty routine
183      WRITE(*,*) 'bdy_vol: You should not have seen this print! error?', kt
184   END SUBROUTINE bdy_vol
185#endif
186
187   !!======================================================================
188END MODULE bdyvol
Note: See TracBrowser for help on using the repository browser.