Changeset 13337 for NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_oce_sponge.F90
- Timestamp:
- 2020-07-24T16:01:24+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_oce_sponge.F90
r13312 r13337 55 55 Agrif_SpecialValue = 0._wp 56 56 Agrif_UseSpecialValue = .TRUE. 57 l_vremap = ln_vremap 57 58 tabspongedone_tsn = .FALSE. 58 59 ! … … 60 61 ! 61 62 Agrif_UseSpecialValue = .FALSE. 63 l_vremap = .FALSE. 62 64 #endif 63 65 ! … … 80 82 Agrif_SpecialValue = 0._wp 81 83 Agrif_UseSpecialValue = ln_spc_dyn 84 l_vremap = ln_vremap 82 85 use_sign_north = .TRUE. 83 86 sign_north = -1._wp … … 93 96 Agrif_UseSpecialValue = .FALSE. 94 97 use_sign_north = .FALSE. 98 l_vremap = .FALSE. 95 99 #endif 96 100 ! … … 109 113 REAL(wp) :: z1_ispongearea, z1_jspongearea 110 114 REAL(wp), DIMENSION(jpi,jpj) :: ztabramp 111 #if defined key_vertical112 115 REAL(wp), DIMENSION(jpi,jpj) :: ztabrampu 113 116 REAL(wp), DIMENSION(jpi,jpj) :: ztabrampv 114 #endif115 REAL(wp), DIMENSION(jpjmax) :: zmskwest, zmskeast116 REAL(wp), DIMENSION(jpimax) :: zmsknorth, zmsksouth117 117 !!---------------------------------------------------------------------- 118 118 ! … … 130 130 #if defined SPONGE || defined SPONGE_TOP 131 131 IF (( .NOT. spongedoneT ).OR.( .NOT. spongedoneU )) THEN 132 !133 ! Retrieve masks at open boundaries:134 135 IF( lk_west ) THEN ! --- West --- !136 ztabramp(:,:) = 0._wp137 ind1 = nn_hls + 1 + nbghostcells ! halo + land + nbghostcells138 DO ji = mi0(ind1), mi1(ind1)139 ztabramp(ji,:) = ssumask(ji,:)140 END DO141 zmskwest( 1:jpj ) = MAXVAL(ztabramp(:,:), dim=1)142 zmskwest(jpj+1:jpjmax) = 0._wp143 ENDIF144 IF( lk_east ) THEN ! --- East --- !145 ztabramp(:,:) = 0._wp146 ind1 = jpiglo - ( nn_hls + nbghostcells + 1) ! halo + land + nbghostcells147 DO ji = mi0(ind1), mi1(ind1)148 ztabramp(ji,:) = ssumask(ji,:)149 END DO150 zmskeast( 1:jpj ) = MAXVAL(ztabramp(:,:), dim=1)151 zmskeast(jpj+1:jpjmax) = 0._wp152 ENDIF153 IF( lk_south ) THEN ! --- South --- !154 ztabramp(:,:) = 0._wp155 ind1 = nn_hls + 1 + nbghostcells ! halo + land + nbghostcells156 DO jj = mj0(ind1), mj1(ind1)157 ztabramp(:,jj) = ssvmask(:,jj)158 END DO159 zmsksouth( 1:jpi ) = MAXVAL(ztabramp(:,:), dim=2)160 zmsksouth(jpi+1:jpimax) = 0._wp161 ENDIF162 IF( lk_north ) THEN ! --- North --- !163 ztabramp(:,:) = 0._wp164 ind1 = jpjglo - ( nn_hls + nbghostcells + 1) ! halo + land + nbghostcells165 DO jj = mj0(ind1), mj1(ind1)166 ztabramp(:,jj) = ssvmask(:,jj)167 END DO168 zmsknorth( 1:jpi ) = MAXVAL(ztabramp(:,:), dim=2)169 zmsknorth(jpi+1:jpimax) = 0._wp170 ENDIF171 172 ! JC: SPONGE MASKING TO BE SORTED OUT:173 zmskwest(:) = 1._wp174 zmskeast(:) = 1._wp175 zmsksouth(:) = 1._wp176 zmsknorth(:) = 1._wp177 #if defined key_mpp_mpi178 ! CALL mpp_max( 'AGRIF_sponge', zmskwest(:) , jpjmax )179 ! CALL mpp_max( 'AGRIF_Sponge', zmskeast(:) , jpjmax )180 ! CALL mpp_max( 'AGRIF_Sponge', zmsksouth(:), jpimax )181 ! CALL mpp_max( 'AGRIF_Sponge', zmsknorth(:), jpimax )182 #endif183 184 132 ! Define ramp from boundaries towards domain interior at T-points 185 133 ! Store it in ztabramp … … 201 149 DO ji = mi0(ind1), mi1(ind2) 202 150 DO jj = 1, jpj 203 ztabramp(ji,jj) = REAL(ind2 - mig(ji), wp) * z1_ispongearea * zmskwest(jj)151 ztabramp(ji,jj) = REAL(ind2 - mig(ji), wp) * z1_ispongearea 204 152 END DO 205 153 END DO … … 209 157 DO ji = mi0(ind1), mi1(ind2) 210 158 DO jj = 1, jpj 211 ztabramp(ji,jj) = zmskwest(jj)159 ztabramp(ji,jj) = 1._wp 212 160 END DO 213 161 END DO … … 218 166 DO ji = mi0(ind1), mi1(ind2) 219 167 DO jj = 1, jpj 220 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mig(ji) - ind1, wp) * z1_ispongearea ) * zmskeast(jj)168 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mig(ji) - ind1, wp) * z1_ispongearea ) 221 169 END DO 222 170 END DO … … 226 174 DO ji = mi0(ind1), mi1(ind2) 227 175 DO jj = 1, jpj 228 ztabramp(ji,jj) = zmskeast(jj)176 ztabramp(ji,jj) = 1._wp 229 177 END DO 230 178 END DO … … 235 183 DO jj = mj0(ind1), mj1(ind2) 236 184 DO ji = 1, jpi 237 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(ind2 - mjg(jj), wp) * z1_jspongearea ) * zmsksouth(ji)185 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(ind2 - mjg(jj), wp) * z1_jspongearea ) 238 186 END DO 239 187 END DO … … 243 191 DO jj = mj0(ind1), mj1(ind2) 244 192 DO ji = 1, jpi 245 ztabramp(ji,jj) = zmsksouth(ji)193 ztabramp(ji,jj) = 1._wp 246 194 END DO 247 195 END DO … … 252 200 DO jj = mj0(ind1), mj1(ind2) 253 201 DO ji = 1, jpi 254 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mjg(jj) - ind1, wp) * z1_jspongearea ) * zmsknorth(ji)202 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mjg(jj) - ind1, wp) * z1_jspongearea ) 255 203 END DO 256 204 END DO … … 260 208 DO jj = mj0(ind1), mj1(ind2) 261 209 DO ji = 1, jpi 262 ztabramp(ji,jj) = zmsknorth(ji)210 ztabramp(ji,jj) = 1._wp 263 211 END DO 264 212 END DO … … 303 251 ENDIF 304 252 305 #if defined key_vertical306 253 ! Remove vertical interpolation where not needed: 307 254 DO_2D( 0, 0, 0, 0 ) … … 327 274 mbku_parent(:,:) = NINT( ztabrampu(:,:) ) 328 275 mbkv_parent(:,:) = NINT( ztabrampv(:,:) ) 329 #endif330 276 ! 331 277 #endif … … 366 312 END DO 367 313 368 # if defined key_vertical 369 ! Interpolate thicknesses 370 ! Warning: these are masked, hence extrapolated prior interpolation. 371 DO jk=k1,k2 372 DO jj=j1,j2 373 DO ji=i1,i2 374 tabres(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kbb_a) 314 IF ( l_vremap ) THEN 315 316 ! Interpolate thicknesses 317 ! Warning: these are masked, hence extrapolated prior interpolation. 318 DO jk=k1,k2 319 DO jj=j1,j2 320 DO ji=i1,i2 321 tabres(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kbb_a) 322 END DO 375 323 END DO 376 324 END DO 377 END DO 378 379 ! Extrapolate thicknesses in partial bottom cells: 380 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 381 IF (ln_zps) THEN 382 DO jj=j1,j2 383 DO ji=i1,i2 384 jk = mbkt(ji,jj) 385 tabres(ji,jj,jk,jpts+1) = 0._wp 386 END DO 387 END DO 325 326 ! Extrapolate thicknesses in partial bottom cells: 327 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 328 IF (ln_zps) THEN 329 DO jj=j1,j2 330 DO ji=i1,i2 331 jk = mbkt(ji,jj) 332 tabres(ji,jj,jk,jpts+1) = 0._wp 333 END DO 334 END DO 335 END IF 336 337 ! Save ssh at last level: 338 IF (.NOT.ln_linssh) THEN 339 tabres(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kbb_a)*tmask(i1:i2,j1:j2,1) 340 ELSE 341 tabres(i1:i2,j1:j2,k2,jpts+1) = 0._wp 342 END IF 388 343 END IF 389 390 ! Save ssh at last level:391 IF (.NOT.ln_linssh) THEN392 tabres(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kbb_a)*tmask(i1:i2,j1:j2,1)393 ELSE394 tabres(i1:i2,j1:j2,k2,jpts+1) = 0._wp395 END IF396 # endif397 344 398 345 ELSE 399 346 ! 400 # if defined key_vertical 401 402 IF (ln_linssh) tabres(i1:i2,j1:j2,k2,n2) = 0._wp 403 404 DO jj=j1,j2 405 DO ji=i1,i2 406 tabres_child(ji,jj,:,:) = 0._wp 407 N_in = mbkt_parent(ji,jj) 408 zhtot = 0._wp 409 DO jk=1,N_in !k2 = jpk of parent grid 410 IF (jk==N_in) THEN 411 h_in(jk) = ht0_parent(ji,jj) + tabres(ji,jj,k2,n2) - zhtot 347 IF ( l_vremap ) THEN 348 349 IF (ln_linssh) tabres(i1:i2,j1:j2,k2,n2) = 0._wp 350 351 DO jj=j1,j2 352 DO ji=i1,i2 353 tabres_child(ji,jj,:,:) = 0._wp 354 N_in = mbkt_parent(ji,jj) 355 zhtot = 0._wp 356 DO jk=1,N_in !k2 = jpk of parent grid 357 IF (jk==N_in) THEN 358 h_in(jk) = ht0_parent(ji,jj) + tabres(ji,jj,k2,n2) - zhtot 359 ELSE 360 h_in(jk) = tabres(ji,jj,jk,n2) 361 END IF 362 zhtot = zhtot + h_in(jk) 363 tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1) 364 END DO 365 N_out = 0 366 DO jk=1,jpk ! jpk of child grid 367 IF (tmask(ji,jj,jk) == 0) EXIT 368 N_out = N_out + 1 369 h_out(jk) = e3t(ji,jj,jk,Kbb_a) !Child grid scale factors. Could multiply by e1e2t here instead of division above 370 END DO 371 372 ! Account for small differences in the free-surface 373 IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN 374 h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) ) 412 375 ELSE 413 h_in(jk) = tabres(ji,jj,jk,n2) 376 h_in(1) = h_in(1) - (sum(h_in(1:N_in))-sum(h_out(1:N_out)) ) 377 END IF 378 IF (N_in*N_out > 0) THEN 379 CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),tabres_child(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts) 414 380 ENDIF 415 zhtot = zhtot + h_in(jk) 416 tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1) 417 END DO 418 N_out = 0 419 DO jk=1,jpk ! jpk of child grid 420 IF (tmask(ji,jj,jk) == 0) EXIT 421 N_out = N_out + 1 422 h_out(jk) = e3t(ji,jj,jk,Kbb_a) !Child grid scale factors. Could multiply by e1e2t here instead of division above 423 END DO 424 425 ! Account for small differences in free-surface 426 IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN 427 h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) ) 428 ELSE 429 h_in(1) = h_in(1) - (sum(h_in(1:N_in))-sum(h_out(1:N_out)) ) 430 ENDIF 431 IF (N_in*N_out > 0) THEN 432 CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),tabres_child(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts) 433 ENDIF 434 END DO 435 END DO 436 # endif 437 438 DO jj=j1,j2 439 DO ji=i1,i2 440 DO jk=1,jpkm1 441 # if defined key_vertical 442 tsbdiff(ji,jj,jk,1:jpts) = (ts(ji,jj,jk,1:jpts,Kbb_a) - tabres_child(ji,jj,jk,1:jpts)) * tmask(ji,jj,jk) 443 # else 444 tsbdiff(ji,jj,jk,1:jpts) = (ts(ji,jj,jk,1:jpts,Kbb_a) - tabres(ji,jj,jk,1:jpts))*tmask(ji,jj,jk) 445 # endif 446 END DO 447 END DO 448 END DO 381 END DO 382 END DO 383 384 DO jj=j1,j2 385 DO ji=i1,i2 386 DO jk=1,jpkm1 387 tsbdiff(ji,jj,jk,1:jpts) = (ts(ji,jj,jk,1:jpts,Kbb_a) - tabres_child(ji,jj,jk,1:jpts)) * tmask(ji,jj,jk) 388 END DO 389 END DO 390 END DO 391 392 ELSE 393 394 DO jj=j1,j2 395 DO ji=i1,i2 396 DO jk=1,jpkm1 397 tsbdiff(ji,jj,jk,1:jpts) = (ts(ji,jj,jk,1:jpts,Kbb_a) - tabres(ji,jj,jk,1:jpts))*tmask(ji,jj,jk) 398 END DO 399 END DO 400 END DO 401 END IF 449 402 450 403 DO jn = 1, jpts … … 528 481 DO ji=i1,i2 529 482 tabres(ji,jj,jk,m1) = uu(ji,jj,jk,Kbb_a) 530 # if defined key_vertical 531 tabres(ji,jj,jk,m2) = e3u(ji,jj,jk,Kbb_a)*umask(ji,jj,jk) 532 # endif 533 END DO 534 END DO 535 END DO 536 537 # if defined key_vertical 538 ! Extrapolate thicknesses in partial bottom cells: 539 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 540 IF (ln_zps) THEN 483 END DO 484 END DO 485 END DO 486 487 IF ( l_vremap ) THEN 488 489 DO jk=k1,k2 490 DO jj=j1,j2 491 DO ji=i1,i2 492 tabres(ji,jj,jk,m2) = e3u(ji,jj,jk,Kbb_a)*umask(ji,jj,jk) 493 END DO 494 END DO 495 END DO 496 497 ! Extrapolate thicknesses in partial bottom cells: 498 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 499 IF (ln_zps) THEN 500 DO jj=j1,j2 501 DO ji=i1,i2 502 jk = mbku(ji,jj) 503 tabres(ji,jj,jk,m2) = 0._wp 504 END DO 505 END DO 506 END IF 507 ! Save ssh at last level: 508 tabres(i1:i2,j1:j2,k2,m2) = 0._wp 509 IF (.NOT.ln_linssh) THEN 510 ! This vertical sum below should be replaced by the sea-level at U-points (optimization): 511 DO jk=1,jpk 512 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) + e3u(i1:i2,j1:j2,jk,Kbb_a) * umask(i1:i2,j1:j2,jk) 513 END DO 514 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) - hu_0(i1:i2,j1:j2) 515 END IF 516 END IF 517 518 ELSE 519 520 IF ( l_vremap ) THEN 521 522 IF (ln_linssh) tabres(i1:i2,j1:j2,k2,m2) = 0._wp 523 541 524 DO jj=j1,j2 542 525 DO ji=i1,i2 543 jk = mbku(ji,jj) 544 tabres(ji,jj,jk,m2) = 0._wp 545 END DO 546 END DO 547 END IF 548 ! Save ssh at last level: 549 tabres(i1:i2,j1:j2,k2,m2) = 0._wp 550 IF (.NOT.ln_linssh) THEN 551 ! This vertical sum below should be replaced by the sea-level at U-points (optimization): 552 DO jk=1,jpk 553 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) + e3u(i1:i2,j1:j2,jk,Kbb_a) * umask(i1:i2,j1:j2,jk) 554 END DO 555 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) - hu_0(i1:i2,j1:j2) 556 END IF 557 # endif 558 559 ELSE 560 561 # if defined key_vertical 562 IF (ln_linssh) tabres(i1:i2,j1:j2,k2,m2) = 0._wp 563 564 DO jj=j1,j2 565 DO ji=i1,i2 566 tabres_child(ji,jj,:) = 0._wp 567 N_in = mbku_parent(ji,jj) 568 zhtot = 0._wp 569 DO jk=1,N_in 570 IF (jk==N_in) THEN 571 h_in(jk) = hu0_parent(ji,jj) + tabres(ji,jj,k2,m2) - zhtot 526 tabres_child(ji,jj,:) = 0._wp 527 N_in = mbku_parent(ji,jj) 528 zhtot = 0._wp 529 DO jk=1,N_in 530 IF (jk==N_in) THEN 531 h_in(jk) = hu0_parent(ji,jj) + tabres(ji,jj,k2,m2) - zhtot 532 ELSE 533 h_in(jk) = tabres(ji,jj,jk,m2) 534 ENDIF 535 zhtot = zhtot + h_in(jk) 536 tabin(jk) = tabres(ji,jj,jk,m1) 537 END DO 538 ! 539 N_out = 0 540 DO jk=1,jpk 541 IF (umask(ji,jj,jk) == 0) EXIT 542 N_out = N_out + 1 543 h_out(N_out) = e3u(ji,jj,jk,Kbb_a) 544 END DO 545 546 ! Account for small differences in free-surface 547 IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN 548 h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) ) 572 549 ELSE 573 h_in( jk) = tabres(ji,jj,jk,m2)550 h_in(1) = h_in(1) - (sum(h_in(1:N_in))-sum(h_out(1:N_out)) ) 574 551 ENDIF 575 zhtot = zhtot + h_in(jk)576 tabin(jk) = tabres(ji,jj,jk,m1)577 END DO578 !579 N_out = 0580 DO jk=1,jpk581 IF (umask(ji,jj,jk) == 0) EXIT582 N_out = N_out + 1583 h_out(N_out) = e3u(ji,jj,jk,Kbb_a)584 END DO585 586 ! Account for small differences in free-surface587 IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN588 h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) )589 ELSE590 h_in(1) = h_in(1) - (sum(h_in(1:N_in))-sum(h_out(1:N_out)) )591 ENDIF592 552 593 IF (N_in * N_out > 0) THEN 594 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 595 ENDIF 596 END DO 597 END DO 598 599 ubdiff(i1:i2,j1:j2,:) = (uu(i1:i2,j1:j2,:,Kbb_a) - tabres_child(i1:i2,j1:j2,:))*umask(i1:i2,j1:j2,:) 600 #else 601 ubdiff(i1:i2,j1:j2,:) = (uu(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*umask(i1:i2,j1:j2,:) 602 #endif 553 IF (N_in * N_out > 0) THEN 554 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 555 ENDIF 556 END DO 557 END DO 558 559 ubdiff(i1:i2,j1:j2,:) = (uu(i1:i2,j1:j2,:,Kbb_a) - tabres_child(i1:i2,j1:j2,:))*umask(i1:i2,j1:j2,:) 560 561 ELSE 562 563 ubdiff(i1:i2,j1:j2,:) = (uu(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*umask(i1:i2,j1:j2,:) 564 565 ENDIF 603 566 ! 604 567 DO jk = 1, jpkm1 ! Horizontal slab … … 705 668 DO ji=i1,i2 706 669 tabres(ji,jj,jk,m1) = vv(ji,jj,jk,Kbb_a) 707 # if defined key_vertical 708 tabres(ji,jj,jk,m2) = vmask(ji,jj,jk) * e3v(ji,jj,jk,Kbb_a) 709 # endif 710 END DO 711 END DO 712 END DO 713 714 # if defined key_vertical 715 ! Extrapolate thicknesses in partial bottom cells: 716 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 717 IF (ln_zps) THEN 670 END DO 671 END DO 672 END DO 673 674 IF ( l_vremap ) THEN 675 676 DO jk=k1,k2 677 DO jj=j1,j2 678 DO ji=i1,i2 679 tabres(ji,jj,jk,m2) = vmask(ji,jj,jk) * e3v(ji,jj,jk,Kbb_a) 680 END DO 681 END DO 682 END DO 683 ! Extrapolate thicknesses in partial bottom cells: 684 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 685 IF (ln_zps) THEN 686 DO jj=j1,j2 687 DO ji=i1,i2 688 jk = mbkv(ji,jj) 689 tabres(ji,jj,jk,m2) = 0._wp 690 END DO 691 END DO 692 END IF 693 ! Save ssh at last level: 694 tabres(i1:i2,j1:j2,k2,m2) = 0._wp 695 IF (.NOT.ln_linssh) THEN 696 ! This vertical sum below should be replaced by the sea-level at V-points (optimization): 697 DO jk=1,jpk 698 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) + e3v(i1:i2,j1:j2,jk,Kbb_a) * vmask(i1:i2,j1:j2,jk) 699 END DO 700 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) - hv_0(i1:i2,j1:j2) 701 END IF 702 703 END IF 704 705 ELSE 706 707 IF ( l_vremap ) THEN 708 IF (ln_linssh) tabres(i1:i2,j1:j2,k2,m2) = 0._wp 718 709 DO jj=j1,j2 719 710 DO ji=i1,i2 720 jk = mbkv(ji,jj) 721 tabres(ji,jj,jk,m2) = 0._wp 722 END DO 723 END DO 724 END IF 725 ! Save ssh at last level: 726 tabres(i1:i2,j1:j2,k2,m2) = 0._wp 727 IF (.NOT.ln_linssh) THEN 728 ! This vertical sum below should be replaced by the sea-level at V-points (optimization): 729 DO jk=1,jpk 730 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) + e3v(i1:i2,j1:j2,jk,Kbb_a) * vmask(i1:i2,j1:j2,jk) 731 END DO 732 tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) - hv_0(i1:i2,j1:j2) 733 END IF 734 # endif 735 736 ELSE 737 738 # if defined key_vertical 739 IF (ln_linssh) tabres(i1:i2,j1:j2,k2,m2) = 0._wp 740 DO jj=j1,j2 741 DO ji=i1,i2 742 tabres_child(ji,jj,:) = 0._wp 743 N_in = mbkv_parent(ji,jj) 744 zhtot = 0._wp 745 DO jk=1,N_in 746 IF (jk==N_in) THEN 747 h_in(jk) = hv0_parent(ji,jj) + tabres(ji,jj,k2,m2) - zhtot 711 tabres_child(ji,jj,:) = 0._wp 712 N_in = mbkv_parent(ji,jj) 713 zhtot = 0._wp 714 DO jk=1,N_in 715 IF (jk==N_in) THEN 716 h_in(jk) = hv0_parent(ji,jj) + tabres(ji,jj,k2,m2) - zhtot 717 ELSE 718 h_in(jk) = tabres(ji,jj,jk,m2) 719 ENDIF 720 zhtot = zhtot + h_in(jk) 721 tabin(jk) = tabres(ji,jj,jk,m1) 722 END DO 723 ! 724 N_out = 0 725 DO jk=1,jpk 726 IF (vmask(ji,jj,jk) == 0) EXIT 727 N_out = N_out + 1 728 h_out(N_out) = e3v(ji,jj,jk,Kbb_a) 729 END DO 730 731 ! Account for small differences in free-surface 732 IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN 733 h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) ) 748 734 ELSE 749 h_in( jk) = tabres(ji,jj,jk,m2)735 h_in(1) = h_in(1) - ( sum(h_in(1:N_in))-sum(h_out(1:N_out)) ) 750 736 ENDIF 751 zhtot = zhtot + h_in(jk)752 tabin(jk) = tabres(ji,jj,jk,m1)753 END DO754 !755 N_out = 0756 DO jk=1,jpk757 IF (vmask(ji,jj,jk) == 0) EXIT758 N_out = N_out + 1759 h_out(N_out) = e3v(ji,jj,jk,Kbb_a)760 END DO761 762 ! Account for small differences in free-surface763 IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN764 h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) )765 ELSE766 h_in(1) = h_in(1) - ( sum(h_in(1:N_in))-sum(h_out(1:N_out)) )767 ENDIF768 737 769 IF (N_in * N_out > 0) THEN 770 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 771 ENDIF 772 END DO 773 END DO 774 775 vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres_child(i1:i2,j1:j2,:))*vmask(i1:i2,j1:j2,:) 776 # else 777 vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*vmask(i1:i2,j1:j2,:) 778 # endif 738 IF (N_in * N_out > 0) THEN 739 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 740 ENDIF 741 END DO 742 END DO 743 744 vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres_child(i1:i2,j1:j2,:))*vmask(i1:i2,j1:j2,:) 745 746 ELSE 747 748 vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*vmask(i1:i2,j1:j2,:) 749 750 ENDIF 779 751 ! 780 752 DO jk = 1, jpkm1 ! Horizontal slab
Note: See TracChangeset
for help on using the changeset viewer.
