- Timestamp:
- 2012-10-08T16:27:20+02:00 (12 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2012/dev_r3327_MERCATOR1_BDY/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90
r3367 r3490 11 11 !! 3.3 ! 2010-09 (D.Storkey) add ice boundary conditions 12 12 !! 3.4 ! 2011 (D. Storkey) rewrite in preparation for OBC-BDY merge 13 !! 3.4 ! 2012 (J. Chanut) straight open boundary case update 13 14 !!---------------------------------------------------------------------- 14 15 #if defined key_bdy … … 26 27 USE lib_mpp ! for mpp_sum 27 28 USE iom ! I/O 29 USE sbctide, ONLY: lk_tide ! Tidal forcing or not 30 USE phycst, ONLY: rday 28 31 29 32 IMPLICIT NONE … … 32 35 PUBLIC bdy_init ! routine called in nemo_init 33 36 37 INTEGER, PARAMETER :: jp_nseg = 100 38 INTEGER, PARAMETER :: nrimmax = 20 ! maximum rimwidth in structured 39 ! open boundary data files 40 ! Straight open boundary segment parameters: 41 INTEGER :: nbdysege, nbdysegw, nbdysegn, nbdysegs 42 INTEGER, DIMENSION(jp_nseg) :: jpieob, jpjedt, jpjeft, npckge 43 INTEGER, DIMENSION(jp_nseg) :: jpiwob, jpjwdt, jpjwft, npckgw 44 INTEGER, DIMENSION(jp_nseg) :: jpjnob, jpindt, jpinft, npckgn 45 INTEGER, DIMENSION(jp_nseg) :: jpjsob, jpisdt, jpisft, npckgs 34 46 !!---------------------------------------------------------------------- 35 47 !! NEMO/OPA 4.0 , NEMO Consortium (2011) … … 53 65 ! namelist variables 54 66 !------------------- 55 INTEGER, PARAMETER :: jp_nseg = 100 56 INTEGER :: nbdysege, nbdysegw, nbdysegn, nbdysegs 57 INTEGER, DIMENSION(jp_nseg) :: jpieob, jpjedt, jpjeft 58 INTEGER, DIMENSION(jp_nseg) :: jpiwob, jpjwdt, jpjwft 59 INTEGER, DIMENSION(jp_nseg) :: jpjnob, jpindt, jpinft 60 INTEGER, DIMENSION(jp_nseg) :: jpjsob, jpisdt, jpisft 67 CHARACTER(LEN=80),DIMENSION(jpbgrd) :: clfile 68 CHARACTER(LEN=1) :: ctypebdy 69 INTEGER :: nbdyind, nbdybeg, nbdyend 61 70 62 71 ! local variables … … 66 75 INTEGER :: iw, ie, is, in, inum, id_dummy ! - - 67 76 INTEGER :: igrd_start, igrd_end, jpbdta ! - - 77 INTEGER :: jpbdtau, jpbdtas ! - - 78 INTEGER :: ib_bdy1, ib_bdy2, ib1, ib2 ! - - 68 79 INTEGER, POINTER :: nbi, nbj, nbr ! short cuts 69 80 REAL , POINTER :: flagu, flagv ! - - 70 81 REAL(wp) :: zefl, zwfl, znfl, zsfl ! local scalars 71 INTEGER, DIMENSION (2) :: kdimsz82 INTEGER, DIMENSION (2) :: kdimsz 72 83 INTEGER, DIMENSION(jpbgrd,jp_bdy) :: nblendta ! Length of index arrays 73 84 INTEGER, ALLOCATABLE, DIMENSION(:,:,:) :: nbidta, nbjdta ! Index arrays: i and j indices of bdy dta 74 85 INTEGER, ALLOCATABLE, DIMENSION(:,:,:) :: nbrdta ! Discrete distance from rim points 75 REAL(wp), DIMENSION(jpidta,jpjdta) :: zmask ! global domain mask 76 CHARACTER(LEN=80),DIMENSION(jpbgrd) :: clfile 77 CHARACTER(LEN=1),DIMENSION(jpbgrd) :: cgrid 86 CHARACTER(LEN=1),DIMENSION(jpbgrd) :: cgrid 78 87 !! 79 88 NAMELIST/nambdy/ nb_bdy, ln_coords_file, cn_coords_file, & … … 86 95 & ln_vol, nn_volctl, nn_rimwidth 87 96 !! 88 NAMELIST/nambdy_index/ nbdysege, jpieob, jpjedt, jpjeft, & 89 nbdysegw, jpiwob, jpjwdt, jpjwft, & 90 nbdysegn, jpjnob, jpindt, jpinft, & 91 nbdysegs, jpjsob, jpisdt, jpisft 97 NAMELIST/nambdy_index/ ctypebdy, nbdyind, nbdybeg, nbdyend 92 98 93 99 !!---------------------------------------------------------------------- … … 106 112 107 113 cgrid= (/'t','u','v'/) 108 114 109 115 ! ----------------------------------------- 110 116 ! Initialise and read namelist parameters … … 139 145 ! Check and write out namelist parameters 140 146 ! ----------------------------------------- 141 142 ln_tra_dmp_tot=.false.143 ln_dyn3d_dmp_tot=.false.144 147 ! ! control prints 145 148 IF(lwp) WRITE(numout,*) ' nambdy' … … 164 167 IF(lwp) WRITE(numout,*) 'Boundary conditions for barotropic solution: ' 165 168 SELECT CASE( nn_dyn2d(ib_bdy) ) 166 CASE( 0 ); IF(lwp) WRITE(numout,*) ' no open boundary condition'167 CASE( 1 ); IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme'168 CASE( 2) ; IF(lwp) WRITE(numout,*) ' Flather radiation condition'169 CASE(jp_none) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 170 CASE(jp_frs) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 171 CASE(jp_flather) ; IF(lwp) WRITE(numout,*) ' Flather radiation condition' 169 172 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_dyn2d' ) 170 173 END SELECT … … 177 180 CASE DEFAULT ; CALL ctl_stop( 'nn_dyn2d_dta must be between 0 and 3' ) 178 181 END SELECT 182 IF (( nn_dyn2d_dta(ib_bdy) .ge. 2 ).AND.(.NOT.lk_tide)) THEN 183 CALL ctl_stop( 'You must activate key_tide to add tidal forcing at open boundaries' ) 184 ENDIF 179 185 ENDIF 180 186 IF(lwp) WRITE(numout,*) … … 182 188 IF(lwp) WRITE(numout,*) 'Boundary conditions for baroclinic velocities: ' 183 189 SELECT CASE( nn_dyn3d(ib_bdy) ) 184 CASE( 0 ); IF(lwp) WRITE(numout,*) ' no open boundary condition'185 CASE( 1 ); IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme'190 CASE(jp_none) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 191 CASE(jp_frs) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 186 192 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Specified value' 187 193 CASE( 3 ) ; IF(lwp) WRITE(numout,*) ' Zero baroclinic velocities (runoff case)' … … 195 201 END SELECT 196 202 ENDIF 197 IF(lwp) WRITE(numout,*)198 203 199 204 IF ( ln_dyn3d_dmp(ib_bdy) ) THEN 200 205 IF ( nn_dyn3d(ib_bdy).EQ.0 ) THEN 201 IF(lwp) WRITE(numout,*) 'No open boundary condition for the baroclinic velocitues : we force ln_dyn3d_dmp=.false.'206 IF(lwp) WRITE(numout,*) 'No open boundary condition for baroclinic velocities: ln_dyn3d_dmp is set to .false.' 202 207 ln_dyn3d_dmp(ib_bdy)=.false. 208 ELSEIF ( nn_dyn3d(ib_bdy).EQ.1 ) THEN 209 CALL ctl_stop( 'Use FRS OR relaxation' ) 203 210 ELSE 204 IF(lwp) WRITE(numout,*) ' Damping of the baroclinic velocities along the boundaries'205 IF(lwp) WRITE(numout,*) ' Time damping :',rn_time_dmp(ib_bdy),' days'206 ln_dyn3d_dmp_tot=.true.211 IF(lwp) WRITE(numout,*) ' + baroclinic velocities relaxation zone' 212 IF(lwp) WRITE(numout,*) ' Damping time scale: ',rn_time_dmp(ib_bdy),' days' 213 IF((lwp).AND.rn_time_dmp(ib_bdy)<0) CALL ctl_stop( 'Time scale must be positive' ) 207 214 ENDIF 215 ELSE 216 IF(lwp) WRITE(numout,*) ' NO relaxation on baroclinic velocities' 208 217 ENDIF 209 218 IF(lwp) WRITE(numout,*) … … 211 220 IF(lwp) WRITE(numout,*) 'Boundary conditions for temperature and salinity: ' 212 221 SELECT CASE( nn_tra(ib_bdy) ) 213 CASE( 0 ); IF(lwp) WRITE(numout,*) ' no open boundary condition'214 CASE( 1 ); IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme'222 CASE(jp_none) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 223 CASE(jp_frs) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 215 224 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Specified value' 216 225 CASE( 3 ) ; IF(lwp) WRITE(numout,*) ' Neumann conditions' … … 225 234 END SELECT 226 235 ENDIF 227 IF(lwp) WRITE(numout,*)228 236 229 237 IF ( ln_tra_dmp(ib_bdy) ) THEN 230 238 IF ( nn_tra(ib_bdy).EQ.0 ) THEN 231 IF(lwp) WRITE(numout,*) 'No open boundary condition for t he tracer : we force ln_tra_dmp=.false.'239 IF(lwp) WRITE(numout,*) 'No open boundary condition for tracers: ln_tra_dmp is set to .false.' 232 240 ln_tra_dmp(ib_bdy)=.false. 241 ELSEIF ( nn_tra(ib_bdy).EQ.1 ) THEN 242 CALL ctl_stop( 'Use FRS OR relaxation' ) 233 243 ELSE 234 IF(lwp) WRITE(numout,*) ' Damping of the scalars along the boundaries'235 IF(lwp) WRITE(numout,*) ' Time dumping :',rn_time_dmp(ib_bdy),' days'236 ln_tra_dmp_tot=.true.244 IF(lwp) WRITE(numout,*) ' + T/S relaxation zone' 245 IF(lwp) WRITE(numout,*) ' Damping time scale: ',rn_time_dmp(ib_bdy),' days' 246 IF((lwp).AND.rn_time_dmp(ib_bdy)<0) CALL ctl_stop( 'Time scale must be positive' ) 237 247 ENDIF 248 ELSE 249 IF(lwp) WRITE(numout,*) ' NO T/S relaxation' 238 250 ENDIF 239 251 IF(lwp) WRITE(numout,*) … … 256 268 #endif 257 269 258 IF(lwp) WRITE(numout,*) ' Boundary rim width for the FRS scheme = ', nn_rimwidth(ib_bdy)270 IF(lwp) WRITE(numout,*) ' Width of relaxation zone = ', nn_rimwidth(ib_bdy) 259 271 IF(lwp) WRITE(numout,*) 260 272 261 273 ENDDO 262 274 263 IF( ln_vol ) THEN ! check volume conservation (nn_volctl value) 264 IF(lwp) WRITE(numout,*) 'Volume correction applied at open boundaries' 265 IF(lwp) WRITE(numout,*) 266 SELECT CASE ( nn_volctl ) 267 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' The total volume will be constant' 268 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' The total volume will vary according to the surface E-P flux' 269 CASE DEFAULT ; CALL ctl_stop( 'nn_volctl must be 0 or 1' ) 270 END SELECT 271 IF(lwp) WRITE(numout,*) 272 ELSE 273 IF(lwp) WRITE(numout,*) 'No volume correction applied at open boundaries' 274 IF(lwp) WRITE(numout,*) 275 IF (nb_bdy .gt. 0) THEN 276 IF( ln_vol ) THEN ! check volume conservation (nn_volctl value) 277 IF(lwp) WRITE(numout,*) 'Volume correction applied at open boundaries' 278 IF(lwp) WRITE(numout,*) 279 SELECT CASE ( nn_volctl ) 280 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' The total volume will be constant' 281 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' The total volume will vary according to the surface E-P flux' 282 CASE DEFAULT ; CALL ctl_stop( 'nn_volctl must be 0 or 1' ) 283 END SELECT 284 IF(lwp) WRITE(numout,*) 285 ELSE 286 IF(lwp) WRITE(numout,*) 'No volume correction applied at open boundaries' 287 IF(lwp) WRITE(numout,*) 288 ENDIF 275 289 ENDIF 276 290 … … 282 296 ! --------------------------------------------- 283 297 REWIND( numnam ) 284 jpbdta = 1298 285 299 nblendta(:,:) = 0 300 nbdysege = 0 301 nbdysegw = 0 302 nbdysegn = 0 303 nbdysegs = 0 304 icount = 0 ! count user defined segments 305 ! Dimensions below are used to allocate arrays to read external data 306 jpbdtas = 1 ! Maximum size of boundary data (structured case) 307 jpbdtau = 1 ! Maximum size of boundary data (unstructured case) 308 286 309 DO ib_bdy = 1, nb_bdy 287 310 288 311 IF( .NOT. ln_coords_file(ib_bdy) ) THEN ! Work out size of global arrays from namelist parameters 289 312 313 icount = icount + 1 290 314 ! No REWIND here because may need to read more than one nambdy_index namelist. 291 315 READ ( numnam, nambdy_index ) 292 316 293 ! Automatic boundary definition: if nbdysegX = -1294 ! set boundary to whole side of model domain.295 IF( nbdysege == -1 ) THEN296 nbdysege = 1297 jpieob(1) = jpiglo - 1298 jpjedt(1) = 2299 jpjeft(1) = jpjglo - 1300 ENDIF301 IF( nbdysegw == -1 ) THEN302 nbdysegw = 1303 jpiwob(1) = 2304 jpjwdt(1) = 2305 jpjwft(1) = jpjglo - 1306 ENDIF307 IF( nbdysegn == -1 ) THEN308 nbdysegn = 1309 jpjnob(1) = jpjglo - 1310 jpindt(1) = 2311 jpinft(1) = jpiglo - 1312 ENDIF313 IF( nbdysegs == -1 ) THEN314 nbdysegs = 1315 jpjsob(1) = 2316 jpisdt(1) = 2317 jpisft(1) = jpiglo - 1318 ENDIF319 320 DO iseg = 1, nbdysege321 igrd = 1322 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpjeft(iseg) - jpjedt(iseg) + 1323 igrd = 2324 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpjeft(iseg) - jpjedt(iseg) + 1325 igrd = 3326 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpjeft(iseg) - jpjedt(iseg)327 ENDDO328 DO iseg = 1, nbdysegw329 igrd = 1330 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpjwft(iseg) - jpjwdt(iseg) + 1331 igrd = 2332 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpjwft(iseg) - jpjwdt(iseg) + 1333 igrd = 3334 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpjwft(iseg) - jpjwdt(iseg)335 ENDDO336 DO iseg = 1, nbdysegn337 igrd = 1338 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpinft(iseg) - jpindt(iseg) + 1339 igrd = 2340 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpinft(iseg) - jpindt(iseg)341 igrd = 3342 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpinft(iseg) - jpindt(iseg) + 1343 END DO344 DO iseg = 1, nbdysegs 345 igrd = 1346 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpisft(iseg) - jpisdt(iseg) + 1347 igrd = 2348 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpisft(iseg) - jpisdt(iseg)349 igrd = 3 350 nblendta(igrd,ib_bdy) = nblendta(igrd,ib_bdy) + jpisft(iseg) - jpisdt(iseg) + 1351 ENDDO352 353 nblendta(:,ib_bdy) = nblendta(:,ib_bdy) * nn_rimwidth(ib_bdy)317 SELECT CASE ( TRIM(ctypebdy) ) 318 CASE( 'N' ) 319 IF( nbdyind == -1 ) THEN ! Automatic boundary definition: if nbdysegX = -1 320 nbdyind = jpjglo - 2 ! set boundary to whole side of model domain. 321 nbdybeg = 2 322 nbdyend = jpiglo - 1 323 ENDIF 324 nbdysegn = nbdysegn + 1 325 npckgn(nbdysegn) = ib_bdy ! Save bdy package number 326 jpjnob(nbdysegn) = nbdyind 327 jpindt(nbdysegn) = nbdybeg 328 jpinft(nbdysegn) = nbdyend 329 ! 330 CASE( 'S' ) 331 IF( nbdyind == -1 ) THEN ! Automatic boundary definition: if nbdysegX = -1 332 nbdyind = 2 ! set boundary to whole side of model domain. 333 nbdybeg = 2 334 nbdyend = jpiglo - 1 335 ENDIF 336 nbdysegs = nbdysegs + 1 337 npckgs(nbdysegs) = ib_bdy ! Save bdy package number 338 jpjsob(nbdysegs) = nbdyind 339 jpisdt(nbdysegs) = nbdybeg 340 jpisft(nbdysegs) = nbdyend 341 ! 342 CASE( 'E' ) 343 IF( nbdyind == -1 ) THEN ! Automatic boundary definition: if nbdysegX = -1 344 nbdyind = jpiglo - 2 ! set boundary to whole side of model domain. 345 nbdybeg = 2 346 nbdyend = jpjglo - 1 347 ENDIF 348 nbdysege = nbdysege + 1 349 npckge(nbdysege) = ib_bdy ! Save bdy package number 350 jpieob(nbdysege) = nbdyind 351 jpjedt(nbdysege) = nbdybeg 352 jpjeft(nbdysege) = nbdyend 353 ! 354 CASE( 'W' ) 355 IF( nbdyind == -1 ) THEN ! Automatic boundary definition: if nbdysegX = -1 356 nbdyind = 2 ! set boundary to whole side of model domain. 357 nbdybeg = 2 358 nbdyend = jpjglo - 1 359 ENDIF 360 nbdysegw = nbdysegw + 1 361 npckgw(nbdysegw) = ib_bdy ! Save bdy package number 362 jpiwob(nbdysegw) = nbdyind 363 jpjwdt(nbdysegw) = nbdybeg 364 jpjwft(nbdysegw) = nbdyend 365 ! 366 CASE DEFAULT ; CALL ctl_stop( 'ctypebdy must be N, S, E or W' ) 367 END SELECT 368 369 ! For simplicity we assume that in case of straight bdy, arrays have the same length 370 ! (even if it is true that last tangential velocity points 371 ! are useless). This simplifies a little bit boundary data format (and agrees with format 372 ! used so far in obc package) 373 374 nblendta(1:jpbgrd,ib_bdy) = (nbdyend - nbdybeg + 1) * nn_rimwidth(ib_bdy) 375 jpbdtas = MAX(jpbdtas, (nbdyend - nbdybeg + 1)) 376 IF (lwp.and.(nn_rimwidth(ib_bdy)>nrimmax)) & 377 & CALL ctl_stop( 'rimwidth must be lower than nrimmax' ) 354 378 355 379 ELSE ! Read size of arrays in boundary coordinates file. 356 357 358 380 CALL iom_open( cn_coords_file(ib_bdy), inum ) 359 381 DO igrd = 1, jpbgrd 360 382 id_dummy = iom_varid( inum, 'nbi'//cgrid(igrd), kdimsz=kdimsz ) 361 383 nblendta(igrd,ib_bdy) = kdimsz(1) 362 ENDDO 384 jpbdtau = MAX(jpbdtau, kdimsz(1)) 385 ENDDO 386 CALL iom_close( inum ) 363 387 364 388 ENDIF … … 366 390 ENDDO ! ib_bdy 367 391 368 jpbdta = MAXVAL(nblendta(:,:)) 369 370 ! Allocate arrays 371 !--------------- 372 ALLOCATE( nbidta(jpbdta, jpbgrd, nb_bdy), nbjdta(jpbdta, jpbgrd, nb_bdy), & 373 & nbrdta(jpbdta, jpbgrd, nb_bdy) ) 374 375 ALLOCATE( dta_global(jpbdta, 1, jpk) ) 392 IF (nb_bdy>0) THEN 393 jpbdta = MAXVAL(nblendta(1:jpbgrd,1:nb_bdy)) 394 395 ! Allocate arrays 396 !--------------- 397 ALLOCATE( nbidta(jpbdta, jpbgrd, nb_bdy), nbjdta(jpbdta, jpbgrd, nb_bdy), & 398 & nbrdta(jpbdta, jpbgrd, nb_bdy) ) 399 400 ALLOCATE( dta_global(jpbdtau, 1, jpk) ) 401 IF ( icount>0 ) ALLOCATE( dta_global2(jpbdtas, nrimmax, jpk) ) 402 ! 403 ENDIF 404 405 ! Now look for crossings in user (namelist) defined open boundary segments: 406 !-------------------------------------------------------------------------- 407 IF ( icount>0 ) CALL bdy_ctl_seg 376 408 377 409 ! Calculate global boundary index arrays or read in from file 378 !------------------------------------------------------------ 379 REWIND( numnam )410 !------------------------------------------------------------ 411 ! 1. Read global index arrays from boundary coordinates file. 380 412 DO ib_bdy = 1, nb_bdy 381 413 382 IF( .NOT. ln_coords_file(ib_bdy) ) THEN ! Calculate global index arrays from namelist parameters 383 384 ! No REWIND here because may need to read more than one nambdy_index namelist. 385 READ ( numnam, nambdy_index ) 386 387 ! Automatic boundary definition: if nbdysegX = -1 388 ! set boundary to whole side of model domain. 389 IF( nbdysege == -1 ) THEN 390 nbdysege = 1 391 jpieob(1) = jpiglo - 1 392 jpjedt(1) = 2 393 jpjeft(1) = jpjglo - 1 394 ENDIF 395 IF( nbdysegw == -1 ) THEN 396 nbdysegw = 1 397 jpiwob(1) = 2 398 jpjwdt(1) = 2 399 jpjwft(1) = jpjglo - 1 400 ENDIF 401 IF( nbdysegn == -1 ) THEN 402 nbdysegn = 1 403 jpjnob(1) = jpjglo - 1 404 jpindt(1) = 2 405 jpinft(1) = jpiglo - 1 406 ENDIF 407 IF( nbdysegs == -1 ) THEN 408 nbdysegs = 1 409 jpjsob(1) = 2 410 jpisdt(1) = 2 411 jpisft(1) = jpiglo - 1 412 ENDIF 413 414 ! ------------ T points ------------- 415 igrd = 1 416 icount = 0 417 DO ir = 1, nn_rimwidth(ib_bdy) 418 ! east 419 DO iseg = 1, nbdysege 420 DO ij = jpjedt(iseg), jpjeft(iseg) 421 icount = icount + 1 422 nbidta(icount, igrd, ib_bdy) = jpieob(iseg) - ir + 1 423 nbjdta(icount, igrd, ib_bdy) = ij 424 nbrdta(icount, igrd, ib_bdy) = ir 425 ENDDO 426 ENDDO 427 ! west 428 DO iseg = 1, nbdysegw 429 DO ij = jpjwdt(iseg), jpjwft(iseg) 430 icount = icount + 1 431 nbidta(icount, igrd, ib_bdy) = jpiwob(iseg) + ir - 1 432 nbjdta(icount, igrd, ib_bdy) = ij 433 nbrdta(icount, igrd, ib_bdy) = ir 434 ENDDO 435 ENDDO 436 ! north 437 DO iseg = 1, nbdysegn 438 DO ii = jpindt(iseg), jpinft(iseg) 439 icount = icount + 1 440 nbidta(icount, igrd, ib_bdy) = ii 441 nbjdta(icount, igrd, ib_bdy) = jpjnob(iseg) - ir + 1 442 nbrdta(icount, igrd, ib_bdy) = ir 443 ENDDO 444 ENDDO 445 ! south 446 DO iseg = 1, nbdysegs 447 DO ii = jpisdt(iseg), jpisft(iseg) 448 icount = icount + 1 449 nbidta(icount, igrd, ib_bdy) = ii 450 nbjdta(icount, igrd, ib_bdy) = jpjsob(iseg) + ir - 1 451 nbrdta(icount, igrd, ib_bdy) = ir 452 ENDDO 453 ENDDO 454 ENDDO 455 456 ! ------------ U points ------------- 457 igrd = 2 458 icount = 0 459 DO ir = 1, nn_rimwidth(ib_bdy) 460 ! east 461 DO iseg = 1, nbdysege 462 DO ij = jpjedt(iseg), jpjeft(iseg) 463 icount = icount + 1 464 nbidta(icount, igrd, ib_bdy) = jpieob(iseg) - ir 465 nbjdta(icount, igrd, ib_bdy) = ij 466 nbrdta(icount, igrd, ib_bdy) = ir 467 ENDDO 468 ENDDO 469 ! west 470 DO iseg = 1, nbdysegw 471 DO ij = jpjwdt(iseg), jpjwft(iseg) 472 icount = icount + 1 473 nbidta(icount, igrd, ib_bdy) = jpiwob(iseg) + ir - 1 474 nbjdta(icount, igrd, ib_bdy) = ij 475 nbrdta(icount, igrd, ib_bdy) = ir 476 ENDDO 477 ENDDO 478 ! north 479 DO iseg = 1, nbdysegn 480 DO ii = jpindt(iseg), jpinft(iseg) - 1 481 icount = icount + 1 482 nbidta(icount, igrd, ib_bdy) = ii 483 nbjdta(icount, igrd, ib_bdy) = jpjnob(iseg) - ir + 1 484 nbrdta(icount, igrd, ib_bdy) = ir 485 ENDDO 486 ENDDO 487 ! south 488 DO iseg = 1, nbdysegs 489 DO ii = jpisdt(iseg), jpisft(iseg) - 1 490 icount = icount + 1 491 nbidta(icount, igrd, ib_bdy) = ii 492 nbjdta(icount, igrd, ib_bdy) = jpjsob(iseg) + ir - 1 493 nbrdta(icount, igrd, ib_bdy) = ir 494 ENDDO 495 ENDDO 496 ENDDO 497 498 ! ------------ V points ------------- 499 igrd = 3 500 icount = 0 501 DO ir = 1, nn_rimwidth(ib_bdy) 502 ! east 503 DO iseg = 1, nbdysege 504 DO ij = jpjedt(iseg), jpjeft(iseg) - 1 505 icount = icount + 1 506 nbidta(icount, igrd, ib_bdy) = jpieob(iseg) - ir + 1 507 nbjdta(icount, igrd, ib_bdy) = ij 508 nbrdta(icount, igrd, ib_bdy) = ir 509 ENDDO 510 ENDDO 511 ! west 512 DO iseg = 1, nbdysegw 513 DO ij = jpjwdt(iseg), jpjwft(iseg) - 1 514 icount = icount + 1 515 nbidta(icount, igrd, ib_bdy) = jpiwob(iseg) + ir - 1 516 nbjdta(icount, igrd, ib_bdy) = ij 517 nbrdta(icount, igrd, ib_bdy) = ir 518 ENDDO 519 ENDDO 520 ! north 521 DO iseg = 1, nbdysegn 522 DO ii = jpindt(iseg), jpinft(iseg) 523 icount = icount + 1 524 nbidta(icount, igrd, ib_bdy) = ii 525 nbjdta(icount, igrd, ib_bdy) = jpjnob(iseg) - ir 526 nbrdta(icount, igrd, ib_bdy) = ir 527 ENDDO 528 ENDDO 529 ! south 530 DO iseg = 1, nbdysegs 531 DO ii = jpisdt(iseg), jpisft(iseg) 532 icount = icount + 1 533 nbidta(icount, igrd, ib_bdy) = ii 534 nbjdta(icount, igrd, ib_bdy) = jpjsob(iseg) + ir - 1 535 nbrdta(icount, igrd, ib_bdy) = ir 536 ENDDO 537 ENDDO 538 ENDDO 539 540 ELSE ! Read global index arrays from boundary coordinates file. 541 414 IF( ln_coords_file(ib_bdy) ) THEN 415 416 CALL iom_open( cn_coords_file(ib_bdy), inum ) 542 417 DO igrd = 1, jpbgrd 543 418 CALL iom_get( inum, jpdom_unknown, 'nbi'//cgrid(igrd), dta_global(1:nblendta(igrd,ib_bdy),:,1) ) … … 560 435 IF (ibr_max < nn_rimwidth(ib_bdy)) & 561 436 CALL ctl_stop( 'nn_rimwidth is larger than maximum rimwidth in file',cn_coords_file(ib_bdy) ) 562 563 437 END DO 564 438 CALL iom_close( inum ) … … 566 440 ENDIF 567 441 568 ENDDO 442 ENDDO 443 444 ! 2. Now fill indices corresponding to straight open boundary arrays: 445 ! East 446 !----- 447 DO iseg = 1, nbdysege 448 ib_bdy = npckge(iseg) 449 ! 450 ! ------------ T points ------------- 451 igrd=1 452 icount=0 453 DO ir = 1, nn_rimwidth(ib_bdy) 454 DO ij = jpjedt(iseg), jpjeft(iseg) 455 icount = icount + 1 456 nbidta(icount, igrd, ib_bdy) = jpieob(iseg) + 2 - ir 457 nbjdta(icount, igrd, ib_bdy) = ij 458 nbrdta(icount, igrd, ib_bdy) = ir 459 ENDDO 460 ENDDO 461 ! 462 ! ------------ U points ------------- 463 igrd=2 464 icount=0 465 DO ir = 1, nn_rimwidth(ib_bdy) 466 DO ij = jpjedt(iseg), jpjeft(iseg) 467 icount = icount + 1 468 nbidta(icount, igrd, ib_bdy) = jpieob(iseg) + 1 - ir 469 nbjdta(icount, igrd, ib_bdy) = ij 470 nbrdta(icount, igrd, ib_bdy) = ir 471 ENDDO 472 ENDDO 473 ! 474 ! ------------ V points ------------- 475 igrd=3 476 icount=0 477 DO ir = 1, nn_rimwidth(ib_bdy) 478 ! DO ij = jpjedt(iseg), jpjeft(iseg) - 1 479 DO ij = jpjedt(iseg), jpjeft(iseg) 480 icount = icount + 1 481 nbidta(icount, igrd, ib_bdy) = jpieob(iseg) + 2 - ir 482 nbjdta(icount, igrd, ib_bdy) = ij 483 nbrdta(icount, igrd, ib_bdy) = ir 484 ENDDO 485 nbidta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 486 nbjdta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 487 ENDDO 488 ENDDO 489 ! 490 ! West 491 !----- 492 DO iseg = 1, nbdysegw 493 ib_bdy = npckgw(iseg) 494 ! 495 ! ------------ T points ------------- 496 igrd=1 497 icount=0 498 DO ir = 1, nn_rimwidth(ib_bdy) 499 DO ij = jpjwdt(iseg), jpjwft(iseg) 500 icount = icount + 1 501 nbidta(icount, igrd, ib_bdy) = jpiwob(iseg) + ir - 1 502 nbjdta(icount, igrd, ib_bdy) = ij 503 nbrdta(icount, igrd, ib_bdy) = ir 504 ENDDO 505 ENDDO 506 ! 507 ! ------------ U points ------------- 508 igrd=2 509 icount=0 510 DO ir = 1, nn_rimwidth(ib_bdy) 511 DO ij = jpjwdt(iseg), jpjwft(iseg) 512 icount = icount + 1 513 nbidta(icount, igrd, ib_bdy) = jpiwob(iseg) + ir - 1 514 nbjdta(icount, igrd, ib_bdy) = ij 515 nbrdta(icount, igrd, ib_bdy) = ir 516 ENDDO 517 ENDDO 518 ! 519 ! ------------ V points ------------- 520 igrd=3 521 icount=0 522 DO ir = 1, nn_rimwidth(ib_bdy) 523 ! DO ij = jpjwdt(iseg), jpjwft(iseg) - 1 524 DO ij = jpjwdt(iseg), jpjwft(iseg) 525 icount = icount + 1 526 nbidta(icount, igrd, ib_bdy) = jpiwob(iseg) + ir - 1 527 nbjdta(icount, igrd, ib_bdy) = ij 528 nbrdta(icount, igrd, ib_bdy) = ir 529 ENDDO 530 nbidta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 531 nbjdta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 532 ENDDO 533 ENDDO 534 ! 535 ! North 536 !----- 537 DO iseg = 1, nbdysegn 538 ib_bdy = npckgn(iseg) 539 ! 540 ! ------------ T points ------------- 541 igrd=1 542 icount=0 543 DO ir = 1, nn_rimwidth(ib_bdy) 544 DO ii = jpindt(iseg), jpinft(iseg) 545 icount = icount + 1 546 nbidta(icount, igrd, ib_bdy) = ii 547 nbjdta(icount, igrd, ib_bdy) = jpjnob(iseg) + 2 - ir 548 nbrdta(icount, igrd, ib_bdy) = ir 549 ENDDO 550 ENDDO 551 ! 552 ! ------------ U points ------------- 553 igrd=2 554 icount=0 555 DO ir = 1, nn_rimwidth(ib_bdy) 556 ! DO ii = jpindt(iseg), jpinft(iseg) - 1 557 DO ii = jpindt(iseg), jpinft(iseg) 558 icount = icount + 1 559 nbidta(icount, igrd, ib_bdy) = ii 560 nbjdta(icount, igrd, ib_bdy) = jpjnob(iseg) + 2 - ir 561 nbrdta(icount, igrd, ib_bdy) = ir 562 ENDDO 563 nbidta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 564 nbjdta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 565 ENDDO 566 ! 567 ! ------------ V points ------------- 568 igrd=3 569 icount=0 570 DO ir = 1, nn_rimwidth(ib_bdy) 571 DO ii = jpindt(iseg), jpinft(iseg) 572 icount = icount + 1 573 nbidta(icount, igrd, ib_bdy) = ii 574 nbjdta(icount, igrd, ib_bdy) = jpjnob(iseg) + 1 - ir 575 nbrdta(icount, igrd, ib_bdy) = ir 576 ENDDO 577 ENDDO 578 ENDDO 579 ! 580 ! South 581 !----- 582 DO iseg = 1, nbdysegs 583 ib_bdy = npckgs(iseg) 584 ! 585 ! ------------ T points ------------- 586 igrd=1 587 icount=0 588 DO ir = 1, nn_rimwidth(ib_bdy) 589 DO ii = jpisdt(iseg), jpisft(iseg) 590 icount = icount + 1 591 nbidta(icount, igrd, ib_bdy) = ii 592 nbjdta(icount, igrd, ib_bdy) = jpjsob(iseg) + ir - 1 593 nbrdta(icount, igrd, ib_bdy) = ir 594 ENDDO 595 ENDDO 596 ! 597 ! ------------ U points ------------- 598 igrd=2 599 icount=0 600 DO ir = 1, nn_rimwidth(ib_bdy) 601 ! DO ii = jpisdt(iseg), jpisft(iseg) - 1 602 DO ii = jpisdt(iseg), jpisft(iseg) 603 icount = icount + 1 604 nbidta(icount, igrd, ib_bdy) = ii 605 nbjdta(icount, igrd, ib_bdy) = jpjsob(iseg) + ir - 1 606 nbrdta(icount, igrd, ib_bdy) = ir 607 ENDDO 608 nbidta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 609 nbjdta(icount, igrd, ib_bdy) = -ib_bdy ! Discount this point 610 ENDDO 611 ! 612 ! ------------ V points ------------- 613 igrd=3 614 icount=0 615 DO ir = 1, nn_rimwidth(ib_bdy) 616 DO ii = jpisdt(iseg), jpisft(iseg) 617 icount = icount + 1 618 nbidta(icount, igrd, ib_bdy) = ii 619 nbjdta(icount, igrd, ib_bdy) = jpjsob(iseg) + ir - 1 620 nbrdta(icount, igrd, ib_bdy) = ir 621 ENDDO 622 ENDDO 623 ENDDO 624 625 ! Deal with duplicated points 626 !----------------------------- 627 ! We assign negative indices to duplicated points (to remove them from bdy points to be updated) 628 ! if their distance to the bdy is greater than the other 629 ! If their distance are the same, just keep only one to avoid updating a point twice 630 DO igrd = 1, jpbgrd 631 DO ib_bdy1 = 1, nb_bdy 632 DO ib_bdy2 = 1, nb_bdy 633 IF (ib_bdy1/=ib_bdy2) THEN 634 DO ib1 = 1, nblendta(igrd,ib_bdy1) 635 DO ib2 = 1, nblendta(igrd,ib_bdy2) 636 IF ((nbidta(ib1, igrd, ib_bdy1)==nbidta(ib2, igrd, ib_bdy2)).AND. & 637 & (nbjdta(ib1, igrd, ib_bdy1)==nbjdta(ib2, igrd, ib_bdy2))) THEN 638 ! IF ((lwp).AND.(igrd==1)) WRITE(numout,*) ' found coincident point ji, jj:', & 639 ! & nbidta(ib1, igrd, ib_bdy1), & 640 ! & nbjdta(ib2, igrd, ib_bdy2) 641 ! keep only points with the lowest distance to boundary: 642 IF (nbrdta(ib1, igrd, ib_bdy1)<nbrdta(ib2, igrd, ib_bdy2)) THEN 643 nbidta(ib2, igrd, ib_bdy2) =-ib_bdy2 644 nbjdta(ib2, igrd, ib_bdy2) =-ib_bdy2 645 ELSEIF (nbrdta(ib1, igrd, ib_bdy1)>nbrdta(ib2, igrd, ib_bdy2)) THEN 646 nbidta(ib1, igrd, ib_bdy1) =-ib_bdy1 647 nbjdta(ib1, igrd, ib_bdy1) =-ib_bdy1 648 ! Arbitrary choice if distances are the same: 649 ELSE 650 nbidta(ib1, igrd, ib_bdy1) =-ib_bdy1 651 nbjdta(ib1, igrd, ib_bdy1) =-ib_bdy1 652 ENDIF 653 END IF 654 END DO 655 END DO 656 ENDIF 657 END DO 658 END DO 659 END DO 569 660 570 661 ! Work out dimensions of boundary data on each processor 571 662 ! ------------------------------------------------------ 572 573 iw = mig(1) + 1 ! if monotasking and no zoom, iw=2 574 ie = mig(1) + nlci-1 - 1 ! if monotasking and no zoom, ie=jpim1 575 is = mjg(1) + 1 ! if monotasking and no zoom, is=2 576 in = mjg(1) + nlcj-1 - 1 ! if monotasking and no zoom, in=jpjm1 663 664 ! Rather assume that boundary data indices are given on global domain 665 ! TO BE DISCUSSED ? 666 ! iw = mig(1) + 1 ! if monotasking and no zoom, iw=2 667 ! ie = mig(1) + nlci-1 - 1 ! if monotasking and no zoom, ie=jpim1 668 ! is = mjg(1) + 1 ! if monotasking and no zoom, is=2 669 ! in = mjg(1) + nlcj-1 - 1 ! if monotasking and no zoom, in=jpjm1 670 iw = mig(1) - jpizoom + 2 ! if monotasking and no zoom, iw=2 671 ie = mig(1) + nlci - jpizoom - 1 ! if monotasking and no zoom, ie=jpim1 672 is = mjg(1) - jpjzoom + 2 ! if monotasking and no zoom, is=2 673 in = mjg(1) + nlcj - jpjzoom - 1 ! if monotasking and no zoom, in=jpjm1 577 674 578 675 DO ib_bdy = 1, nb_bdy … … 610 707 ALLOCATE( idx_bdy(ib_bdy)%nbj(ilen1,jpbgrd) ) 611 708 ALLOCATE( idx_bdy(ib_bdy)%nbr(ilen1,jpbgrd) ) 709 ALLOCATE( idx_bdy(ib_bdy)%nbd(ilen1,jpbgrd) ) 612 710 ALLOCATE( idx_bdy(ib_bdy)%nbmap(ilen1,jpbgrd) ) 613 711 ALLOCATE( idx_bdy(ib_bdy)%nbw(ilen1,jpbgrd) ) … … 629 727 ! 630 728 icount = icount + 1 631 idx_bdy(ib_bdy)%nbi(icount,igrd) = nbidta(ib,igrd,ib_bdy)- mig(1)+1 632 idx_bdy(ib_bdy)%nbj(icount,igrd) = nbjdta(ib,igrd,ib_bdy)- mjg(1)+1 729 730 ! Rather assume that boundary data indices are given on global domain 731 ! TO BE DISCUSSED ? 732 ! idx_bdy(ib_bdy)%nbi(icount,igrd) = nbidta(ib,igrd,ib_bdy)- mig(1)+1 733 ! idx_bdy(ib_bdy)%nbj(icount,igrd) = nbjdta(ib,igrd,ib_bdy)- mjg(1)+1 734 idx_bdy(ib_bdy)%nbi(icount,igrd) = nbidta(ib,igrd,ib_bdy)- mig(1)+jpizoom 735 idx_bdy(ib_bdy)%nbj(icount,igrd) = nbjdta(ib,igrd,ib_bdy)- mjg(1)+jpjzoom 633 736 idx_bdy(ib_bdy)%nbr(icount,igrd) = nbrdta(ib,igrd,ib_bdy) 634 737 idx_bdy(ib_bdy)%nbmap(icount,igrd) = ib … … 643 746 DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 644 747 nbr => idx_bdy(ib_bdy)%nbr(ib,igrd) 645 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( FLOAT( nbr - 1 ) *0.5 ) ! tanh formulation 646 idx_bdy(ib_bdy)%nbw(ib,igrd) = (FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2. ! quadratic 647 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = FLOAT(nn_rimwidth+1-nbr)/FLOAT(nn_rimwidth) ! linear 748 idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( FLOAT( nbr - 1 ) *0.5 ) ! tanh formulation 749 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = (FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2. ! quadratic 750 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)) ! linear 751 END DO 752 END DO 753 754 ! Compute damping coefficients 755 ! ---------------------------- 756 DO igrd = 1, jpbgrd 757 DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 758 nbr => idx_bdy(ib_bdy)%nbr(ib,igrd) 759 idx_bdy(ib_bdy)%nbd(ib,igrd) = 1. / ( rn_time_dmp(ib_bdy) * rday ) & 760 & *(FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2. ! quadratic 648 761 END DO 649 762 END DO … … 660 773 ! = 0 elsewhere 661 774 662 IF( cp_cfg == "eel" .AND. jp_cfg == 5 ) THEN ! EEL configuration at 5km resolution 663 zmask( : ,:) = 0.e0 664 zmask(jpizoom+1:jpizoom+jpiglo-2,:) = 1.e0 665 ELSE IF( ln_mask_file ) THEN 775 IF( ln_mask_file ) THEN 666 776 CALL iom_open( cn_mask_file, inum ) 667 CALL iom_get ( inum, jpdom_data, 'bdy_msk', zmask(:,:) )777 CALL iom_get ( inum, jpdom_data, 'bdy_msk', bdytmask(:,:) ) 668 778 CALL iom_close( inum ) 669 ELSE 670 zmask(:,:) = 1.e0 671 ENDIF 672 673 DO ij = 1, nlcj ! Save mask over local domain 674 DO ii = 1, nlci 675 bdytmask(ii,ij) = zmask( mig(ii), mjg(ij) ) 779 780 ! Derive mask on U and V grid from mask on T grid 781 bdyumask(:,:) = 0.e0 782 bdyvmask(:,:) = 0.e0 783 DO ij=1, jpjm1 784 DO ii=1, jpim1 785 bdyumask(ii,ij)=bdytmask(ii,ij)*bdytmask(ii+1, ij ) 786 bdyvmask(ii,ij)=bdytmask(ii,ij)*bdytmask(ii ,ij+1) 787 END DO 676 788 END DO 677 END DO 678 679 ! Derive mask on U and V grid from mask on T grid 680 bdyumask(:,:) = 0.e0 681 bdyvmask(:,:) = 0.e0 682 DO ij=1, jpjm1 683 DO ii=1, jpim1 684 bdyumask(ii,ij)=bdytmask(ii,ij)*bdytmask(ii+1, ij ) 685 bdyvmask(ii,ij)=bdytmask(ii,ij)*bdytmask(ii ,ij+1) 789 CALL lbc_lnk( bdyumask(:,:), 'U', 1. ) ; CALL lbc_lnk( bdyvmask(:,:), 'V', 1. ) ! Lateral boundary cond. 790 791 792 ! Mask corrections 793 ! ---------------- 794 DO ik = 1, jpkm1 795 DO ij = 1, jpj 796 DO ii = 1, jpi 797 tmask(ii,ij,ik) = tmask(ii,ij,ik) * bdytmask(ii,ij) 798 umask(ii,ij,ik) = umask(ii,ij,ik) * bdyumask(ii,ij) 799 vmask(ii,ij,ik) = vmask(ii,ij,ik) * bdyvmask(ii,ij) 800 bmask(ii,ij) = bmask(ii,ij) * bdytmask(ii,ij) 801 END DO 802 END DO 686 803 END DO 687 END DO 688 CALL lbc_lnk( bdyumask(:,:), 'U', 1. ) ; CALL lbc_lnk( bdyvmask(:,:), 'V', 1. ) ! Lateral boundary cond. 689 690 691 ! Mask corrections 692 ! ---------------- 693 DO ik = 1, jpkm1 694 DO ij = 1, jpj 695 DO ii = 1, jpi 696 tmask(ii,ij,ik) = tmask(ii,ij,ik) * bdytmask(ii,ij) 697 umask(ii,ij,ik) = umask(ii,ij,ik) * bdyumask(ii,ij) 698 vmask(ii,ij,ik) = vmask(ii,ij,ik) * bdyvmask(ii,ij) 699 bmask(ii,ij) = bmask(ii,ij) * bdytmask(ii,ij) 700 END DO 804 805 DO ik = 1, jpkm1 806 DO ij = 2, jpjm1 807 DO ii = 2, jpim1 808 fmask(ii,ij,ik) = fmask(ii,ij,ik) * bdytmask(ii,ij ) * bdytmask(ii+1,ij ) & 809 & * bdytmask(ii,ij+1) * bdytmask(ii+1,ij+1) 810 END DO 811 END DO 701 812 END DO 702 END DO 703 704 DO ik = 1, jpkm1 705 DO ij = 2, jpjm1 706 DO ii = 2, jpim1 707 fmask(ii,ij,ik) = fmask(ii,ij,ik) * bdytmask(ii,ij ) * bdytmask(ii+1,ij ) & 708 & * bdytmask(ii,ij+1) * bdytmask(ii+1,ij+1) 709 END DO 710 END DO 711 END DO 712 713 tmask_i (:,:) = tmask(:,:,1) * tmask_i(:,:) 813 814 tmask_i (:,:) = tmask(:,:,1) * tmask_i(:,:) 815 816 ENDIF ! ln_mask_file=.TRUE. 817 714 818 bdytmask(:,:) = tmask(:,:,1) 715 819 … … 800 904 ! Compute total lateral surface for volume correction: 801 905 ! ---------------------------------------------------- 906 ! JC: this must be done at each time step with key_vvl 802 907 bdysurftot = 0.e0 803 908 IF( ln_vol ) THEN … … 806 911 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 807 912 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 808 nbj => idx_bdy(ib_bdy)%nb i(ib,igrd)913 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 809 914 flagu => idx_bdy(ib_bdy)%flagu(ib) 810 915 bdysurftot = bdysurftot + hu (nbi , nbj) & … … 819 924 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 820 925 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 821 nbj => idx_bdy(ib_bdy)%nb i(ib,igrd)926 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 822 927 flagv => idx_bdy(ib_bdy)%flagv(ib) 823 928 bdysurftot = bdysurftot + hv (nbi, nbj ) & … … 833 938 ! Tidy up 834 939 !-------- 835 DEALLOCATE(nbidta, nbjdta, nbrdta) 940 IF (nb_bdy>0) THEN 941 DEALLOCATE(nbidta, nbjdta, nbrdta) 942 ENDIF 836 943 837 944 IF( nn_timing == 1 ) CALL timing_stop('bdy_init') 838 945 839 946 END SUBROUTINE bdy_init 947 948 SUBROUTINE bdy_ctl_seg 949 !!---------------------------------------------------------------------- 950 !! *** ROUTINE bdy_ctl_seg *** 951 !! 952 !! ** Purpose : Check straight open boundary segments location 953 !! 954 !! ** Method : - Look for open boundary corners 955 !! - Check that segments start or end on land 956 !!---------------------------------------------------------------------- 957 INTEGER :: ib, ib1, ib2, ji ,jj, itest 958 INTEGER, DIMENSION(jp_nseg,2) :: icorne, icornw, icornn, icorns 959 REAL(wp), DIMENSION(2) :: ztestmask 960 !!---------------------------------------------------------------------- 961 ! 962 IF (lwp) WRITE(numout,*) ' ' 963 IF (lwp) WRITE(numout,*) 'bdy_ctl_seg: Check analytical segments' 964 IF (lwp) WRITE(numout,*) '~~~~~~~~~~~~' 965 ! 966 IF(lwp) WRITE(numout,*) 'Number of east segments : ', nbdysege 967 IF(lwp) WRITE(numout,*) 'Number of west segments : ', nbdysegw 968 IF(lwp) WRITE(numout,*) 'Number of north segments : ', nbdysegn 969 IF(lwp) WRITE(numout,*) 'Number of south segments : ', nbdysegs 970 ! 1. Check bounds 971 !---------------- 972 DO ib = 1, nbdysegn 973 IF (lwp) WRITE(numout,*) '**check north seg bounds pckg: ', npckgn(ib) 974 IF ((jpjnob(ib).ge.jpjglo-1).or.& 975 &(jpjnob(ib).le.1)) CALL ctl_stop( 'nbdyind out of domain' ) 976 IF (jpindt(ib).ge.jpinft(ib)) CALL ctl_stop( 'Bdy start index is greater than end index' ) 977 IF (jpindt(ib).le.1 ) CALL ctl_stop( 'Start index out of domain' ) 978 IF (jpinft(ib).ge.jpiglo) CALL ctl_stop( 'End index out of domain' ) 979 END DO 980 ! 981 DO ib = 1, nbdysegs 982 IF (lwp) WRITE(numout,*) '**check south seg bounds pckg: ', npckgs(ib) 983 IF ((jpjsob(ib).ge.jpjglo-1).or.& 984 &(jpjsob(ib).le.1)) CALL ctl_stop( 'nbdyind out of domain' ) 985 IF (jpisdt(ib).ge.jpisft(ib)) CALL ctl_stop( 'Bdy start index is greater than end index' ) 986 IF (jpisdt(ib).le.1 ) CALL ctl_stop( 'Start index out of domain' ) 987 IF (jpisft(ib).ge.jpiglo) CALL ctl_stop( 'End index out of domain' ) 988 END DO 989 ! 990 DO ib = 1, nbdysege 991 IF (lwp) WRITE(numout,*) '**check east seg bounds pckg: ', npckge(ib) 992 IF ((jpieob(ib).ge.jpiglo-1).or.& 993 &(jpieob(ib).le.1)) CALL ctl_stop( 'nbdyind out of domain' ) 994 IF (jpjedt(ib).ge.jpjeft(ib)) CALL ctl_stop( 'Bdy start index is greater than end index' ) 995 IF (jpjedt(ib).le.1 ) CALL ctl_stop( 'Start index out of domain' ) 996 IF (jpjeft(ib).ge.jpjglo) CALL ctl_stop( 'End index out of domain' ) 997 END DO 998 ! 999 DO ib = 1, nbdysegw 1000 IF (lwp) WRITE(numout,*) '**check west seg bounds pckg: ', npckgw(ib) 1001 IF ((jpiwob(ib).ge.jpiglo-1).or.& 1002 &(jpiwob(ib).le.1)) CALL ctl_stop( 'nbdyind out of domain' ) 1003 IF (jpjwdt(ib).ge.jpjwft(ib)) CALL ctl_stop( 'Bdy start index is greater than end index' ) 1004 IF (jpjwdt(ib).le.1 ) CALL ctl_stop( 'Start index out of domain' ) 1005 IF (jpjwft(ib).ge.jpjglo) CALL ctl_stop( 'End index out of domain' ) 1006 ENDDO 1007 ! 1008 ! 1009 ! 2. Look for segment crossings 1010 !------------------------------ 1011 IF (lwp) WRITE(numout,*) '**Look for segments corners :' 1012 ! 1013 itest = 0 ! corner number 1014 ! 1015 ! flag to detect if start or end of open boundary belongs to a corner 1016 ! if not (=0), it must be on land. 1017 ! if a corner is detected, save bdy package number for further tests 1018 icorne(:,:)=0. ; icornw(:,:)=0. ; icornn(:,:)=0. ; icorns(:,:)=0. 1019 ! South/West crossings 1020 IF ((nbdysegw > 0).AND.(nbdysegs > 0)) THEN 1021 DO ib1 = 1, nbdysegw 1022 DO ib2 = 1, nbdysegs 1023 IF (( jpisdt(ib2)<=jpiwob(ib1)).AND. & 1024 & ( jpisft(ib2)>=jpiwob(ib1)).AND. & 1025 & ( jpjwdt(ib1)<=jpjsob(ib2)).AND. & 1026 & ( jpjwft(ib1)>=jpjsob(ib2))) THEN 1027 IF ((jpjwdt(ib1)==jpjsob(ib2)).AND.(jpisdt(ib2)==jpiwob(ib1))) THEN 1028 ! We have a possible South-West corner 1029 ! WRITE(numout,*) ' Found a South-West corner at (i,j): ', jpisdt(ib2), jpjwdt(ib1) 1030 ! WRITE(numout,*) ' between segments: ', npckgw(ib1), npckgs(ib2) 1031 icornw(ib1,1) = npckgs(ib2) 1032 icorns(ib2,1) = npckgw(ib1) 1033 ELSEIF ((jpisft(ib2)==jpiwob(ib1)).AND.(jpjwft(ib1)==jpjsob(ib2))) THEN 1034 IF(lwp) WRITE(numout,*) 1035 IF(lwp) WRITE(numout,*) ' E R R O R : Found an acute open boundary corner at point (i,j)= ', & 1036 & jpisft(ib2), jpjwft(ib1) 1037 IF(lwp) WRITE(numout,*) ' ========== Not allowed yet' 1038 IF(lwp) WRITE(numout,*) ' Crossing problem with West segment: ',npckgw(ib1), & 1039 & ' and South segment: ',npckgs(ib2) 1040 IF(lwp) WRITE(numout,*) 1041 nstop = nstop + 1 1042 ELSE 1043 IF(lwp) WRITE(numout,*) 1044 IF(lwp) WRITE(numout,*) ' E R R O R : Check South and West Open boundary indices' 1045 IF(lwp) WRITE(numout,*) ' ========== Crossing problem with West segment: ',npckgw(ib1) , & 1046 & ' and South segment: ',npckgs(ib2) 1047 IF(lwp) WRITE(numout,*) 1048 nstop = nstop+1 1049 END IF 1050 END IF 1051 END DO 1052 END DO 1053 END IF 1054 ! 1055 ! South/East crossings 1056 IF ((nbdysege > 0).AND.(nbdysegs > 0)) THEN 1057 DO ib1 = 1, nbdysege 1058 DO ib2 = 1, nbdysegs 1059 IF (( jpisdt(ib2)<=jpieob(ib1)+1).AND. & 1060 & ( jpisft(ib2)>=jpieob(ib1)+1).AND. & 1061 & ( jpjedt(ib1)<=jpjsob(ib2) ).AND. & 1062 & ( jpjeft(ib1)>=jpjsob(ib2) )) THEN 1063 IF ((jpjedt(ib1)==jpjsob(ib2)).AND.(jpisft(ib2)==jpieob(ib1)+1)) THEN 1064 ! We have a possible South-East corner 1065 ! WRITE(numout,*) ' Found a South-East corner at (i,j): ', jpisft(ib2), jpjedt(ib1) 1066 ! WRITE(numout,*) ' between segments: ', npckge(ib1), npckgs(ib2) 1067 icorne(ib1,1) = npckgs(ib2) 1068 icorns(ib2,2) = npckge(ib1) 1069 ELSEIF ((jpjeft(ib1)==jpjsob(ib2)).AND.(jpisdt(ib2)==jpieob(ib1)+1)) THEN 1070 IF(lwp) WRITE(numout,*) 1071 IF(lwp) WRITE(numout,*) ' E R R O R : Found an acute open boundary corner at point (i,j)= ', & 1072 & jpisdt(ib2), jpjeft(ib1) 1073 IF(lwp) WRITE(numout,*) ' ========== Not allowed yet' 1074 IF(lwp) WRITE(numout,*) ' Crossing problem with East segment: ',npckge(ib1), & 1075 & ' and South segment: ',npckgs(ib2) 1076 IF(lwp) WRITE(numout,*) 1077 nstop = nstop + 1 1078 ELSE 1079 IF(lwp) WRITE(numout,*) 1080 IF(lwp) WRITE(numout,*) ' E R R O R : Check South and East Open boundary indices' 1081 IF(lwp) WRITE(numout,*) ' ========== Crossing problem with East segment: ',npckge(ib1), & 1082 & ' and South segment: ',npckgs(ib2) 1083 IF(lwp) WRITE(numout,*) 1084 nstop = nstop + 1 1085 END IF 1086 END IF 1087 END DO 1088 END DO 1089 END IF 1090 ! 1091 ! North/West crossings 1092 IF ((nbdysegn > 0).AND.(nbdysegw > 0)) THEN 1093 DO ib1 = 1, nbdysegw 1094 DO ib2 = 1, nbdysegn 1095 IF (( jpindt(ib2)<=jpiwob(ib1) ).AND. & 1096 & ( jpinft(ib2)>=jpiwob(ib1) ).AND. & 1097 & ( jpjwdt(ib1)<=jpjnob(ib2)+1).AND. & 1098 & ( jpjwft(ib1)>=jpjnob(ib2)+1)) THEN 1099 IF ((jpjwft(ib1)==jpjnob(ib2)+1).AND.(jpindt(ib2)==jpiwob(ib1))) THEN 1100 ! We have a possible North-West corner 1101 ! WRITE(numout,*) ' Found a North-West corner at (i,j): ', jpindt(ib2), jpjwft(ib1) 1102 ! WRITE(numout,*) ' between segments: ', npckgw(ib1), npckgn(ib2) 1103 icornw(ib1,2) = npckgn(ib2) 1104 icornn(ib2,1) = npckgw(ib1) 1105 ELSEIF ((jpjwdt(ib1)==jpjnob(ib2)+1).AND.(jpinft(ib2)==jpiwob(ib1))) THEN 1106 IF(lwp) WRITE(numout,*) 1107 IF(lwp) WRITE(numout,*) ' E R R O R : Found an acute open boundary corner at point (i,j)= ', & 1108 & jpinft(ib2), jpjwdt(ib1) 1109 IF(lwp) WRITE(numout,*) ' ========== Not allowed yet' 1110 IF(lwp) WRITE(numout,*) ' Crossing problem with West segment: ',npckgw(ib1), & 1111 & ' and North segment: ',npckgn(ib2) 1112 IF(lwp) WRITE(numout,*) 1113 nstop = nstop + 1 1114 ELSE 1115 IF(lwp) WRITE(numout,*) 1116 IF(lwp) WRITE(numout,*) ' E R R O R : Check North and West Open boundary indices' 1117 IF(lwp) WRITE(numout,*) ' ========== Crossing problem with West segment: ',npckgw(ib1), & 1118 & ' and North segment: ',npckgn(ib2) 1119 IF(lwp) WRITE(numout,*) 1120 nstop = nstop + 1 1121 END IF 1122 END IF 1123 END DO 1124 END DO 1125 END IF 1126 ! 1127 ! North/East crossings 1128 IF ((nbdysegn > 0).AND.(nbdysege > 0)) THEN 1129 DO ib1 = 1, nbdysege 1130 DO ib2 = 1, nbdysegn 1131 IF (( jpindt(ib2)<=jpieob(ib1)+1).AND. & 1132 & ( jpinft(ib2)>=jpieob(ib1)+1).AND. & 1133 & ( jpjedt(ib1)<=jpjnob(ib2)+1).AND. & 1134 & ( jpjeft(ib1)>=jpjnob(ib2)+1)) THEN 1135 IF ((jpjeft(ib1)==jpjnob(ib2)+1).AND.(jpinft(ib2)==jpieob(ib1)+1)) THEN 1136 ! We have a possible North-East corner 1137 ! WRITE(numout,*) ' Found a North-East corner at (i,j): ', jpinft(ib2), jpjeft(ib1) 1138 ! WRITE(numout,*) ' between segments: ', npckge(ib1), npckgn(ib2) 1139 icorne(ib1,2) = npckgn(ib2) 1140 icornn(ib2,2) = npckge(ib1) 1141 ELSEIF ((jpjedt(ib1)==jpjnob(ib2)+1).AND.(jpindt(ib2)==jpieob(ib1)+1)) THEN 1142 IF(lwp) WRITE(numout,*) 1143 IF(lwp) WRITE(numout,*) ' E R R O R : Found an acute open boundary corner at point (i,j)= ', & 1144 & jpindt(ib2), jpjedt(ib1) 1145 IF(lwp) WRITE(numout,*) ' ========== Not allowed yet' 1146 IF(lwp) WRITE(numout,*) ' Crossing problem with East segment: ',npckge(ib1), & 1147 & ' and North segment: ',npckgn(ib2) 1148 IF(lwp) WRITE(numout,*) 1149 nstop = nstop + 1 1150 ELSE 1151 IF(lwp) WRITE(numout,*) 1152 IF(lwp) WRITE(numout,*) ' E R R O R : Check North and East Open boundary indices' 1153 IF(lwp) WRITE(numout,*) ' ========== Crossing problem with East segment: ',npckge(ib1), & 1154 & ' and North segment: ',npckgn(ib2) 1155 IF(lwp) WRITE(numout,*) 1156 nstop = nstop + 1 1157 END IF 1158 END IF 1159 END DO 1160 END DO 1161 END IF 1162 ! 1163 ! 3. Check if segment extremities are on land 1164 !-------------------------------------------- 1165 ! 1166 ! West segments 1167 DO ib = 1, nbdysegw 1168 ! get mask at boundary extremities: 1169 ztestmask(1:2)=0. 1170 DO ji = 1, jpi 1171 DO jj = 1, jpj 1172 IF (((ji + nimpp - 1) == jpiwob(ib)).AND. & 1173 & ((jj + njmpp - 1) == jpjwdt(ib))) ztestmask(1)=tmask(ji,jj,1) 1174 IF (((ji + nimpp - 1) == jpiwob(ib)).AND. & 1175 & ((jj + njmpp - 1) == jpjwft(ib))) ztestmask(2)=tmask(ji,jj,1) 1176 END DO 1177 END DO 1178 IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain 1179 1180 IF (ztestmask(1)==1) THEN 1181 IF (icornw(ib,1)==0) THEN 1182 IF(lwp) WRITE(numout,*) 1183 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckgw(ib) 1184 IF(lwp) WRITE(numout,*) ' ========== does not start on land or on a corner' 1185 IF(lwp) WRITE(numout,*) 1186 nstop = nstop + 1 1187 ELSE 1188 ! This is a corner 1189 WRITE(numout,*) 'Found a South-West corner at (i,j): ', jpiwob(ib), jpjwdt(ib) 1190 CALL bdy_ctl_corn(npckgw(ib), icornw(ib,1)) 1191 itest=itest+1 1192 ENDIF 1193 ENDIF 1194 IF (ztestmask(2)==1) THEN 1195 IF (icornw(ib,2)==0) THEN 1196 IF(lwp) WRITE(numout,*) 1197 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckgw(ib) 1198 IF(lwp) WRITE(numout,*) ' ========== does not end on land or on a corner' 1199 IF(lwp) WRITE(numout,*) 1200 nstop = nstop + 1 1201 ELSE 1202 ! This is a corner 1203 WRITE(numout,*) 'Found a North-West corner at (i,j): ', jpiwob(ib), jpjwft(ib) 1204 CALL bdy_ctl_corn(npckgw(ib), icornw(ib,2)) 1205 itest=itest+1 1206 ENDIF 1207 ENDIF 1208 END DO 1209 ! 1210 ! East segments 1211 DO ib = 1, nbdysege 1212 ! get mask at boundary extremities: 1213 ztestmask(1:2)=0. 1214 DO ji = 1, jpi 1215 DO jj = 1, jpj 1216 IF (((ji + nimpp - 1) == jpieob(ib)+1).AND. & 1217 & ((jj + njmpp - 1) == jpjedt(ib))) ztestmask(1)=tmask(ji,jj,1) 1218 IF (((ji + nimpp - 1) == jpieob(ib)+1).AND. & 1219 & ((jj + njmpp - 1) == jpjeft(ib))) ztestmask(2)=tmask(ji,jj,1) 1220 END DO 1221 END DO 1222 IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain 1223 1224 IF (ztestmask(1)==1) THEN 1225 IF (icorne(ib,1)==0) THEN 1226 IF(lwp) WRITE(numout,*) 1227 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckge(ib) 1228 IF(lwp) WRITE(numout,*) ' ========== does not start on land or on a corner' 1229 IF(lwp) WRITE(numout,*) 1230 nstop = nstop + 1 1231 ELSE 1232 ! This is a corner 1233 WRITE(numout,*) 'Found a South-East corner at (i,j): ', jpieob(ib)+1, jpjedt(ib) 1234 CALL bdy_ctl_corn(npckge(ib), icorne(ib,1)) 1235 itest=itest+1 1236 ENDIF 1237 ENDIF 1238 IF (ztestmask(2)==1) THEN 1239 IF (icorne(ib,2)==0) THEN 1240 IF(lwp) WRITE(numout,*) 1241 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckge(ib) 1242 IF(lwp) WRITE(numout,*) ' ========== does not end on land or on a corner' 1243 IF(lwp) WRITE(numout,*) 1244 nstop = nstop + 1 1245 ELSE 1246 ! This is a corner 1247 WRITE(numout,*) 'Found a North-East corner at (i,j): ', jpieob(ib)+1, jpjeft(ib) 1248 CALL bdy_ctl_corn(npckge(ib), icorne(ib,2)) 1249 itest=itest+1 1250 ENDIF 1251 ENDIF 1252 END DO 1253 ! 1254 ! South segments 1255 DO ib = 1, nbdysegs 1256 ! get mask at boundary extremities: 1257 ztestmask(1:2)=0. 1258 DO ji = 1, jpi 1259 DO jj = 1, jpj 1260 IF (((jj + njmpp - 1) == jpjsob(ib)).AND. & 1261 & ((ji + nimpp - 1) == jpisdt(ib))) ztestmask(1)=tmask(ji,jj,1) 1262 IF (((jj + njmpp - 1) == jpjsob(ib)).AND. & 1263 & ((ji + nimpp - 1) == jpisft(ib))) ztestmask(2)=tmask(ji,jj,1) 1264 END DO 1265 END DO 1266 IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain 1267 1268 IF ((ztestmask(1)==1).AND.(icorns(ib,1)==0)) THEN 1269 IF(lwp) WRITE(numout,*) 1270 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckgs(ib) 1271 IF(lwp) WRITE(numout,*) ' ========== does not start on land or on a corner' 1272 IF(lwp) WRITE(numout,*) 1273 nstop = nstop + 1 1274 ENDIF 1275 IF ((ztestmask(2)==1).AND.(icorns(ib,2)==0)) THEN 1276 IF(lwp) WRITE(numout,*) 1277 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckgs(ib) 1278 IF(lwp) WRITE(numout,*) ' ========== does not end on land or on a corner' 1279 IF(lwp) WRITE(numout,*) 1280 nstop = nstop + 1 1281 ENDIF 1282 END DO 1283 ! 1284 ! North segments 1285 DO ib = 1, nbdysegn 1286 ! get mask at boundary extremities: 1287 ztestmask(1:2)=0. 1288 DO ji = 1, jpi 1289 DO jj = 1, jpj 1290 IF (((jj + njmpp - 1) == jpjnob(ib)+1).AND. & 1291 & ((ji + nimpp - 1) == jpindt(ib))) ztestmask(1)=tmask(ji,jj,1) 1292 IF (((jj + njmpp - 1) == jpjnob(ib)+1).AND. & 1293 & ((ji + nimpp - 1) == jpinft(ib))) ztestmask(2)=tmask(ji,jj,1) 1294 END DO 1295 END DO 1296 IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain 1297 1298 IF ((ztestmask(1)==1).AND.(icornn(ib,1)==0)) THEN 1299 IF(lwp) WRITE(numout,*) 1300 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckgn(ib) 1301 IF(lwp) WRITE(numout,*) ' ========== does not start on land' 1302 IF(lwp) WRITE(numout,*) 1303 nstop = nstop + 1 1304 ENDIF 1305 IF ((ztestmask(2)==1).AND.(icornn(ib,2)==0)) THEN 1306 IF(lwp) WRITE(numout,*) 1307 IF(lwp) WRITE(numout,*) ' E R R O R : Open boundary segment ', npckgn(ib) 1308 IF(lwp) WRITE(numout,*) ' ========== does not end on land' 1309 IF(lwp) WRITE(numout,*) 1310 nstop = nstop + 1 1311 ENDIF 1312 END DO 1313 ! 1314 IF ((itest==0).AND.(lwp)) WRITE(numout,*) 'NO open boundary corner found' 1315 ! 1316 ! Other tests TBD: 1317 ! segments completly on land 1318 ! optimized open boundary array length according to landmask 1319 ! Nudging layers that overlap with interior domain 1320 ! 1321 END SUBROUTINE bdy_ctl_seg 1322 1323 SUBROUTINE bdy_ctl_corn( ib1, ib2 ) 1324 !!---------------------------------------------------------------------- 1325 !! *** ROUTINE bdy_ctl_corn *** 1326 !! 1327 !! ** Purpose : Check numerical schemes consistency between 1328 !! segments having a common corner 1329 !! 1330 !! ** Method : 1331 !!---------------------------------------------------------------------- 1332 INTEGER, INTENT(in) :: ib1, ib2 1333 INTEGER :: itest 1334 !!---------------------------------------------------------------------- 1335 itest = 0 1336 1337 IF (nn_dyn2d(ib1)/=nn_dyn2d(ib2)) itest = itest + 1 1338 IF (nn_dyn3d(ib1)/=nn_dyn3d(ib2)) itest = itest + 1 1339 IF (nn_tra(ib1)/=nn_tra(ib2)) itest = itest + 1 1340 ! 1341 IF (nn_dyn2d_dta(ib1)/=nn_dyn2d_dta(ib2)) itest = itest + 1 1342 IF (nn_dyn3d_dta(ib1)/=nn_dyn3d_dta(ib2)) itest = itest + 1 1343 IF (nn_tra_dta(ib1)/=nn_tra_dta(ib2)) itest = itest + 1 1344 ! 1345 IF (nn_rimwidth(ib1)/=nn_rimwidth(ib2)) itest = itest + 1 1346 ! 1347 IF ( itest>0 ) THEN 1348 IF(lwp) WRITE(numout,*) ' E R R O R : Segments ', ib1, 'and ', ib2 1349 IF(lwp) WRITE(numout,*) ' ========== have different open bdy schemes' 1350 IF(lwp) WRITE(numout,*) 1351 nstop = nstop + 1 1352 ENDIF 1353 ! 1354 END SUBROUTINE bdy_ctl_corn 840 1355 841 1356 #else
Note: See TracChangeset
for help on using the changeset viewer.