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

Last change on this file since 2814 was 2814, checked in by davestorkey, 13 years ago
  1. Implement tidal harmonics forcing (UKMO version) in new structure.
  2. Other bug fixes and updates.
File size: 8.9 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 in_out_manager  !
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC   obc_dyn2d     ! routine called in dynspg_ts and dyn_nxt and dynspg_flt
32
33   !!----------------------------------------------------------------------
34   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
35   !! $Id: obcdyn.F90 2528 2010-12-27 17:33:53Z rblod $
36   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
37   !!----------------------------------------------------------------------
38CONTAINS
39
40   SUBROUTINE obc_dyn2d( kt )
41      !!----------------------------------------------------------------------
42      !!                  ***  SUBROUTINE obc_dyn2d  ***
43      !!
44      !! ** Purpose : - Apply open boundary conditions for barotropic variables
45      !!
46      !!----------------------------------------------------------------------
47      INTEGER,                      INTENT(in) ::   kt   ! Main time step counter
48      !!
49      INTEGER                                  ::   ib_obc ! Loop counter
50
51      DO ib_obc=1, nb_obc
52
53         SELECT CASE( nn_dyn2d(ib_obc) )
54         CASE(jp_none)
55            CYCLE
56         CASE(jp_frs)
57            CALL obc_dyn2d_frs( idx_obc(ib_obc), dta_obc(ib_obc), kt )
58         CASE(jp_flather)
59            CALL obc_dyn2d_fla( idx_obc(ib_obc), dta_obc(ib_obc) )
60         CASE DEFAULT
61            CALL ctl_stop( 'obc_dyn3d : unrecognised option for open boundaries for barotropic variables' )
62         END SELECT
63      ENDDO
64
65   END SUBROUTINE obc_dyn2d
66
67   SUBROUTINE obc_dyn2d_frs( idx, dta, kt )
68      !!----------------------------------------------------------------------
69      !!                  ***  SUBROUTINE obc_dyn2d_frs  ***
70      !!
71      !! ** Purpose : - Apply the Flow Relaxation Scheme for barotropic velocities
72      !!                at open boundaries.
73      !!
74      !! References :- Engedahl H., 1995: Use of the flow relaxation scheme in
75      !!               a three-dimensional baroclinic ocean model with realistic
76      !!               topography. Tellus, 365-382.
77      !!----------------------------------------------------------------------
78      INTEGER,         INTENT(in) ::   kt
79      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
80      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
81      !!
82      INTEGER  ::   jb, jk         ! dummy loop indices
83      INTEGER  ::   ii, ij, igrd   ! local integers
84      REAL(wp) ::   zwgt           ! boundary weight
85      !!----------------------------------------------------------------------
86      !
87      !
88      igrd = 2                      ! Relaxation of zonal velocity
89      DO jb = 1, idx%nblen(igrd)
90         ii   = idx%nbi(jb,igrd)
91         ij   = idx%nbj(jb,igrd)
92         zwgt = idx%nbw(jb,igrd)
93         pu2d(ii,ij) = ( pu2d(ii,ij) + zwgt * ( dta%u2d(jb) - pu2d(ii,ij) ) ) * umask(ii,ij,1)
94      END DO
95      !
96      igrd = 3                      ! Relaxation of meridional velocity
97      DO jb = 1, idx%nblen(igrd)
98         ii   = idx%nbi(jb,igrd)
99         ij   = idx%nbj(jb,igrd)
100         zwgt = idx%nbw(jb,igrd)
101         pv2d(ii,ij) = ( pv2d(ii,ij) + zwgt * ( dta%v2d(jb) - pv2d(ii,ij) ) ) * vmask(ii,ij,1)
102      END DO
103      CALL lbc_lnk( pu2d, 'U', -1. ) 
104      CALL lbc_lnk( pv2d, 'V', -1. )   ! Boundary points should be updated
105      !
106
107   END SUBROUTINE obc_dyn2d_frs
108
109
110   SUBROUTINE obc_dyn2d_fla( idx, dta )
111      !!----------------------------------------------------------------------
112      !!                 ***  SUBROUTINE obc_dyn2d_fla  ***
113      !!             
114      !!              - Apply Flather boundary conditions on normal barotropic velocities
115      !!
116      !! ** WARNINGS about FLATHER implementation:
117      !!1. According to Palma and Matano, 1998 "after ssh" is used.
118      !!   In ROMS and POM implementations, it is "now ssh". In the current
119      !!   implementation (tested only in the EEL-R5 conf.), both cases were unstable.
120      !!   So I use "before ssh" in the following.
121      !!
122      !!2. We assume that the normal ssh gradient at the obc is zero. As a matter of
123      !!   fact, the model ssh just inside the dynamical boundary is used (the outside 
124      !!   ssh in the code is not updated).
125      !!
126      !! References:  Flather, R. A., 1976: A tidal model of the northwest European
127      !!              continental shelf. Mem. Soc. R. Sci. Liege, Ser. 6,10, 141-164.     
128      !!----------------------------------------------------------------------
129      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
130      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data
131      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh
132
133      INTEGER  ::   jb, igrd                         ! dummy loop indices
134      INTEGER  ::   ii, ij, iim1, iip1, ijm1, ijp1   ! 2D addresses
135      REAL(wp) ::   zcorr                            ! Flather correction
136      REAL(wp) ::   zforc                            ! temporary scalar
137      !!----------------------------------------------------------------------
138
139      ! ---------------------------------!
140      ! Flather boundary conditions     :!
141      ! ---------------------------------!
142     
143!!! REPLACE spgu with nemo_wrk work space
144
145      ! Fill temporary array with ssh data (here spgu):
146      igrd = 1
147      spgu(:,:) = 0.0
148      DO jb = 1, idx%nblenrim(igrd)
149         ii = idx%nbi(jb,igrd)
150         ij = idx%nbj(jb,igrd)
151         spgu(ii, ij) = dta%ssh(jb)
152      END DO
153      !
154      igrd = 2      ! Flather bc on u-velocity;
155      !             ! remember that flagu=-1 if normal velocity direction is outward
156      !             ! I think we should rather use after ssh ?
157      DO jb = 1, idx%nblenrim(igrd)
158         ii  = idx%nbi(jb,igrd)
159         ij  = idx%nbj(jb,igrd) 
160         iim1 = ii + MAX( 0, INT( idx%flagu(jb) ) )   ! T pts i-indice inside the boundary
161         iip1 = ii - MIN( 0, INT( idx%flagu(jb) ) )   ! T pts i-indice outside the boundary
162         !
163         zcorr = - idx%flagu(jb) * SQRT( grav * phur(ii, ij) ) * ( pssh(iim1, ij) - spgu(iip1,ij) )
164         zforc = dta%u2d(jb)
165         pu2d(ii,ij) = zforc + zcorr * umask(ii,ij,1) 
166      END DO
167      !
168      igrd = 3      ! Flather bc on v-velocity
169      !             ! remember that flagv=-1 if normal velocity direction is outward
170      DO jb = 1, idx%nblenrim(igrd)
171         ii  = idx%nbi(jb,igrd)
172         ij  = idx%nbj(jb,igrd) 
173         ijm1 = ij + MAX( 0, INT( idx%flagv(jb) ) )   ! T pts j-indice inside the boundary
174         ijp1 = ij - MIN( 0, INT( idx%flagv(jb) ) )   ! T pts j-indice outside the boundary
175         !
176         zcorr = - idx%flagv(jb) * SQRT( grav * phvr(ii, ij) ) * ( pssh(ii, ijm1) - spgu(ii,ijp1) )
177         zforc = dta%v2d(jb)
178         pv2d(ii,ij) = zforc + zcorr * vmask(ii,ij,1)
179      END DO
180      CALL lbc_lnk( pu2d, 'U', -1. )   ! Boundary points should be updated
181      CALL lbc_lnk( pv2d, 'V', -1. )   !
182      !
183      !
184   END SUBROUTINE obc_dyn2d_fla
185#else
186   !!----------------------------------------------------------------------
187   !!   Dummy module                   NO Unstruct Open Boundary Conditions
188   !!----------------------------------------------------------------------
189CONTAINS
190   SUBROUTINE obc_dyn2d( kt )      ! Empty routine
191      WRITE(*,*) 'obc_dyn_frs: You should not have seen this print! error?', kt
192   END SUBROUTINE obc_dyn2d
193#endif
194
195   !!======================================================================
196END MODULE obcdyn2d
Note: See TracBrowser for help on using the repository browser.