Changeset 181
- Timestamp:
- 2004-11-05T15:08:18+01:00 (19 years ago)
- Location:
- trunk/NEMO/OPA_SRC
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/DOM/domain.F90
r112 r181 344 344 IF(lwp)WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors' 345 345 IF(lwp)WRITE(numout,*) '~~~~~~~' 346 ze1min = MINVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 347 ze2min = MINVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 348 ze1max = MAXVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 349 ze2max = MAXVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 350 351 IF( lk_mpp ) CALL mpp_min( ze1min ) ! min over the global domain 352 IF( lk_mpp ) CALL mpp_min( ze2min ) 353 IF( lk_mpp ) CALL mpp_max( ze1max ) ! max over the global domain 354 IF( lk_mpp ) CALL mpp_max( ze2max ) 355 356 iloc = MINLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 357 iimi1 = iloc(1) + nimpp - 1 358 ijmi1 = iloc(2) + njmpp - 1 359 iloc = MINLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 360 iimi2 = iloc(1) + nimpp - 1 361 ijmi2 = iloc(2) + njmpp - 1 362 iloc = MAXLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 363 iima1 = iloc(1) + nimpp - 1 364 ijma1 = iloc(2) + njmpp - 1 365 iloc = MAXLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 366 iima2 = iloc(1) + nimpp - 1 367 ijma2 = iloc(2) + njmpp - 1 368 369 IF( lk_mpp ) THEN 370 !CT bug CALL mpp_isl( iimi1 ) 371 !CT bug CALL mpp_isl( ijmi1 ) 372 !CT bug CALL mpp_isl( iimi2 ) 373 !CT bug CALL mpp_isl( ijmi2 ) 374 !CT bug CALL mpp_isl( iima1 ) 375 !CT bug CALL mpp_isl( ijma1 ) 376 !CT bug CALL mpp_isl( iima2 ) 377 !CT bug CALL mpp_isl( ijma2 ) 378 ENDIF 379 380 IF(lwp) THEN 381 IF(lk_mpp) THEN 382 WRITE(numout,cform_war) 383 WRITE(numout,*)' Min(Max) of e1t, e2t are those of the first proc only' 384 WRITE(numout,*) 385 END IF 386 WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i3,' j= ',i3)") ze1max, iima1, ijma1 387 WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i3,' j= ',i3)") ze1min, iimi1, ijmi1 388 WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i3,' j= ',i3)") ze2max, iima2, ijma2 389 WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i3,' j= ',i3)") ze2min, iimi2, ijmi2 346 347 IF (lk_mpp) THEN 348 CALL mpp_minloc( e1t(:,:), tmask(:,:,1), ze1min, iimi1,ijmi1 ) 349 CALL mpp_minloc( e2t(:,:), tmask(:,:,1), ze2min, iimi2,ijmi2 ) 350 CALL mpp_maxloc( e1t(:,:), tmask(:,:,1), ze1max, iima1,ijma1 ) 351 CALL mpp_maxloc( e2t(:,:), tmask(:,:,1), ze2max, iima2,ijma2 ) 352 ELSE 353 ze1min = MINVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 354 ze2min = MINVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 355 ze1max = MAXVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 356 ze2max = MAXVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 357 358 iloc = MINLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 359 iimi1 = iloc(1) + nimpp - 1 360 ijmi1 = iloc(2) + njmpp - 1 361 iloc = MINLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 362 iimi2 = iloc(1) + nimpp - 1 363 ijmi2 = iloc(2) + njmpp - 1 364 iloc = MAXLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 365 iima1 = iloc(1) + nimpp - 1 366 ijma1 = iloc(2) + njmpp - 1 367 iloc = MAXLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 368 iima2 = iloc(1) + nimpp - 1 369 ijma2 = iloc(2) + njmpp - 1 370 ENDIF 371 372 IF(lwp) THEN 373 WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1 374 WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1 375 WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2 376 WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2 390 377 ENDIF 391 378 -
trunk/NEMO/OPA_SRC/lib_mpp.F90
r51 r181 67 67 MODULE PROCEDURE mpp_lbc_north_3d, mpp_lbc_north_2d 68 68 END INTERFACE 69 INTERFACE mpp_minloc 70 MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d 71 END INTERFACE 72 INTERFACE mpp_maxloc 73 MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d 74 END INTERFACE 75 69 76 70 77 !! * Share module variables 71 78 LOGICAL, PUBLIC, PARAMETER :: lk_mpp = .TRUE. !: mpp flag 79 LOGICAL, PUBLIC :: lk_bsend = .FALSE. !: mpp_bsend flag 80 LOGICAL, PUBLIC :: lk_isend = .FALSE. !: mpp_isend flag 72 81 73 82 … … 248 257 ! Enroll in MPI 249 258 ! ------------- 250 !!! CALL mpi_init_opa( ierr ) 251 CALL mpi_init( ierr ) 259 # if defined key_mpi_bsend 260 lk_bsend = .TRUE. !: mpp_bsend flag 261 # endif 262 # if defined key_mpi_isend 263 lk_isend = .TRUE. !: mpp_isend flag 264 # endif 265 266 IF(lk_bsend) THEN 267 CALL mpi_init_opa( ierr ) 268 ELSE 269 CALL mpi_init( ierr ) 270 ENDIF 252 271 CALL mpi_comm_rank( mpi_comm_world, rank, ierr ) 253 272 CALL mpi_comm_size( mpi_comm_world, size, ierr ) … … 504 523 INTEGER :: ji, jk, jl ! dummy loop indices 505 524 INTEGER :: imigr, iihom, ijhom, iloc, ijt, iju ! temporary integers 525 INTEGER :: ml_req1, ml_req2, ml_err ! for key_mpi_isend 526 INTEGER :: ml_stat(MPI_STATUS_SIZE) ! for key_mpi_isend 506 527 !!---------------------------------------------------------------------- 507 528 … … 577 598 SELECT CASE ( nbondi ) 578 599 CASE ( -1 ) 579 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea )600 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea, ml_req1 ) 580 601 CALL mpprecv( 1, t3ew(1,1,1,2), imigr ) 602 IF(lk_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err) 581 603 CASE ( 0 ) 582 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe )583 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea )604 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe, ml_req1 ) 605 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea, ml_req2 ) 584 606 CALL mpprecv( 1, t3ew(1,1,1,2), imigr ) 585 607 CALL mpprecv( 2, t3we(1,1,1,2), imigr ) 608 IF(lk_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err) 609 IF(lk_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err) 586 610 CASE ( 1 ) 587 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe )611 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe, ml_req1 ) 588 612 CALL mpprecv( 2, t3we(1,1,1,2), imigr ) 613 IF(lk_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err) 589 614 END SELECT 590 615 #endif … … 651 676 SELECT CASE ( nbondj ) 652 677 CASE ( -1 ) 653 CALL mppsend( 4, t3sn(1,1,1,1), imigr, nono )678 CALL mppsend( 4, t3sn(1,1,1,1), imigr, nono, ml_req1 ) 654 679 CALL mpprecv( 3, t3ns(1,1,1,2), imigr ) 680 IF(lk_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err) 655 681 CASE ( 0 ) 656 CALL mppsend( 3, t3ns(1,1,1,1), imigr, noso )657 CALL mppsend( 4, t3sn(1,1,1,1), imigr, nono )682 CALL mppsend( 3, t3ns(1,1,1,1), imigr, noso, ml_req1 ) 683 CALL mppsend( 4, t3sn(1,1,1,1), imigr, nono, ml_req2 ) 658 684 CALL mpprecv( 3, t3ns(1,1,1,2), imigr ) 659 685 CALL mpprecv( 4, t3sn(1,1,1,2), imigr ) 686 IF(lk_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err) 687 IF(lk_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err) 660 688 CASE ( 1 ) 661 CALL mppsend( 3, t3ns(1,1,1,1), imigr, noso )689 CALL mppsend( 3, t3ns(1,1,1,1), imigr, noso, ml_req1 ) 662 690 CALL mpprecv( 4, t3sn(1,1,1,2), imigr ) 691 IF(lk_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err) 663 692 END SELECT 664 693 … … 849 878 SELECT CASE ( nbondi ) 850 879 CASE ( -1 ) 851 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea )880 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea, ml_req1 ) 852 881 CALL mpprecv( 1, t3ew(1,1,1,2), imigr ) 882 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 853 883 CASE ( 0 ) 854 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe )855 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea )884 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe, ml_req1 ) 885 CALL mppsend( 2, t3we(1,1,1,1), imigr, noea, ml_req2 ) 856 886 CALL mpprecv( 1, t3ew(1,1,1,2), imigr ) 857 887 CALL mpprecv( 2, t3we(1,1,1,2), imigr ) 888 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 889 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 858 890 CASE ( 1 ) 859 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe )891 CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe, ml_req1 ) 860 892 CALL mpprecv( 2, t3we(1,1,1,2), imigr ) 893 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 861 894 END SELECT 862 895 #endif … … 925 958 imigr, iihom, ijhom, & ! temporary integers 926 959 iloc, ijt, iju ! " " 960 INTEGER :: ml_req1, ml_req2, ml_err ! for key_mpi_isend 961 INTEGER :: ml_stat(MPI_STATUS_SIZE) ! for key_mpi_isend 927 962 !!---------------------------------------------------------------------- 928 963 … … 999 1034 SELECT CASE ( nbondi ) 1000 1035 CASE ( -1 ) 1001 CALL mppsend( 2, t2we(1,1,1), imigr, noea )1036 CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req1 ) 1002 1037 CALL mpprecv( 1, t2ew(1,1,2), imigr ) 1038 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1003 1039 CASE ( 0 ) 1004 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe )1005 CALL mppsend( 2, t2we(1,1,1), imigr, noea )1040 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 ) 1041 CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req2 ) 1006 1042 CALL mpprecv( 1, t2ew(1,1,2), imigr ) 1007 1043 CALL mpprecv( 2, t2we(1,1,2), imigr ) 1044 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1045 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 1008 1046 CASE ( 1 ) 1009 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe )1047 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 ) 1010 1048 CALL mpprecv( 2, t2we(1,1,2), imigr ) 1049 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1011 1050 END SELECT 1012 1051 … … 1072 1111 SELECT CASE ( nbondj ) 1073 1112 CASE ( -1 ) 1074 CALL mppsend( 4, t2sn(1,1,1), imigr, nono )1113 CALL mppsend( 4, t2sn(1,1,1), imigr, nono, ml_req1 ) 1075 1114 CALL mpprecv( 3, t2ns(1,1,2), imigr ) 1115 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1076 1116 CASE ( 0 ) 1077 CALL mppsend( 3, t2ns(1,1,1), imigr, noso )1078 CALL mppsend( 4, t2sn(1,1,1), imigr, nono )1117 CALL mppsend( 3, t2ns(1,1,1), imigr, noso, ml_req1 ) 1118 CALL mppsend( 4, t2sn(1,1,1), imigr, nono, ml_req2 ) 1079 1119 CALL mpprecv( 3, t2ns(1,1,2), imigr ) 1080 1120 CALL mpprecv( 4, t2sn(1,1,2), imigr ) 1121 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1122 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 1081 1123 CASE ( 1 ) 1082 CALL mppsend( 3, t2ns(1,1,1), imigr, noso )1124 CALL mppsend( 3, t2ns(1,1,1), imigr, noso, ml_req1 ) 1083 1125 CALL mpprecv( 4, t2sn(1,1,2), imigr ) 1126 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1084 1127 END SELECT 1085 1128 … … 1268 1311 SELECT CASE ( nbondi ) 1269 1312 CASE ( -1 ) 1270 CALL mppsend( 2, t2we(1,1,1), imigr, noea )1313 CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req1 ) 1271 1314 CALL mpprecv( 1, t2ew(1,1,2), imigr ) 1315 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1272 1316 CASE ( 0 ) 1273 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe )1274 CALL mppsend( 2, t2we(1,1,1), imigr, noea )1317 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 ) 1318 CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req2 ) 1275 1319 CALL mpprecv( 1, t2ew(1,1,2), imigr ) 1276 1320 CALL mpprecv( 2, t2we(1,1,2), imigr ) 1321 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1322 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 1277 1323 CASE ( 1 ) 1278 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe )1324 CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 ) 1279 1325 CALL mpprecv( 2, t2we(1,1,2), imigr ) 1326 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1280 1327 END SELECT 1281 1328 #endif … … 1333 1380 INTEGER :: & 1334 1381 imigr, iihom, ijhom ! temporary integers 1382 INTEGER :: ml_req1, ml_req2, ml_err ! for key_mpi_isend 1383 INTEGER :: ml_stat(MPI_STATUS_SIZE) ! for key_mpi_isend 1335 1384 !!---------------------------------------------------------------------- 1336 1385 … … 1369 1418 imigr=jprecj*jpi 1370 1419 1371 CALL mppsend(3,t2p1(1,1,1),imigr,nono )1420 CALL mppsend(3,t2p1(1,1,1),imigr,nono, ml_req1) 1372 1421 CALL mpprecv(3,t2p1(1,1,2),imigr) 1422 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1373 1423 1374 1424 #endif … … 1397 1447 imigr=jprecj*jpi 1398 1448 1399 CALL mppsend(3,t2p1(1,1,1),imigr,nono )1449 CALL mppsend(3,t2p1(1,1,1),imigr,nono, ml_req1) 1400 1450 CALL mpprecv(3,t2p1(1,1,2),imigr) 1451 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1401 1452 1402 1453 #endif … … 1459 1510 1460 1511 CASE ( -1 ) 1461 CALL mppsend(2,t2we(1,1,1),imigr,noea )1512 CALL mppsend(2,t2we(1,1,1),imigr,noea, ml_req1) 1462 1513 CALL mpprecv(1,t2ew(1,1,2),imigr) 1463 1514 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1464 1515 CASE ( 0 ) 1465 CALL mppsend(1,t2ew(1,1,1),imigr,nowe )1466 CALL mppsend(2,t2we(1,1,1),imigr,noea )1516 CALL mppsend(1,t2ew(1,1,1),imigr,nowe, ml_req1) 1517 CALL mppsend(2,t2we(1,1,1),imigr,noea, ml_req2) 1467 1518 CALL mpprecv(1,t2ew(1,1,2),imigr) 1468 1519 CALL mpprecv(2,t2we(1,1,2),imigr) 1520 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1521 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 1469 1522 1470 1523 CASE ( 1 ) 1471 CALL mppsend(1,t2ew(1,1,1),imigr,nowe )1524 CALL mppsend(1,t2ew(1,1,1),imigr,nowe, ml_req1) 1472 1525 CALL mpprecv(2,t2we(1,1,2),imigr) 1526 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1473 1527 1474 1528 END SELECT … … 1549 1603 1550 1604 CASE ( -1 ) 1551 CALL mppsend(4,t2sn(1,1,1),imigr,nono )1605 CALL mppsend(4,t2sn(1,1,1),imigr,nono, ml_req1) 1552 1606 CALL mpprecv(3,t2ns(1,1,2),imigr) 1607 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1553 1608 1554 1609 CASE ( 0 ) 1555 CALL mppsend(3,t2ns(1,1,1),imigr,noso )1556 CALL mppsend(4,t2sn(1,1,1),imigr,nono )1610 CALL mppsend(3,t2ns(1,1,1),imigr,noso, ml_req1) 1611 CALL mppsend(4,t2sn(1,1,1),imigr,nono, ml_req2) 1557 1612 CALL mpprecv(3,t2ns(1,1,2),imigr) 1558 1613 CALL mpprecv(4,t2sn(1,1,2),imigr) 1614 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1615 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 1559 1616 1560 1617 CASE ( 1 ) 1561 CALL mppsend(3,t2ns(1,1,1),imigr,noso )1618 CALL mppsend(3,t2ns(1,1,1),imigr,noso, ml_req1) 1562 1619 CALL mpprecv(4,t2sn(1,1,2),imigr) 1620 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 1563 1621 END SELECT 1564 1622 … … 1592 1650 1593 1651 1594 SUBROUTINE mppsend( ktyp, pmess, kbytes, kdest 1652 SUBROUTINE mppsend( ktyp, pmess, kbytes, kdest, md_req) 1595 1653 !!---------------------------------------------------------------------- 1596 1654 !! *** routine mppsend *** … … 1603 1661 INTEGER , INTENT( in ) :: kbytes, & ! size of the array pmess 1604 1662 & kdest , & ! receive process number 1605 & ktyp ! Tag of the message 1663 & ktyp, & ! Tag of the message 1664 & md_req ! Argument for isend 1606 1665 !!---------------------------------------------------------------------- 1607 1666 #if defined key_mpp_shmem … … 1612 1671 INTEGER :: iflag 1613 1672 1614 CALL mpi_send( pmess, kbytes, mpi_real8, kdest, ktyp, & 1673 IF(lk_bsend) THEN 1674 CALL mpi_bsend( pmess, kbytes, mpi_double_precision, kdest, ktyp, & 1615 1675 & mpi_comm_world, iflag ) 1676 ELSEIF (lk_isend) THEN 1677 ! Carefull here : one more argument for mpi_isend : the mpi request identifier.. 1678 CALL mpi_isend( pmess, kbytes, mpi_double_precision, kdest, ktyp, & 1679 & mpi_comm_world, md_req, iflag ) 1680 ELSE 1681 CALL mpi_send( pmess, kbytes, mpi_double_precision, kdest, ktyp, & 1682 & mpi_comm_world, iflag ) 1683 ENDIF 1616 1684 #endif 1617 1685 … … 1639 1707 INTEGER :: iflag 1640 1708 1641 CALL mpi_recv( pmess, kbytes, mpi_ real8, mpi_any_source, ktyp, &1709 CALL mpi_recv( pmess, kbytes, mpi_double_precision, mpi_any_source, ktyp, & 1642 1710 & mpi_comm_world, istatus, iflag ) 1643 1711 #endif … … 1673 1741 1674 1742 itaille=jpi*jpj 1675 CALL mpi_gather( ptab, itaille, mpi_ real8, pio, itaille, &1676 & mpi_ real8, kp , mpi_comm_world, ierror )1743 CALL mpi_gather( ptab, itaille, mpi_double_precision, pio, itaille, & 1744 & mpi_double_precision, kp , mpi_comm_world, ierror ) 1677 1745 #endif 1678 1746 … … 1706 1774 1707 1775 itaille=jpi*jpj 1708 1709 CALL mpi_scatter( pio, itaille, mpi_ real8, ptab, itaille, &1710 & mpi_ real8, kp, mpi_comm_world, ierror )1776 1777 CALL mpi_scatter( pio, itaille, mpi_double_precision, ptab, itaille, & 1778 & mpi_double_precision, kp, mpi_comm_world, ierror ) 1711 1779 #endif 1712 1780 … … 2098 2166 2099 2167 CALL mpi_op_create(lc_isl,lcommute,mpi_isl,ierror) 2100 CALL mpi_allreduce(ptab, zwork,kdim,mpi_ real8&2168 CALL mpi_allreduce(ptab, zwork,kdim,mpi_double_precision & 2101 2169 ,mpi_isl,mpi_comm_world,ierror) 2102 2170 ptab(:) = zwork(:) … … 2156 2224 2157 2225 CALL mpi_op_create( lc_isl, lcommute, mpi_isl, ierror ) 2158 CALL mpi_allreduce( ptab, zwork, 1, mpi_ real8, &2226 CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, & 2159 2227 & mpi_isl , mpi_comm_world, ierror ) 2160 2228 ptab = zwork … … 2225 2293 REAL(wp), DIMENSION(kdim) :: zwork 2226 2294 2227 CALL mpi_allreduce(ptab, zwork,kdim,mpi_ real8&2295 CALL mpi_allreduce(ptab, zwork,kdim,mpi_double_precision & 2228 2296 ,mpi_max,mpi_comm_world,ierror) 2229 2297 ptab(:) = zwork(:) … … 2269 2337 REAL(wp) :: zwork 2270 2338 2271 CALL mpi_allreduce( ptab, zwork , 1 , mpi_ real8, &2339 CALL mpi_allreduce( ptab, zwork , 1 , mpi_double_precision, & 2272 2340 & mpi_max, mpi_comm_world, ierror ) 2273 2341 ptab = zwork … … 2325 2393 REAL(wp), DIMENSION(kdim) :: zwork 2326 2394 2327 CALL mpi_allreduce(ptab, zwork,kdim,mpi_ real8&2395 CALL mpi_allreduce(ptab, zwork,kdim,mpi_double_precision & 2328 2396 ,mpi_min,mpi_comm_world,ierror) 2329 2397 ptab(:) = zwork(:) … … 2370 2438 REAL(wp) :: zwork 2371 2439 2372 CALL mpi_allreduce( ptab, zwork, 1,mpi_ real8&2440 CALL mpi_allreduce( ptab, zwork, 1,mpi_double_precision & 2373 2441 & ,mpi_min,mpi_comm_world,ierror) 2374 2442 ptab = zwork … … 2426 2494 REAL(wp), DIMENSION(kdim) :: zwork ! temporary workspace 2427 2495 2428 CALL mpi_allreduce(ptab, zwork,kdim,mpi_ real8&2496 CALL mpi_allreduce(ptab, zwork,kdim,mpi_double_precision & 2429 2497 & ,mpi_sum,mpi_comm_world,ierror) 2430 2498 ptab(:) = zwork(:) … … 2470 2538 REAL(wp) :: zwork 2471 2539 2472 CALL mpi_allreduce(ptab, zwork, 1,mpi_ real8&2540 CALL mpi_allreduce(ptab, zwork, 1,mpi_double_precision & 2473 2541 & ,mpi_sum,mpi_comm_world,ierror) 2474 2542 ptab = zwork … … 2478 2546 END SUBROUTINE mppsum_real 2479 2547 2548 SUBROUTINE mpp_minloc2d(ptab, pmask, pmin, ki,kj ) 2549 !!------------------------------------------------------------------------ 2550 !! *** routine mpp_minloc *** 2551 !! 2552 !! ** Purpose : Compute the global minimum of an array ptab 2553 !! and also give its global position 2554 !! 2555 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 2556 !! 2557 !! ** Arguments : I : ptab =local 2D array 2558 !! O : pmin = global minimum 2559 !! O : ki,kj = global position of minimum 2560 !! 2561 !! ** Author : J.M. Molines 10/10/2004 2562 !!-------------------------------------------------------------------------- 2563 #ifdef key_mpp_shmem 2564 IF (lwp) THEN 2565 WRITE(numout,*) ' mpp_minloc not yet available in SHMEM' 2566 STOP 2567 ENDIF 2568 # elif key_mpp_mpi 2569 !! * Arguments 2570 REAL(wp), DIMENSION (jpi,jpj), INTENT (in) :: ptab ,& ! Local 2D array 2571 & pmask ! Local mask 2572 REAL(wp) , INTENT (out) :: pmin ! Global minimum of ptab 2573 INTEGER , INTENT (out) :: ki,kj ! index of minimum in global frame 2574 2575 !! * Local variables 2576 REAL(wp) :: zmin ! local minimum 2577 REAL(wp) ,DIMENSION(2,1) :: zain, zaout 2578 INTEGER, DIMENSION (2) :: ilocs 2579 INTEGER :: ierror 2580 2581 2582 zmin = MINVAL( ptab(:,:) , mask= pmask == 1.e0 ) 2583 ilocs = MINLOC( ptab(:,:) , mask= pmask == 1.e0 ) 2584 2585 ki = ilocs(1) + nimpp - 1 2586 kj = ilocs(2) + njmpp - 1 2587 2588 zain(1,:)=zmin 2589 zain(2,:)=ki+10000.*kj 2590 2591 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_WORLD,ierror) 2592 2593 pmin=zaout(1,1) 2594 kj= INT(zaout(2,1)/10000.) 2595 ki= INT(zaout(2,1) - 10000.*kj ) 2596 #endif 2597 2598 END SUBROUTINE mpp_minloc2d 2599 2600 2601 SUBROUTINE mpp_minloc3d(ptab, pmask, pmin, ki,kj ,kk) 2602 !!------------------------------------------------------------------------ 2603 !! *** routine mpp_minloc *** 2604 !! 2605 !! ** Purpose : Compute the global minimum of an array ptab 2606 !! and also give its global position 2607 !! 2608 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 2609 !! 2610 !! ** Arguments : I : ptab =local 2D array 2611 !! O : pmin = global minimum 2612 !! O : ki,kj = global position of minimum 2613 !! 2614 !! ** Author : J.M. Molines 10/10/2004 2615 !!-------------------------------------------------------------------------- 2616 #ifdef key_mpp_shmem 2617 IF (lwp) THEN 2618 WRITE(numout,*) ' mpp_minloc not yet available in SHMEM' 2619 STOP 2620 ENDIF 2621 # elif key_mpp_mpi 2622 !! * Arguments 2623 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT (in) :: ptab ,& ! Local 2D array 2624 & pmask ! Local mask 2625 REAL(wp) , INTENT (out) :: pmin ! Global minimum of ptab 2626 INTEGER , INTENT (out) :: ki,kj,kk ! index of minimum in global frame 2627 2628 !! * Local variables 2629 REAL(wp) :: zmin ! local minimum 2630 REAL(wp) ,DIMENSION(2,1) :: zain, zaout 2631 INTEGER, DIMENSION (3) :: ilocs 2632 INTEGER :: ierror 2633 2634 2635 zmin = MINVAL( ptab(:,:,:) , mask= pmask == 1.e0 ) 2636 ilocs = MINLOC( ptab(:,:,:) , mask= pmask == 1.e0 ) 2637 2638 ki = ilocs(1) + nimpp - 1 2639 kj = ilocs(2) + njmpp - 1 2640 kk = ilocs(3) 2641 2642 zain(1,:)=zmin 2643 zain(2,:)=ki+10000.*kj+100000000.*kk 2644 2645 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_WORLD,ierror) 2646 2647 pmin=zaout(1,1) 2648 kk= INT(zaout(2,1)/100000000.) 2649 kj= INT(zaout(2,1) - kk * 100000000. )/10000 2650 ki= INT(zaout(2,1) - kk * 100000000. -kj * 10000. ) 2651 #endif 2652 2653 END SUBROUTINE mpp_minloc3d 2654 2655 2656 SUBROUTINE mpp_maxloc2d(ptab, pmask, pmax, ki,kj ) 2657 !!------------------------------------------------------------------------ 2658 !! *** routine mpp_maxloc *** 2659 !! 2660 !! ** Purpose : Compute the global maximum of an array ptab 2661 !! and also give its global position 2662 !! 2663 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 2664 !! 2665 !! ** Arguments : I : ptab =local 2D array 2666 !! O : pmax = global maximum 2667 !! O : ki,kj = global position of maximum 2668 !! 2669 !! ** Author : J.M. Molines 10/10/2004 2670 !!-------------------------------------------------------------------------- 2671 #ifdef key_mpp_shmem 2672 IF (lwp) THEN 2673 WRITE(numout,*) ' mpp_maxloc not yet available in SHMEM' 2674 STOP 2675 ENDIF 2676 # elif key_mpp_mpi 2677 !! * Arguments 2678 REAL(wp), DIMENSION (jpi,jpj), INTENT (in) :: ptab ,& ! Local 2D array 2679 & pmask ! Local mask 2680 REAL(wp) , INTENT (out) :: pmax ! Global maximum of ptab 2681 INTEGER , INTENT (out) :: ki,kj ! index of maximum in global frame 2682 2683 !! * Local variables 2684 REAL(wp) :: zmax ! local maximum 2685 REAL(wp) ,DIMENSION(2,1) :: zain, zaout 2686 INTEGER, DIMENSION (2) :: ilocs 2687 INTEGER :: ierror 2688 2689 2690 zmax = MAXVAL( ptab(:,:) , mask= pmask == 1.e0 ) 2691 ilocs = MAXLOC( ptab(:,:) , mask= pmask == 1.e0 ) 2692 2693 ki = ilocs(1) + nimpp - 1 2694 kj = ilocs(2) + njmpp - 1 2695 2696 zain(1,:)=zmax 2697 zain(2,:)=ki+10000.*kj 2698 2699 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_WORLD,ierror) 2700 2701 pmax=zaout(1,1) 2702 kj= INT(zaout(2,1)/10000.) 2703 ki= INT(zaout(2,1) - 10000.*kj ) 2704 #endif 2705 2706 END SUBROUTINE mpp_maxloc2d 2707 2708 SUBROUTINE mpp_maxloc3d(ptab, pmask, pmax, ki,kj,kk ) 2709 !!------------------------------------------------------------------------ 2710 !! *** routine mpp_maxloc *** 2711 !! 2712 !! ** Purpose : Compute the global maximum of an array ptab 2713 !! and also give its global position 2714 !! 2715 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 2716 !! 2717 !! ** Arguments : I : ptab =local 2D array 2718 !! O : pmax = global maximum 2719 !! O : ki,kj = global position of maximum 2720 !! 2721 !! ** Author : J.M. Molines 10/10/2004 2722 !!-------------------------------------------------------------------------- 2723 #ifdef key_mpp_shmem 2724 IF (lwp) THEN 2725 WRITE(numout,*) ' mpp_maxloc not yet available in SHMEM' 2726 STOP 2727 ENDIF 2728 # elif key_mpp_mpi 2729 !! * Arguments 2730 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT (in) :: ptab ,& ! Local 2D array 2731 & pmask ! Local mask 2732 REAL(wp) , INTENT (out) :: pmax ! Global maximum of ptab 2733 INTEGER , INTENT (out) :: ki,kj,kk ! index of maximum in global frame 2734 2735 !! * Local variables 2736 REAL(wp) :: zmax ! local maximum 2737 REAL(wp) ,DIMENSION(2,1) :: zain, zaout 2738 INTEGER, DIMENSION (3) :: ilocs 2739 INTEGER :: ierror 2740 2741 2742 zmax = MAXVAL( ptab(:,:,:) , mask= pmask == 1.e0 ) 2743 ilocs = MAXLOC( ptab(:,:,:) , mask= pmask == 1.e0 ) 2744 2745 ki = ilocs(1) + nimpp - 1 2746 kj = ilocs(2) + njmpp - 1 2747 kk = ilocs(3) 2748 2749 zain(1,:)=zmax 2750 zain(2,:)=ki+10000.*kj+100000000.*kk 2751 2752 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_WORLD,ierror) 2753 2754 pmax=zaout(1,1) 2755 kk= INT(zaout(2,1)/100000000.) 2756 kj= INT(zaout(2,1) - kk * 100000000. )/10000 2757 ki= INT(zaout(2,1) - kk * 100000000. -kj * 10000. ) 2758 #endif 2759 2760 END SUBROUTINE mpp_maxloc3d 2480 2761 2481 2762 SUBROUTINE mppsync() … … 2598 2879 ijpt0, ijpt1, & ! " " 2599 2880 imigr, iihom, ijhom ! " " 2881 INTEGER :: ml_req1, ml_req2, ml_err ! for key_mpi_isend 2882 INTEGER :: ml_stat(MPI_STATUS_SIZE) ! for key_mpi_isend 2600 2883 REAL(wp), DIMENSION(jpi,jpj) :: & 2601 2884 ztab ! temporary workspace … … 2678 2961 2679 2962 IF( nbondi == -1 ) THEN 2680 CALL mppsend(2,t2we(1,1,1),imigr,noea )2963 CALL mppsend(2,t2we(1,1,1),imigr,noea, ml_req1) 2681 2964 CALL mpprecv(1,t2ew(1,1,2),imigr) 2965 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 2682 2966 ELSEIF( nbondi == 0 ) THEN 2683 CALL mppsend(1,t2ew(1,1,1),imigr,nowe )2684 CALL mppsend(2,t2we(1,1,1),imigr,noea )2967 CALL mppsend(1,t2ew(1,1,1),imigr,nowe, ml_req1) 2968 CALL mppsend(2,t2we(1,1,1),imigr,noea, ml_req2) 2685 2969 CALL mpprecv(1,t2ew(1,1,2),imigr) 2686 2970 CALL mpprecv(2,t2we(1,1,2),imigr) 2971 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 2972 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 2687 2973 ELSEIF( nbondi == 1 ) THEN 2688 CALL mppsend(1,t2ew(1,1,1),imigr,nowe )2974 CALL mppsend(1,t2ew(1,1,1),imigr,nowe, ml_req1) 2689 2975 CALL mpprecv(2,t2we(1,1,2),imigr) 2976 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 2690 2977 ENDIF 2691 2978 #endif … … 2745 3032 2746 3033 IF( nbondj == -1 ) THEN 2747 CALL mppsend(4,t2sn(1,1,1),imigr,nono )3034 CALL mppsend(4,t2sn(1,1,1),imigr,nono, ml_req1) 2748 3035 CALL mpprecv(3,t2ns(1,1,2),imigr) 3036 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 2749 3037 ELSEIF( nbondj == 0 ) THEN 2750 CALL mppsend(3,t2ns(1,1,1),imigr,noso )2751 CALL mppsend(4,t2sn(1,1,1),imigr,nono )3038 CALL mppsend(3,t2ns(1,1,1),imigr,noso, ml_req1) 3039 CALL mppsend(4,t2sn(1,1,1),imigr,nono, ml_req2) 2752 3040 CALL mpprecv(3,t2ns(1,1,2),imigr) 2753 3041 CALL mpprecv(4,t2sn(1,1,2),imigr) 3042 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 3043 IF(lk_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err) 2754 3044 ELSEIF( nbondj == 1 ) THEN 2755 CALL mppsend(3,t2ns(1,1,1),imigr,noso )3045 CALL mppsend(3,t2ns(1,1,1),imigr,noso, ml_req1) 2756 3046 CALL mpprecv(4,t2sn(1,1,2),imigr) 3047 IF(lk_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err) 2757 3048 ENDIF 2758 3049 … … 2941 3232 #elif defined key_mpp_mpi 2942 3233 itaille=jpi*jpk*ijpj 2943 CALL MPI_GATHER(znorthloc,itaille,MPI_ REAL8,znorthgloio,itaille,MPI_REAL8,0,ncomm_north,ierr)3234 CALL MPI_GATHER(znorthloc,itaille,MPI_DOUBLE_PRECISION,znorthgloio,itaille,MPI_DOUBLE_PRECISION,0,ncomm_north,ierr) 2944 3235 #endif 2945 3236 … … 3102 3393 IF ( npolj /= 0 ) THEN 3103 3394 itaille=jpi*jpk*ijpj 3104 CALL MPI_SCATTER(znorthgloio,itaille,MPI_ REAL8,znorthloc,itaille,MPI_REAL8,0,ncomm_north,ierr)3395 CALL MPI_SCATTER(znorthgloio,itaille,MPI_DOUBLE_PRECISION,znorthloc,itaille,MPI_DOUBLE_PRECISION,0,ncomm_north,ierr) 3105 3396 ENDIF 3106 3397 #endif … … 3185 3476 #elif defined key_mpp_mpi 3186 3477 itaille=jpi*ijpj 3187 CALL MPI_GATHER(znorthloc,itaille,MPI_ REAL8,znorthgloio,itaille,MPI_REAL8,0,ncomm_north,ierr)3478 CALL MPI_GATHER(znorthloc,itaille,MPI_DOUBLE_PRECISION,znorthgloio,itaille,MPI_DOUBLE_PRECISION,0,ncomm_north,ierr) 3188 3479 #endif 3189 3480 ENDIF … … 3344 3635 IF ( npolj /= 0 ) THEN 3345 3636 itaille=jpi*ijpj 3346 CALL MPI_SCATTER(znorthgloio,itaille,MPI_ REAL8,znorthloc,itaille,MPI_REAL8,0,ncomm_north,ierr)3637 CALL MPI_SCATTER(znorthgloio,itaille,MPI_DOUBLE_PRECISION,znorthloc,itaille,MPI_DOUBLE_PRECISION,0,ncomm_north,ierr) 3347 3638 ENDIF 3348 3639 #endif … … 3448 3739 MODULE PROCEDURE mppobc_1d, mppobc_2d, mppobc_3d, mppobc_4d 3449 3740 END INTERFACE 3741 INTERFACE mpp_minloc 3742 MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d 3743 END INTERFACE 3744 INTERFACE mpp_maxloc 3745 MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d 3746 END INTERFACE 3747 3450 3748 3451 3749 LOGICAL, PUBLIC, PARAMETER :: lk_mpp = .FALSE. !: mpp flag … … 3577 3875 END SUBROUTINE mppisl_real 3578 3876 3877 SUBROUTINE mpp_minloc2d ( ptab, pmask, pmin, ki, kj ) 3878 REAL :: pmin 3879 REAL , DIMENSION (:,:) :: ptab, pmask 3880 INTEGER :: ki, kj 3881 WRITE(*,*) 'mppisl_real: You should not have seen this print! error?', pmin, ki, kj 3882 WRITE(*,*) ' " ": " " ', ptab(1,1), pmask(1,1) 3883 END SUBROUTINE mpp_minloc2d 3884 3885 SUBROUTINE mpp_minloc3d ( ptab, pmask, pmin, ki, kj, kk ) 3886 REAL :: pmin 3887 REAL , DIMENSION (:,:,:) :: ptab, pmask 3888 INTEGER :: ki, kj, kk 3889 WRITE(*,*) 'mppisl_real: You should not have seen this print! error?', pmin, ki, kj, kk 3890 WRITE(*,*) ' " ": " " ', ptab(1,1,1), pmask(1,1,1) 3891 END SUBROUTINE mpp_minloc3d 3892 3893 SUBROUTINE mpp_maxloc2d ( ptab, pmask, pmax, ki, kj ) 3894 REAL :: pmax 3895 REAL , DIMENSION (:,:) :: ptab, pmask 3896 INTEGER :: ki, kj 3897 WRITE(*,*) 'mppisl_real: You should not have seen this print! error?', pmax, ki, kj 3898 WRITE(*,*) ' " ": " " ', ptab(1,1), pmask(1,1) 3899 END SUBROUTINE mpp_maxloc2d 3900 3901 SUBROUTINE mpp_maxloc3d ( ptab, pmask, pmax, ki, kj, kk ) 3902 REAL :: pmax 3903 REAL , DIMENSION (:,:,:) :: ptab, pmask 3904 INTEGER :: ki, kj, kk 3905 WRITE(*,*) 'mppisl_real: You should not have seen this print! error?', pmax, ki, kj, kk 3906 WRITE(*,*) ' " ": " " ', ptab(1,1,1), pmask(1,1,1) 3907 END SUBROUTINE mpp_maxloc3d 3908 3579 3909 SUBROUTINE mppstop 3580 3910 WRITE(*,*) 'mppstop: You should not have seen this print! error?' -
trunk/NEMO/OPA_SRC/stpctl.F90
r79 r181 117 117 ENDIF 118 118 IF( zumax > 20.) THEN 119 ilocu = MAXLOC( ABS( un(:,:,:) ) ) 120 ii = ilocu(1) + nimpp - 1 121 ij = ilocu(2) + njmpp - 1 122 ik = ilocu(3) 123 IF( lk_mpp ) THEN 124 CALL mpp_isl( ii ) 125 CALL mpp_isl( ij ) 126 CALL mpp_isl( ik ) 119 IF (lk_mpp ) THEN 120 CALL mpp_maxloc(un,umask,zumax,ii,ij,ik) 121 ELSE 122 ilocu = MAXLOC( ABS( un(:,:,:) ) ) 123 ii = ilocu(1) + nimpp - 1 124 ij = ilocu(2) + njmpp - 1 125 ik = ilocu(3) 127 126 ENDIF 128 127 IF(lwp) THEN … … 156 155 ENDIF 157 156 IF( zsmin < 0.) THEN 158 ilocs = MINLOC( sn(:,:,1), mask = tmask(:,:,1) == 1.e0 ) 159 ii = ilocs(1) + nimpp - 1 160 ij = ilocs(2) + njmpp - 1 161 IF( lk_mpp ) CALL mpp_isl( ii ) 162 IF( lk_mpp ) CALL mpp_isl( ij ) 157 IF (lk_mpp) THEN 158 CALL mpp_minloc ( sn(:,:,1),tmask(:,:,1), zsmin, ii,ij ) 159 ELSE 160 ilocs = MINLOC( sn(:,:,1), mask = tmask(:,:,1) == 1.e0 ) 161 ii = ilocs(1) + nimpp - 1 162 ij = ilocs(2) + njmpp - 1 163 END IF 163 164 164 165 IF(lwp) THEN
Note: See TracChangeset
for help on using the changeset viewer.