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

source: trunk/NEMO/OPA_SRC/OBC/obcfla.F90 @ 475

Last change on this file since 475 was 475, checked in by opalod, 18 years ago

nemo_v1_bugfix_043 : CT : tangential velocities on the east and north boundaries must be computed over ghost cells instead of internal values

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 9.5 KB
Line 
1MODULE obcfla
2#if defined key_obc && defined key_dynspg_ts
3   !!=================================================================================
4   !!                       ***  MODULE  obcfla  ***
5   !! Ocean dynamics:   Flather's algorithm at open boundaries for the time-splitting
6   !!=================================================================================
7
8   !!---------------------------------------------------------------------------------
9   !!   obc_fla_ts        : call the subroutine for each open boundary
10   !!   obc_fla_ts_east   : Flather on the east  open boundary velocities & ssh
11   !!   obc_fla_ts_west   : Flather on the west  open boundary velocities & ssh
12   !!   obc_fla_ts_north  : Flather on the north open boundary velocities & ssh
13   !!   obc_fla_ts_south  : Flather on the south open boundary velocities & ssh
14   !!----------------------------------------------------------------------------------
15
16   !!----------------------------------------------------------------------------------
17   !! * Modules used
18   USE oce             ! ocean dynamics and tracers
19   USE dom_oce         ! ocean space and time domain
20   USE dynspg_oce      ! surface pressure gradient variables
21   USE phycst          ! physical constants
22   USE obc_oce         ! ocean open boundary conditions
23   USE obcdta          ! ocean open boundary conditions: climatology
24
25   IMPLICIT NONE
26   PRIVATE
27
28   !! * Accessibility
29   PUBLIC obc_fla_ts  ! routine called in dynspg_ts (free surface time splitting case)
30
31   !!---------------------------------------------------------------------------------
32   !!  OPA 9.0 , LOCEAN-IPSL (2005)
33   !! $Header$
34   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
35   !!---------------------------------------------------------------------------------
36
37CONTAINS
38
39   SUBROUTINE obc_fla_ts
40      !!------------------------------------------------------------------------------
41      !!                      SUBROUTINE obc_fla_ts
42      !!                     **********************
43      !! ** Purpose :
44      !!      Apply Flather's algorithm at open boundaries for the time-splitting
45      !!      free surface case (barotropic variables)
46      !!
47      !!      This routine is called in dynspg_ts.F90 routine
48      !!
49      !!      The logical variable lp_obc_east, and/or lp_obc_west, and/or lp_obc_north,
50      !!      and/or lp_obc_south allow the user to determine which boundary is an
51      !!      open one (must be done in the obc_par.F90 file).
52      !!
53      !! ** Reference :
54      !!         Flather, R. A., 1976, Mem. Soc. R. Sci. Liege, Ser. 6, 10, 141-164
55      !!
56      !! History :
57      !!   9.0  !  05-12  (V. Garnier) original
58      !!------------------------------------------------------------------------------
59
60      IF( lp_obc_east  )   CALL obc_fla_ts_east 
61      IF( lp_obc_west  )   CALL obc_fla_ts_west 
62      IF( lp_obc_north )   CALL obc_fla_ts_north
63      IF( lp_obc_south )   CALL obc_fla_ts_south
64
65   END SUBROUTINE obc_fla_ts
66
67
68   SUBROUTINE obc_fla_ts_east
69      !!------------------------------------------------------------------------------
70      !!                  ***  SUBROUTINE obc_fla_ts_east  ***
71      !!
72      !! ** Purpose :
73      !!      Apply Flather's algorithm on east OBC velocities ua, va
74      !!      Fix sea surface height (sshn_e) on east open boundary
75      !!
76      !!  History :
77      !!   9.0  !  05-12  (V. Garnier) original
78      !!------------------------------------------------------------------------------
79      !! * Local declaration
80      INTEGER ::   ji, jj, jk ! dummy loop indices
81      !!------------------------------------------------------------------------------
82
83      DO ji = nie0, nie1
84         DO jk = 1, jpkm1
85            DO jj = 1, jpj
86               ua_e(ji,jj) = (  ubtfoe(jj) + sqrt( grav*hu(ji,jj) )           &
87                  &            * ( ( sshn_e(ji,jj) + sshn_e(ji+1,jj) ) * 0.5  &
88                  &            - sshfoe(jj) )  ) * uemsk(jj,jk)
89            END DO
90         END DO
91      END DO
92      DO ji = nie0p1, nie1p1
93         DO jj = 1, jpj
94            sshfoe_b(ji,jj) = sshfoe_b(ji,jj) + sqrt( grav*hur(ji,jj) )     &
95               &             * ( ( sshn_e(ji,jj) + sshn_e(ji+1,jj) ) * 0.5  &
96               &                 - sshfoe(jj) ) * uemsk(jj,1)
97            ssha_e(ji,jj) = ssha_e(ji,jj) * ( 1. - temsk(jj,1) ) &
98               &            + temsk(jj,1) * sshfoe(jj)
99            va_e(ji,jj) = vbtfoe(jj) * uemsk(jj,jk)
100         END DO
101      END DO
102
103   END SUBROUTINE obc_fla_ts_east
104
105
106   SUBROUTINE obc_fla_ts_west
107      !!------------------------------------------------------------------------------
108      !!                  ***  SUBROUTINE obc_fla_ts_west  ***
109      !!
110      !! ** Purpose :
111      !!      Apply Flather's algorithm on west OBC velocities ua, va
112      !!      Fix sea surface height (sshn_e) on west open boundary
113      !!
114      !!  History :
115      !!   9.0  !  05-12  (V. Garnier) original
116      !!------------------------------------------------------------------------------
117      !! * Local declaration
118      INTEGER ::   ji, jj, jk ! dummy loop indices
119      !!------------------------------------------------------------------------------
120
121      DO ji = niw0, niw1
122         DO jk = 1, jpkm1
123            DO jj = 1, jpj
124               ua_e(ji,jj) = ( ubtfow(jj) - sqrt( grav * hu(ji,jj) )          &
125                  &            * ( ( sshn_e(ji,jj) + sshn_e(ji+1,jj) ) * 0.5  &
126                  &                - sshfow(jj) ) ) * uwmsk(jj,jk)
127               va_e(ji,jj) = vbtfow(jj) * uwmsk(jj,jk)
128            END DO
129         END DO
130         DO jj = 1, jpj
131            sshfow_b(ji,jj) = sshfow_b(ji,jj) - sqrt( grav * hur(ji,jj) )     &
132                              * ( ( sshn_e(ji,jj) + sshn_e(ji+1,jj) ) * 0.5   &
133                                 - sshfow(jj) ) * uwmsk(jj,1)
134            ssha_e(ji,jj) = ssha_e(ji,jj) * ( 1. - twmsk(jj,1) ) &
135               &            + twmsk(jj,1)*sshfow(jj)
136         END DO
137      END DO
138
139   END SUBROUTINE obc_fla_ts_west
140
141   SUBROUTINE obc_fla_ts_north
142      !!------------------------------------------------------------------------------
143      !!                     SUBROUTINE obc_fla_ts_north
144      !!                    *************************
145      !! ** Purpose :
146      !!      Apply Flather's algorithm on north OBC velocities ua, va
147      !!      Fix sea surface height (sshn_e) on north open boundary
148      !!
149      !!  History :
150      !!   9.0  !  05-12  (V. Garnier) original
151      !!------------------------------------------------------------------------------
152      !! * Local declaration
153      INTEGER ::   ji, jj, jk ! dummy loop indices
154      !!------------------------------------------------------------------------------
155
156      DO jj = njn0, njn1
157         DO jk = 1, jpkm1
158            DO ji = 1, jpi
159               va_e(ji,jj) = ( vbtfon(ji) + sqrt( grav * hv(ji,jj) )           &
160                  &            * ( ( sshn_e(ji,jj) + sshn_e(ji,jj+1) ) * 0.5   &
161                  &                - sshfon(ji) ) ) * vnmsk(ji,jk)
162            END DO
163         END DO
164      END DO
165      DO jj = njn0p1, njn1p1
166         DO ji = 1, jpi
167            sshfon_b(ji,jj) = sshfon_b(ji,jj) + sqrt( grav * hvr(ji,jj) )  &
168               &              * ( ( sshn_e(ji,jj) + sshn_e(ji,jj+1) ) * 0.5    &
169               &                  - sshfon(ji) ) * vnmsk(ji,1)
170            ssha_e(ji,jj) = ssha_e(ji,jj) * ( 1. - tnmsk(ji,1) ) &
171               &            + sshfon(ji) * tnmsk(ji,1)
172            ua_e(ji,jj) = ubtfon(ji) * vnmsk(ji,jk)
173         END DO
174      END DO
175
176   END SUBROUTINE obc_fla_ts_north
177
178   SUBROUTINE obc_fla_ts_south
179      !!------------------------------------------------------------------------------
180      !!                     SUBROUTINE obc_fla_ts_south
181      !!                    *************************
182      !! ** Purpose :
183      !!      Apply Flather's algorithm on south OBC velocities ua, va
184      !!      Fix sea surface height (sshn_e) on south open boundary
185      !!
186      !!  History :
187      !!   9.0  !  05-12  (V. Garnier) original
188      !!------------------------------------------------------------------------------
189      !! * Local declaration
190      INTEGER ::   ji, jj, jk ! dummy loop indices
191
192      !!------------------------------------------------------------------------------
193
194      DO jj = njs0, njs1
195         DO jk = 1, jpkm1
196            DO ji = 1, jpi
197               va_e(ji,jj) = ( vbtfos(ji) - sqrt( grav * hv(ji,jj) )            &
198                  &            * ( ( sshn_e(ji,jj) + sshn_e(ji,jj+1) ) * 0.5    &
199                  &                - sshfos(ji) ) ) * vsmsk(ji,jk)
200               ua_e(ji,jj) = ubtfos(ji) * vsmsk(ji,jk)
201            END DO
202         END DO
203         DO ji = 1, jpi
204            sshfos_b(ji,jj) = sshfos_b(ji,jj) - sqrt( grav * hvr(ji,jj) )      &
205               &              * ( ( sshn_e(ji,jj) + sshn_e(ji,jj+1) ) * 0.5    &
206               &                  - sshfos(ji) ) * vsmsk(ji,1)
207            ssha_e(ji,jj) = ssha_e(ji,jj) * (1. - tsmsk(ji,1) ) &
208               &            + tsmsk(ji,1) * sshfos(ji)
209         END DO
210      END DO
211
212   END SUBROUTINE obc_fla_ts_south
213#else
214   !!=================================================================================
215   !!                       ***  MODULE  obcfla  ***
216   !! Ocean dynamics:   Flather's algorithm at open boundaries for the time-splitting
217   !!=================================================================================
218CONTAINS
219
220   SUBROUTINE obc_fla_ts
221      WRITE(*,*) 'obc_fla_ts: You should not have seen this print! error?'
222   END SUBROUTINE obc_fla_ts
223#endif
224
225END MODULE obcfla
Note: See TracBrowser for help on using the repository browser.