- Timestamp:
- 2018-09-17T15:16:43+02:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2018/dev_r9759_HPC09_ESIWACE/src/OCE/TRA/traadv_fct.F90
r10103 r10136 27 27 USE lbclnk ! ocean lateral boundary condition (or mpp link) 28 28 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 29 USE timing ! Timing 29 30 30 31 IMPLICIT NONE … … 325 326 326 327 SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt ) 328 #ifdef SCOREP_USER_ENABLE 329 use mpi 330 #include "scorep/SCOREP_User.inc" 331 #endif 327 332 !!--------------------------------------------------------------------- 328 333 !! *** ROUTINE nonosc *** … … 346 351 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zup, zdo ! - - 347 352 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zbetup, zbetdo, zbup, zbdo 353 !dir$ attributes align:64 :: zbetup, zbetdo, zbup, zbdo 354 #ifdef SCOREP_USER_ENABLE 355 integer :: ier 356 SCOREP_USER_REGION_DEFINE( reg_nonosc ) 357 SCOREP_USER_REGION_DEFINE( reg_nonosc_setup ) 358 SCOREP_USER_REGION_DEFINE( reg_nonosc_cb1 ) 359 SCOREP_USER_REGION_DEFINE( reg_nonosc_cb2 ) 360 SCOREP_USER_REGION_DEFINE( reg_nonosc_barrier ) 361 SCOREP_USER_REGION_DEFINE( reg_nonosc_imbalance ) 362 363 SCOREP_USER_REGION_BEGIN( reg_nonosc_barrier, "nonosc barrier", SCOREP_USER_REGION_TYPE_COMMON ) 364 call MPI_Barrier(MPI_COMM_WORLD, ier) 365 SCOREP_USER_REGION_END( reg_nonosc_barrier ) 366 SCOREP_USER_REGION_BEGIN( reg_nonosc, "nonosc", SCOREP_USER_REGION_TYPE_FUNCTION ) 367 SCOREP_USER_REGION_BEGIN( reg_nonosc_setup, "nonosc setup", SCOREP_USER_REGION_TYPE_COMMON ) 368 #endif 369 IF( ln_timing ) CALL timing_start( 'nonosc' ) 348 370 !!---------------------------------------------------------------------- 349 371 ! 350 372 zbig = 1.e+40_wp 351 373 zrtrn = 1.e-15_wp 374 #ifndef BULL_NONOSC_INIT 352 375 zbetup(:,:,:) = 0._wp ; zbetdo(:,:,:) = 0._wp 376 #else 377 zbetup(:,:,jpk) = 0._wp ; zbetdo(:,:,jpk) = 0._wp 378 #endif 353 379 354 380 ! Search local extrema … … 360 386 & paft * tmask + zbig * ( 1._wp - tmask ) ) 361 387 388 #ifdef SCOREP_USER_ENABLE 389 SCOREP_USER_REGION_END( reg_nonosc_setup ) 390 #endif 391 392 #ifndef BULL_ASYNC 393 #ifdef SCOREP_USER_ENABLE 394 SCOREP_USER_REGION_BEGIN( reg_nonosc_cb1, "cb1", SCOREP_USER_REGION_TYPE_LOOP ) 395 #endif 396 ! loads: 397 ! - zbup: ji-1/ji/ji+1, jj-1/jj/jj+1, ji/jk+1/jk-1 398 ! - zbdo: " 399 ! - paa: ji-1/ji 400 ! - pbb: jj-1/jj 401 ! - pcc: ji, jj, jk/jk+1 402 ! - e1e2t, e3t_n, paft (*2): ji,jj,jk 403 ! 404 ! stores: 405 ! - zbetup 406 ! - zbetdo 362 407 DO jk = 1, jpkm1 363 408 ikm1 = MAX(jk-1,1) … … 394 439 END DO 395 440 END DO 441 #ifdef SCOREP_USER_ENABLE 442 SCOREP_USER_REGION_END( reg_nonosc_cb1 ) 443 #endif 396 444 CALL lbc_lnk_multi("traadv_fct",zbetup, 'T', 1. , zbetdo, 'T', 1. ) ! lateral boundary cond. (unchanged sign) 397 445 #else 446 call lbc_lnk_multi_async( "traadv_fct", cb1, zbetup, 'T', 1. , zbetdo, 'T', 1. ) ! lateral boundary cond. (unchanged sign) 447 #endif 448 449 #ifndef BULL_ASYNC 398 450 ! 3. monotonic flux in the i & j direction (paa & pbb) 399 451 ! ---------------------------------------- 452 #ifdef SCOREP_USER_ENABLE 453 SCOREP_USER_REGION_BEGIN( reg_nonosc_cb2, "cb2", SCOREP_USER_REGION_TYPE_LOOP ) 454 #endif 400 455 DO jk = 1, jpkm1 401 456 DO jj = 2, jpjm1 … … 420 475 END DO 421 476 END DO 477 #ifdef SCOREP_USER_ENABLE 478 SCOREP_USER_REGION_END( reg_nonosc_cb2 ) 479 #endif 422 480 CALL lbc_lnk_multi("traadv_fct",paa, 'U', -1. , pbb, 'V', -1. ) ! lateral boundary condition (changed sign) 423 ! 481 #else 482 call lbc_lnk_multi_async( "traadv_fct", cb2, paa, 'U', -1. , pbb, 'V', -1. ) ! lateral boundary condition (changed sign) 483 #endif 484 ! 485 IF( ln_timing ) CALL timing_stop( 'nonosc' ) 486 #ifdef SCOREP_USER_ENABLE 487 SCOREP_USER_REGION_END( reg_nonosc ) 488 SCOREP_USER_REGION_BEGIN( reg_nonosc_imbalance, "nonosc imbalance", SCOREP_USER_REGION_TYPE_COMMON ) 489 call MPI_Barrier(MPI_COMM_WORLD, ier) 490 SCOREP_USER_REGION_END( reg_nonosc_imbalance ) 491 #endif 492 #ifdef BULL_ASYNC 493 contains 494 subroutine cb1(i0, i1, j0, j1, k0, k1, buf) 495 integer, intent(in) :: i0, i1, j0, j1, k0, k1 496 real(wp), dimension(:,:,:,:,:,:), optional, intent(out) :: buf 497 integer jji, jjj, jjk 498 real(wp) :: p2dt_inv 499 !REAL(wp), DIMENSION (40,jpj,jpk) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions 500 !REAL(wp), DIMENSION (40,jpj,jpk) :: e3t_n, paft 501 !REAL(wp), DIMENSION (40,jpj) :: e1e2t 502 !REAL(wp), DIMENSION(40,jpj,jpk) :: zbetup, zbetdo, zbup, zbdo 503 !DIR$ ASSUME_ALIGNED zbup:64 504 !DIR$ ASSUME (jpi .EQ.40) 505 !DIR$ ASSUME (jpj .EQ.42) 506 !DIR$ ASSUME (jpk .EQ.75) 507 508 p2dt_inv = 1._wp * p2dt 509 if(i0 == i1) then 510 ji=i0 511 ! DO jjj = j0, j1, 8 512 DO jk = k0, k1 513 ikm1 = MAX(jk-1,1) 514 !DIR$ vector always 515 DO jj = j0, j1 516 !DO jj = jjj, min(jjj+7, j1) 517 ! search maximum in neighbourhood 518 zup = MAX( zbup(ji ,jj ,jk ), & 519 & zbup(ji-1,jj ,jk ), zbup(ji+1,jj ,jk ), & 520 & zbup(ji ,jj-1,jk ), zbup(ji ,jj+1,jk ), & 521 & zbup(ji ,jj ,ikm1), zbup(ji ,jj ,jk+1) ) 522 523 ! search minimum in neighbourhood 524 zdo = MIN( zbdo(ji ,jj ,jk ), & 525 & zbdo(ji-1,jj ,jk ), zbdo(ji+1,jj ,jk ), & 526 & zbdo(ji ,jj-1,jk ), zbdo(ji ,jj+1,jk ), & 527 & zbdo(ji ,jj ,ikm1), zbdo(ji ,jj ,jk+1) ) 528 529 ! positive part of the flux 530 zpos = MAX( 0., paa(ji-1,jj ,jk ) ) - MIN( 0., paa(ji ,jj ,jk ) ) & 531 & + MAX( 0., pbb(ji ,jj-1,jk ) ) - MIN( 0., pbb(ji ,jj ,jk ) ) & 532 & + MAX( 0., pcc(ji ,jj ,jk+1) ) - MIN( 0., pcc(ji ,jj ,jk ) ) 533 534 ! negative part of the flux 535 zneg = MAX( 0., paa(ji ,jj ,jk ) ) - MIN( 0., paa(ji-1,jj ,jk ) ) & 536 & + MAX( 0., pbb(ji ,jj ,jk ) ) - MIN( 0., pbb(ji ,jj-1,jk ) ) & 537 & + MAX( 0., pcc(ji ,jj ,jk ) ) - MIN( 0., pcc(ji ,jj ,jk+1) ) 538 539 ! up & down beta terms 540 zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) * p2dt_inv 541 zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 542 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt 543 544 #ifdef BULL_CB_WITH_BUF 545 ! zt3ew(:,jh,jk,jl,jf,1) = ARRAY_IN(nn_hls+jh,:,jk,jl,jf 546 buf(jj,1,jk,1,1,1) = zbetup(ji,jj,jk) 547 buf(jj,1,jk,1,2,1) = zbetdo(ji,jj,jk) 548 #endif 549 END DO 550 END DO 551 !end do 552 else 553 DO jk = k0, k1 554 ikm1 = MAX(jk-1,1) 555 DO jj = j0, j1 556 !DIR$ vector always 557 DO ji = i0, i1 558 559 ! search maximum in neighbourhood 560 zup = MAX( zbup(ji ,jj ,jk ), & 561 & zbup(ji-1,jj ,jk ), zbup(ji+1,jj ,jk ), & 562 & zbup(ji ,jj-1,jk ), zbup(ji ,jj+1,jk ), & 563 & zbup(ji ,jj ,ikm1), zbup(ji ,jj ,jk+1) ) 564 565 ! search minimum in neighbourhood 566 zdo = MIN( zbdo(ji ,jj ,jk ), & 567 & zbdo(ji-1,jj ,jk ), zbdo(ji+1,jj ,jk ), & 568 & zbdo(ji ,jj-1,jk ), zbdo(ji ,jj+1,jk ), & 569 & zbdo(ji ,jj ,ikm1), zbdo(ji ,jj ,jk+1) ) 570 571 ! positive part of the flux 572 zpos = MAX( 0., paa(ji-1,jj ,jk ) ) - MIN( 0., paa(ji ,jj ,jk ) ) & 573 & + MAX( 0., pbb(ji ,jj-1,jk ) ) - MIN( 0., pbb(ji ,jj ,jk ) ) & 574 & + MAX( 0., pcc(ji ,jj ,jk+1) ) - MIN( 0., pcc(ji ,jj ,jk ) ) 575 576 ! negative part of the flux 577 zneg = MAX( 0., paa(ji ,jj ,jk ) ) - MIN( 0., paa(ji-1,jj ,jk ) ) & 578 & + MAX( 0., pbb(ji ,jj ,jk ) ) - MIN( 0., pbb(ji ,jj-1,jk ) ) & 579 & + MAX( 0., pcc(ji ,jj ,jk ) ) - MIN( 0., pcc(ji ,jj ,jk+1) ) 580 581 ! up & down beta terms 582 zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) * p2dt_inv 583 zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 584 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt 585 586 END DO 587 END DO 588 END DO 589 endif 590 591 end subroutine 592 subroutine cb2(i0, i1, j0, j1, k0, k1, buf) 593 integer, intent(in) :: i0, i1, j0, j1, k0, k1 594 real(wp), dimension(:,:,:,:,:,:), optional, intent(out) :: buf 595 integer jji, jjj, jjk 596 597 ! 3. monotonic flux in the i & j direction (paa & pbb) 598 if(i0 == i1) then 599 ji=i0 600 do jjj=j0, j1, 8 601 DO jk = k0, k1 602 !DIR$ vector always 603 !DO jj = j0, j1 604 DO jj = jjj, min(jjj+7, j1) 605 zau = MIN( 1._wp, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 606 zbu = MIN( 1._wp, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) 607 zcu = ( 0.5 + SIGN( 0.5 , paa(ji,jj,jk) ) ) 608 paa(ji,jj,jk) = paa(ji,jj,jk) * ( zcu * zau + ( 1._wp - zcu) * zbu ) 609 610 zav = MIN( 1._wp, zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) ) 611 zbv = MIN( 1._wp, zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) ) 612 zcv = ( 0.5 + SIGN( 0.5 , pbb(ji,jj,jk) ) ) 613 pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zcv * zav + ( 1._wp - zcv) * zbv ) 614 615 ! monotonic flux in the k direction, i.e. pcc 616 ! ------------------------------------------- 617 za = MIN( 1., zbetdo(ji,jj,jk+1), zbetup(ji,jj,jk) ) 618 zb = MIN( 1., zbetup(ji,jj,jk+1), zbetdo(ji,jj,jk) ) 619 zc = ( 0.5 + SIGN( 0.5 , pcc(ji,jj,jk+1) ) ) 620 pcc(ji,jj,jk+1) = pcc(ji,jj,jk+1) * ( zc * za + ( 1._wp - zc) * zb ) 621 #ifdef BULL_CB_WITH_BUF 622 ! zt3ew(:,jh,jk,jl,jf,1) = ARRAY_IN(nn_hls+jh,:,jk,jl,jf 623 buf(jj,1,jk,1,1,1) = paa(ji,jj,jk) 624 buf(jj,1,jk,1,2,1) = pbb(ji,jj,jk) 625 #endif 626 END DO 627 END DO 628 end do 629 else 630 DO jk = k0, k1 631 DO jj = j0, j1 632 !DIR$ vector always 633 DO ji = i0, i1 634 zau = MIN( 1._wp, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 635 zbu = MIN( 1._wp, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) 636 zcu = ( 0.5 + SIGN( 0.5 , paa(ji,jj,jk) ) ) 637 paa(ji,jj,jk) = paa(ji,jj,jk) * ( zcu * zau + ( 1._wp - zcu) * zbu ) 638 639 zav = MIN( 1._wp, zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) ) 640 zbv = MIN( 1._wp, zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) ) 641 zcv = ( 0.5 + SIGN( 0.5 , pbb(ji,jj,jk) ) ) 642 pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zcv * zav + ( 1._wp - zcv) * zbv ) 643 644 ! monotonic flux in the k direction, i.e. pcc 645 ! ------------------------------------------- 646 za = MIN( 1., zbetdo(ji,jj,jk+1), zbetup(ji,jj,jk) ) 647 zb = MIN( 1., zbetup(ji,jj,jk+1), zbetdo(ji,jj,jk) ) 648 zc = ( 0.5 + SIGN( 0.5 , pcc(ji,jj,jk+1) ) ) 649 pcc(ji,jj,jk+1) = pcc(ji,jj,jk+1) * ( zc * za + ( 1._wp - zc) * zb ) 650 END DO 651 END DO 652 END DO 653 endif 654 end subroutine 655 #endif 424 656 END SUBROUTINE nonosc 425 657
Note: See TracChangeset
for help on using the changeset viewer.