Changeset 2715 for trunk/NEMOGCM/NEMO/OPA_SRC/OBC/obcdyn_bt.F90
- Timestamp:
- 2011-03-30T17:58:35+02:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/OBC/obcdyn_bt.F90
r2528 r2715 1 1 MODULE obcdyn_bt 2 !!====================================================================== 3 !! *** MODULE obcdyn_bt *** 4 !! Ocean dynamics: Radiation/prescription of sea surface heights on each open boundary 5 !!====================================================================== 6 !! History : 1.0 ! 2005-12 (V. Garnier) original code 7 !!---------------------------------------------------------------------- 2 8 #if ( defined key_dynspg_ts || defined key_dynspg_exp ) && defined key_obc 3 !!================================================================================= 4 !! *** MODULE obcdyn_bt *** 5 !! Ocean dynamics: Radiation/prescription of sea surface heights 6 !! on each open boundary 7 !!================================================================================= 8 9 !!--------------------------------------------------------------------------------- 9 !!---------------------------------------------------------------------- 10 !! 'key_dynspg_ts' OR time spliting free surface 11 !! 'key_dynspg_exp' AND explicit free surface 12 !! 'key_obc' Open Boundary Condition 13 !!---------------------------------------------------------------------- 10 14 !! obc_dyn_bt : call the subroutine for each open boundary 11 15 !! obc_dyn_bt_east : Flather's algorithm at the east open boundary … … 13 17 !! obc_dyn_bt_north : Flather's algorithm at the north open boundary 14 18 !! obc_dyn_bt_south : Flather's algorithm at the south open boundary 15 !!---------------------------------------------------------------------------------- 16 17 !!---------------------------------------------------------------------------------- 18 !! * Modules used 19 !!---------------------------------------------------------------------- 19 20 USE oce ! ocean dynamics and tracers 20 21 USE dom_oce ! ocean space and time domain … … 30 31 PRIVATE 31 32 32 !! * Accessibility 33 PUBLIC obc_dyn_bt ! routine called in dynnxt (explicit free surface case) 34 35 !!--------------------------------------------------------------------------------- 33 PUBLIC obc_dyn_bt ! routine called in dynnxt (explicit free surface case) 34 35 !!---------------------------------------------------------------------- 36 36 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 37 37 !! $Id$ 38 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 39 !!---------------------------------------------------------------------- 40 38 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 39 !!---------------------------------------------------------------------- 41 40 CONTAINS 42 41 43 SUBROUTINE obc_dyn_bt ( kt ) 44 !!------------------------------------------------------------------------------ 45 !! SUBROUTINE obc_dyn_bt 46 !! *********************** 47 !! ** Purpose : 48 !! Apply Flather's algorithm at open boundaries for the explicit 49 !! free surface case and free surface case with time-splitting 42 SUBROUTINE obc_dyn_bt( kt ) 43 !!---------------------------------------------------------------------- 44 !! *** SUBROUTINE obc_dyn_bt *** 45 !! 46 !! ** Purpose : Apply Flather's algorithm at open boundaries for the explicit 47 !! free surface case and free surface case with time-splitting 50 48 !! 51 49 !! This routine is called in dynnxt.F routine and updates ua, va and sshn. … … 55 53 !! open one (must be done in the param_obc.h90 file). 56 54 !! 57 !! ** Reference : 58 !! Flather, R. A., 1976, Mem. Soc. R. Sci. Liege, Ser. 6, 10, 141-164 59 !! 60 !! History : 61 !! 9.0 ! 05-12 (V. Garnier) original 62 !!---------------------------------------------------------------------- 63 !! * Arguments 64 INTEGER, INTENT( in ) :: kt 65 55 !! Reference : Flather, R. A., 1976, Mem. Soc. R. Sci. Liege, Ser. 6, 10, 141-164 56 !!---------------------------------------------------------------------- 57 INTEGER, INTENT(in) :: kt 66 58 !!---------------------------------------------------------------------- 67 59 … … 85 77 86 78 # if defined key_dynspg_exp 79 87 80 SUBROUTINE obc_dyn_bt_east 88 !!---------------------------------------------------------------------- --------81 !!---------------------------------------------------------------------- 89 82 !! *** SUBROUTINE obc_dyn_bt_east *** 90 83 !! … … 93 86 !! Fix sea surface height (sshn) on east open boundary 94 87 !! The logical lfbceast must be .TRUE. 95 !! 96 !! History : 97 !! 9.0 ! 05-12 (V. Garnier) original 98 !!------------------------------------------------------------------------------ 99 !! * Local declaration 100 INTEGER :: ji, jj, jk ! dummy loop indices 101 !!------------------------------------------------------------------------------ 88 !!---------------------------------------------------------------------- 89 INTEGER, INTENT(in) :: kt 90 !!---------------------------------------------------------------------- 91 INTEGER :: ji, jj, jk ! dummy loop indices 92 !!---------------------------------------------------------------------- 102 93 103 94 DO ji = nie0, nie1 … … 120 111 121 112 SUBROUTINE obc_dyn_bt_west 122 !!---------------------------------------------------------------------- --------113 !!---------------------------------------------------------------------- 123 114 !! *** SUBROUTINE obc_dyn_bt_west *** 124 115 !! … … 127 118 !! Fix sea surface height (sshn) on west open boundary 128 119 !! The logical lfbcwest must be .TRUE. 129 !! 130 !! History : 131 !! 9.0 ! 05-12 (V. Garnier) original 132 !!------------------------------------------------------------------------------ 133 !! * Local declaration 134 INTEGER :: ji, jj, jk ! dummy loop indices 135 !!------------------------------------------------------------------------------ 136 120 !!---------------------------------------------------------------------- 121 INTEGER :: ji, jj, jk ! dummy loop indices 122 !!---------------------------------------------------------------------- 123 ! 137 124 DO ji = niw0, niw1 138 125 DO jk = 1, jpkm1 … … 147 134 END DO 148 135 END DO 149 136 ! 150 137 END SUBROUTINE obc_dyn_bt_west 138 151 139 152 140 SUBROUTINE obc_dyn_bt_north … … 158 146 !! Fix sea surface height (sshn) on north open boundary 159 147 !! The logical lfbcnorth must be .TRUE. 160 !! 161 !! History : 162 !! 9.0 ! 05-12 (V. Garnier) original 163 !!------------------------------------------------------------------------------ 164 !! * Local declaration 165 INTEGER :: ji, jj, jk ! dummy loop indices 166 !!------------------------------------------------------------------------------ 167 148 !!---------------------------------------------------------------------- 149 INTEGER :: ji, jj, jk ! dummy loop indices 150 !!---------------------------------------------------------------------- 151 ! 168 152 DO jj = njn0, njn1 169 153 DO jk = 1, jpkm1 … … 180 164 END DO 181 165 END DO 182 166 ! 183 167 END SUBROUTINE obc_dyn_bt_north 184 168 169 185 170 SUBROUTINE obc_dyn_bt_south 186 !!---------------------------------------------------------------------- --------171 !!---------------------------------------------------------------------- 187 172 !! *** SUBROUTINE obc_dyn_bt_south *** 188 173 !! … … 191 176 !! Fix sea surface height (sshn) on south open boundary 192 177 !! The logical lfbcsouth must be .TRUE. 193 !! 194 !! History : 195 !! 9.0 ! 05-12 (V. Garnier) original 196 !!------------------------------------------------------------------------------ 197 !! * Local declaration 198 INTEGER :: ji, jj, jk ! dummy loop indices 199 200 !!------------------------------------------------------------------------------ 201 178 !!---------------------------------------------------------------------- 179 INTEGER :: ji, jj, jk ! dummy loop indices 180 !!---------------------------------------------------------------------- 181 ! 202 182 DO jj = njs0, njs1 203 183 DO jk = 1, jpkm1 … … 212 192 END DO 213 193 END DO 214 194 ! 215 195 END SUBROUTINE obc_dyn_bt_south 216 196 … … 225 205 !! Fix sea surface height (sshn) on east open boundary 226 206 !! The logical lfbceast must be .TRUE. 227 !! 228 !! History : 229 !! 9.0 ! 05-12 (V. Garnier) original 230 !!------------------------------------------------------------------------------ 231 !! * Local declaration 232 INTEGER :: ji, jj, jk ! dummy loop indices 233 !!------------------------------------------------------------------------------ 234 207 !!---------------------------------------------------------------------- 208 INTEGER :: ji, jj, jk ! dummy loop indices 209 !!---------------------------------------------------------------------- 210 ! 235 211 DO ji = nie0, nie1 236 212 DO jk = 1, jpkm1 … … 245 221 END DO 246 222 END DO 247 223 ! 248 224 END SUBROUTINE obc_dyn_bt_east 249 225 226 250 227 SUBROUTINE obc_dyn_bt_west 251 !!--------------------------------------------------------------------- ---------228 !!--------------------------------------------------------------------- 252 229 !! *** SUBROUTINE obc_dyn_bt_west *** 253 230 !! 254 !! ** Purpose : 255 !! ** Purpose : 256 !! Apply Flather algorithm on west OBC velocities ua, va 231 !! ** Purpose : Apply Flather algorithm on west OBC velocities ua, va 257 232 !! Fix sea surface height (sshn) on west open boundary 258 233 !! The logical lfbcwest must be .TRUE. 259 !! 260 !! History : 261 !! 9.0 ! 05-12 (V. Garnier) original 262 !!------------------------------------------------------------------------------ 263 !! * Local declaration 264 INTEGER :: ji, jj, jk ! dummy loop indices 265 !!------------------------------------------------------------------------------ 266 234 !!---------------------------------------------------------------------- 235 INTEGER :: ji, jj, jk ! dummy loop indices 236 !!---------------------------------------------------------------------- 237 ! 267 238 DO ji = niw0, niw1 268 239 DO jk = 1, jpkm1 … … 275 246 END DO 276 247 END DO 277 248 ! 278 249 END SUBROUTINE obc_dyn_bt_west 250 279 251 280 252 SUBROUTINE obc_dyn_bt_north 281 253 !!------------------------------------------------------------------------------ 282 !! SUBROUTINE obc_dyn_bt_north283 !! *************************254 !! *** SUBROUTINE obc_dyn_bt_north *** 255 !! 284 256 !! ** Purpose : 285 257 !! Apply Flather algorithm on north OBC velocities ua, va 286 258 !! Fix sea surface height (sshn) on north open boundary 287 259 !! The logical lfbcnorth must be .TRUE. 288 !! 289 !! History : 290 !! 9.0 ! 05-12 (V. Garnier) original 291 !!------------------------------------------------------------------------------ 292 !! * Local declaration 293 INTEGER :: ji, jj, jk ! dummy loop indices 294 !!------------------------------------------------------------------------------ 295 260 !!---------------------------------------------------------------------- 261 INTEGER :: ji, jj, jk ! dummy loop indices 262 !!---------------------------------------------------------------------- 263 ! 296 264 DO jj = njn0, njn1 297 265 DO jk = 1, jpkm1 … … 306 274 END DO 307 275 END DO 308 276 ! 309 277 END SUBROUTINE obc_dyn_bt_north 278 310 279 311 280 SUBROUTINE obc_dyn_bt_south 312 281 !!------------------------------------------------------------------------------ 313 !! SUBROUTINE obc_dyn_bt_south314 !! *************************282 !! *** SUBROUTINE obc_dyn_bt_south *** 283 !! 315 284 !! ** Purpose : 316 285 !! Apply Flather algorithm on south OBC velocities ua, va 317 286 !! Fix sea surface height (sshn) on south open boundary 318 287 !! The logical lfbcsouth must be .TRUE. 319 !! 320 !! History : 321 !! 9.0 ! 05-12 (V. Garnier) original 322 !!------------------------------------------------------------------------------ 323 !! * Local declaration 324 INTEGER :: ji, jj, jk ! dummy loop indices 325 326 !!------------------------------------------------------------------------------ 327 288 !!---------------------------------------------------------------------- 289 INTEGER :: ji, jj, jk ! dummy loop indices 290 !!---------------------------------------------------------------------- 291 ! 328 292 DO jj = njs0, njs1 329 293 DO jk = 1, jpkm1 … … 336 300 END DO 337 301 END DO 338 302 ! 339 303 END SUBROUTINE obc_dyn_bt_south 340 304 341 305 # endif 306 342 307 #else 343 !!================================================================================= 344 !! *** MODULE obcdyn_bt *** 345 !! Ocean dynamics: Radiation of velocities on each open boundary 346 !!================================================================================= 308 !!---------------------------------------------------------------------- 309 !! Default option No Open Boundaries or not explicit fre surface 310 !!---------------------------------------------------------------------- 347 311 CONTAINS 348 349 SUBROUTINE obc_dyn_bt 350 ! No open boundaries ==> empty routine 312 SUBROUTINE obc_dyn_bt ! Dummy routine 351 313 END SUBROUTINE obc_dyn_bt 352 314 #endif 353 315 316 !!====================================================================== 354 317 END MODULE obcdyn_bt
Note: See TracChangeset
for help on using the changeset viewer.