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.
obcdyn2d.F90 in branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC – NEMO

source: branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcdyn2d.F90 @ 2800

Last change on this file since 2800 was 2800, checked in by davestorkey, 13 years ago
  1. Application of boundary conditions to barotropic and baroclinic velocities clearly separated.
  2. Option to input full velocities in boundary data (default expects barotropic and baroclinic velocities separately).
  3. Option to use initial conditions as boundary conditions coded.
File size: 9.2 KB
Line 
1MODULE obcdyn2d
2   !!======================================================================
3   !!                       ***  MODULE  obcdyn  ***
4   !! Unstructured Open Boundary Cond. :   Flow relaxation scheme on velocities
5   !!======================================================================
6   !! History :  1.0  !  2005-02  (J. Chanut, A. Sellar)  Original code
7   !!             -   !  2007-07  (D. Storkey) Move Flather implementation to separate routine.
8   !!            3.0  !  2008-04  (NEMO team)  add in the reference version
9   !!            3.2  !  2008-04  (R. Benshila) consider velocity instead of transport
10   !!            3.3  !  2010-09  (E.O'Dea) modifications for Shelf configurations
11   !!            3.3  !  2010-09  (D.Storkey) add ice boundary conditions
12   !!----------------------------------------------------------------------
13#if defined key_obc 
14   !!----------------------------------------------------------------------
15   !!   'key_obc' :                    Unstructured Open Boundary Condition
16   !!----------------------------------------------------------------------
17   !!   obc_dyn2d      : Apply open boundary conditions to barotropic variables.
18   !!   obc_dyn2d_fla    : Apply Flather condition
19   !!----------------------------------------------------------------------
20   USE oce             ! ocean dynamics and tracers
21   USE dom_oce         ! ocean space and time domain
22   USE obc_oce         ! ocean open boundary conditions
23   USE dynspg_oce      ! for barotropic variables
24   USE phycst          ! physical constants
25   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
26   USE obctides        ! for tidal harmonic forcing at boundary
27   USE in_out_manager  !
28
29   IMPLICIT NONE
30   PRIVATE
31
32   PUBLIC   obc_dyn2d     ! routine called in dynspg_ts (time splitting case ONLY)
33
34   !!----------------------------------------------------------------------
35   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
36   !! $Id: obcdyn.F90 2528 2010-12-27 17:33:53Z rblod $
37   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
38   !!----------------------------------------------------------------------
39CONTAINS
40
41   SUBROUTINE obc_dyn2d( kt )
42      !!----------------------------------------------------------------------
43      !!                  ***  SUBROUTINE obc_dyn2d  ***
44      !!
45      !! ** Purpose : - Apply open boundary conditions for barotropic variables
46      !!
47      !!----------------------------------------------------------------------
48      INTEGER,                      INTENT(in) ::   kt   ! Main time step counter
49      !!
50      INTEGER                                  ::   ib_obc ! Loop counter
51
52      DO ib_obc=1, nb_obc
53
54         SELECT CASE( nn_dyn2d(ib_obc) )
55         CASE(jp_none)
56            CYCLE
57         CASE(jp_frs)
58            CALL obc_dyn2d_frs( idx_obc(ib_obc), dta_obc(ib_obc), kt )
59         CASE(jp_flather)
60            CALL obc_dyn2d_fla( idx_obc(ib_obc), dta_obc(ib_obc) )
61         CASE DEFAULT
62            CALL ctl_stop( 'obc_dyn3d : unrecognised option for open boundaries for barotropic variables' )
63         END SELECT
64      ENDDO
65
66   END SUBROUTINE obc_dyn2d
67
68   SUBROUTINE obc_dyn2d_frs( idx, dta, kt )
69      !!----------------------------------------------------------------------
70      !!                  ***  SUBROUTINE obc_dyn2d_frs  ***
71      !!
72      !! ** Purpose : - Apply the Flow Relaxation Scheme for barotropic velocities
73      !!                at open boundaries.
74      !!
75      !! References :- Engedahl H., 1995: Use of the flow relaxation scheme in
76      !!               a three-dimensional baroclinic ocean model with realistic
77      !!               topography. Tellus, 365-382.
78      !!----------------------------------------------------------------------
79      INTEGER,         INTENT(in) :: kt
80      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
81      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
82      !!
83      INTEGER  ::   jb, jk         ! dummy loop indices
84      INTEGER  ::   ii, ij, igrd   ! local integers
85      REAL(wp) ::   zwgt           ! boundary weight
86      !!----------------------------------------------------------------------
87      !
88      !
89      igrd = 2                      ! Relaxation of zonal velocity
90      DO jb = 1, idx%nblen(igrd)
91         ii   = idx%nbi(jb,igrd)
92         ij   = idx%nbj(jb,igrd)
93         zwgt = idx%nbw(jb,igrd)
94         pu2d(ii,ij) = ( pu2d(ii,ij) + zwgt * ( dta%u2d(jb) - pu2d(ii,ij) ) ) * umask(ii,ij,1)
95      END DO
96      !
97      igrd = 3                      ! Relaxation of meridional velocity
98      DO jb = 1, idx%nblen(igrd)
99         ii   = idx%nbi(jb,igrd)
100         ij   = idx%nbj(jb,igrd)
101         zwgt = idx%nbw(jb,igrd)
102         pv2d(ii,ij) = ( pv2d(ii,ij) + zwgt * ( dta%v2d(jb) - pv2d(ii,ij) ) ) * vmask(ii,ij,1)
103      END DO
104      CALL lbc_lnk( pu2d, 'U', -1. ) 
105      CALL lbc_lnk( pv2d, 'V', -1. )   ! Boundary points should be updated
106      !
107
108   END SUBROUTINE obc_dyn2d_frs
109
110
111   SUBROUTINE obc_dyn2d_fla( idx, dta )
112      !!----------------------------------------------------------------------
113      !!                 ***  SUBROUTINE obc_dyn2d_fla  ***
114      !!             
115      !!              - Apply Flather boundary conditions on normal barotropic velocities
116      !!                (ln_dyn_fla=.true. or ln_tides=.true.)
117      !!
118      !! ** WARNINGS about FLATHER implementation:
119      !!1. According to Palma and Matano, 1998 "after ssh" is used.
120      !!   In ROMS and POM implementations, it is "now ssh". In the current
121      !!   implementation (tested only in the EEL-R5 conf.), both cases were unstable.
122      !!   So I use "before ssh" in the following.
123      !!
124      !!2. We assume that the normal ssh gradient at the obc is zero. As a matter of
125      !!   fact, the model ssh just inside the dynamical boundary is used (the outside 
126      !!   ssh in the code is not updated).
127      !!
128      !! References:  Flather, R. A., 1976: A tidal model of the northwest European
129      !!              continental shelf. Mem. Soc. R. Sci. Liege, Ser. 6,10, 141-164.     
130      !!----------------------------------------------------------------------
131      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
132      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data
133      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh
134
135      INTEGER  ::   jb, igrd                         ! dummy loop indices
136      INTEGER  ::   ii, ij, iim1, iip1, ijm1, ijp1   ! 2D addresses
137      REAL(wp) ::   zcorr                            ! Flather correction
138      REAL(wp) ::   zforc                            ! temporary scalar
139      !!----------------------------------------------------------------------
140
141      ! ---------------------------------!
142      ! Flather boundary conditions     :!
143      ! ---------------------------------!
144     
145!!!!!!!!!!!! SOME WORK TO DO ON THE TIDES HERE !!!!!!!!!!!!!
146
147!!! REPLACE spgu with nemo_wrk work space
148
149      ! Fill temporary array with ssh data (here spgu):
150      igrd = 1
151      spgu(:,:) = 0.0
152      DO jb = 1, idx%nblenrim(igrd)
153         ii = idx%nbi(jb,igrd)
154         ij = idx%nbj(jb,igrd)
155         spgu(ii, ij) = dta%ssh(jb)
156!!$         IF( ln_tides )   spgu(ii, ij) = spgu(ii, ij) + sshtide(jb)
157      END DO
158      !
159      igrd = 2      ! Flather bc on u-velocity;
160      !             ! remember that flagu=-1 if normal velocity direction is outward
161      !             ! I think we should rather use after ssh ?
162      DO jb = 1, idx%nblenrim(igrd)
163         ii  = idx%nbi(jb,igrd)
164         ij  = idx%nbj(jb,igrd) 
165         iim1 = ii + MAX( 0, INT( idx%flagu(jb) ) )   ! T pts i-indice inside the boundary
166         iip1 = ii - MIN( 0, INT( idx%flagu(jb) ) )   ! T pts i-indice outside the boundary
167         !
168         zcorr = - idx%flagu(jb) * SQRT( grav * phur(ii, ij) ) * ( pssh(iim1, ij) - spgu(iip1,ij) )
169!!$         zforc = dta%u2d(jb) + utide(jb)
170         zforc = dta%u2d(jb)
171         pu2d(ii,ij) = zforc + zcorr * umask(ii,ij,1) 
172      END DO
173      !
174      igrd = 3      ! Flather bc on v-velocity
175      !             ! remember that flagv=-1 if normal velocity direction is outward
176      DO jb = 1, idx%nblenrim(igrd)
177         ii  = idx%nbi(jb,igrd)
178         ij  = idx%nbj(jb,igrd) 
179         ijm1 = ij + MAX( 0, INT( idx%flagv(jb) ) )   ! T pts j-indice inside the boundary
180         ijp1 = ij - MIN( 0, INT( idx%flagv(jb) ) )   ! T pts j-indice outside the boundary
181         !
182         zcorr = - idx%flagv(jb) * SQRT( grav * phvr(ii, ij) ) * ( pssh(ii, ijm1) - spgu(ii,ijp1) )
183!!$         zforc = dta%v2d(jb) + vtide(jb)
184         zforc = dta%v2d(jb)
185         pv2d(ii,ij) = zforc + zcorr * vmask(ii,ij,1)
186      END DO
187      CALL lbc_lnk( pu2d, 'U', -1. )   ! Boundary points should be updated
188      CALL lbc_lnk( pv2d, 'V', -1. )   !
189      !
190      !
191   END SUBROUTINE obc_dyn2d_fla
192#else
193   !!----------------------------------------------------------------------
194   !!   Dummy module                   NO Unstruct Open Boundary Conditions
195   !!----------------------------------------------------------------------
196CONTAINS
197   SUBROUTINE obc_dyn2d( kt )      ! Empty routine
198      WRITE(*,*) 'obc_dyn_frs: You should not have seen this print! error?', kt
199   END SUBROUTINE obc_dyn2d
200#endif
201
202   !!======================================================================
203END MODULE obcdyn2d
Note: See TracBrowser for help on using the repository browser.