Changeset 10425 for NEMO/trunk/src/OCE/LBC/lib_mpp.F90
- Timestamp:
- 2018-12-19T22:54:16+01:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/LBC/lib_mpp.F90
r10068 r10425 64 64 65 65 INTERFACE mpp_nfd 66 MODULE PROCEDURE mpp_nfd_2d , mpp_nfd_3d, mpp_nfd_4d66 MODULE PROCEDURE mpp_nfd_2d , mpp_nfd_3d , mpp_nfd_4d 67 67 MODULE PROCEDURE mpp_nfd_2d_ptr, mpp_nfd_3d_ptr, mpp_nfd_4d_ptr 68 68 END INTERFACE 69 69 70 70 ! Interface associated to the mpp_lnk_... routines is defined in lbclnk 71 PUBLIC mpp_lnk_2d , mpp_lnk_3d, mpp_lnk_4d71 PUBLIC mpp_lnk_2d , mpp_lnk_3d , mpp_lnk_4d 72 72 PUBLIC mpp_lnk_2d_ptr, mpp_lnk_3d_ptr, mpp_lnk_4d_ptr 73 73 ! … … 83 83 PUBLIC mpp_lbc_north_icb 84 84 PUBLIC mpp_min, mpp_max, mpp_sum, mpp_minloc, mpp_maxloc 85 PUBLIC mpp_ max_multiple85 PUBLIC mpp_delay_max, mpp_delay_sum, mpp_delay_rcv 86 86 PUBLIC mppscatter, mppgather 87 PUBLIC mpp_ini_ice, mpp_ini_znl 88 PUBLIC mppsize 87 PUBLIC mpp_ini_znl 89 88 PUBLIC mppsend, mpprecv ! needed by TAM and ICB routines 90 89 PUBLIC mpp_lnk_bdy_2d, mpp_lnk_bdy_3d, mpp_lnk_bdy_4d 91 PUBLIC mpprank92 90 93 91 !! * Interfaces … … 111 109 MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d 112 110 END INTERFACE 113 INTERFACE mpp_max_multiple114 MODULE PROCEDURE mppmax_real_multiple115 END INTERFACE116 111 117 112 !! ========================= !! … … 126 121 INTEGER, PARAMETER :: nprocmax = 2**10 ! maximun dimension (required to be a power of 2) 127 122 128 INTEGER :: mppsize ! number of process129 INTEGER :: mpprank ! process number [ 0 - size-1 ]123 INTEGER, PUBLIC :: mppsize ! number of process 124 INTEGER, PUBLIC :: mpprank ! process number [ 0 - size-1 ] 130 125 !$AGRIF_DO_NOT_TREAT 131 126 INTEGER, PUBLIC :: mpi_comm_oce ! opa local communicator … … 133 128 134 129 INTEGER :: MPI_SUMDD 135 136 ! variables used in case of sea-ice137 INTEGER, PUBLIC :: ncomm_ice !: communicator made by the processors with sea-ice (public so that it can be freed in icethd)138 INTEGER :: ngrp_iworld ! group ID for the world processors (for rheology)139 INTEGER :: ngrp_ice ! group ID for the ice processors (for rheology)140 INTEGER :: ndim_rank_ice ! number of 'ice' processors141 INTEGER :: n_ice_root ! number (in the comm_ice) of proc 0 in the ice comm142 INTEGER, DIMENSION(:), ALLOCATABLE, SAVE :: nrank_ice ! dimension ndim_rank_ice143 130 144 131 ! variables used for zonal integration … … 164 151 INTEGER , PUBLIC :: nn_buffer !: size of the buffer in case of mpi_bsend 165 152 153 ! Communications summary report 154 CHARACTER(len=128), DIMENSION(:), ALLOCATABLE :: crname_lbc !: names of lbc_lnk calling routines 155 CHARACTER(len=128), DIMENSION(:), ALLOCATABLE :: crname_glb !: names of global comm calling routines 156 INTEGER, PUBLIC :: ncom_stp = 0 !: copy of time step # istp 157 INTEGER, PUBLIC :: ncom_fsbc = 1 !: copy of sbc time step # nn_fsbc 158 INTEGER, PUBLIC :: ncom_dttrc = 1 !: copy of top time step # nn_dttrc 159 INTEGER, PUBLIC :: ncom_freq !: frequency of comm diagnostic 160 INTEGER, PUBLIC , DIMENSION(:,:), ALLOCATABLE :: ncomm_sequence !: size of communicated arrays (halos) 161 INTEGER, PARAMETER, PUBLIC :: ncom_rec_max = 2000 !: max number of communication record 162 INTEGER, PUBLIC :: n_sequence_lbc = 0 !: # of communicated arraysvia lbc 163 INTEGER, PUBLIC :: n_sequence_glb = 0 !: # of global communications 164 INTEGER, PUBLIC :: numcom = -1 !: logical unit for communicaton report 165 LOGICAL, PUBLIC :: l_full_nf_update = .TRUE. !: logical for a full (2lines) update of bc at North fold report 166 INTEGER, PARAMETER, PUBLIC :: nbdelay = 2 !: number of delayed operations 167 !: name (used as id) of allreduce-delayed operations 168 CHARACTER(len=32), DIMENSION(nbdelay), PUBLIC :: c_delaylist = (/ 'cflice', 'fwb' /) 169 !: component name where the allreduce-delayed operation is performed 170 CHARACTER(len=3), DIMENSION(nbdelay), PUBLIC :: c_delaycpnt = (/ 'ICE' , 'OCE' /) 171 TYPE, PUBLIC :: DELAYARR 172 REAL( wp), POINTER, DIMENSION(:) :: z1d => NULL() 173 COMPLEX(wp), POINTER, DIMENSION(:) :: y1d => NULL() 174 END TYPE DELAYARR 175 TYPE( DELAYARR ), DIMENSION(nbdelay), PUBLIC :: todelay 176 INTEGER, DIMENSION(nbdelay), PUBLIC :: ndelayid = -1 !: mpi request id of the delayed operations 177 178 ! timing summary report 179 REAL(wp), DIMENSION(2), PUBLIC :: waiting_time = 0._wp 180 REAL(wp) , PUBLIC :: compute_time = 0._wp, elapsed_time = 0._wp 181 166 182 REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE :: tampon ! buffer in case of bsend 167 183 … … 214 230 WRITE(ldtxt(ii),*) ' size exported buffer nn_buffer = ', nn_buffer,' bytes'; ii = ii + 1 215 231 ! 216 #if defined key_agrif217 IF( .NOT. Agrif_Root() ) THEN218 jpni = Agrif_Parent(jpni )219 jpnj = Agrif_Parent(jpnj )220 jpnij = Agrif_Parent(jpnij)221 ENDIF222 #endif223 !224 IF( jpnij < 1 ) THEN ! If jpnij is not specified in namelist then we calculate it225 jpnij = jpni * jpnj ! this means there will be no land cutting out.226 ENDIF227 228 232 IF( jpni < 1 .OR. jpnj < 1 ) THEN 229 233 WRITE(ldtxt(ii),*) ' jpni, jpnj and jpnij will be calculated automatically' ; ii = ii + 1 … … 231 235 WRITE(ldtxt(ii),*) ' processor grid extent in i jpni = ',jpni ; ii = ii + 1 232 236 WRITE(ldtxt(ii),*) ' processor grid extent in j jpnj = ',jpnj ; ii = ii + 1 233 WRITE(ldtxt(ii),*) ' number of local domains jpnij = ',jpnij ; ii = ii + 1234 237 ENDIF 235 238 … … 264 267 ! 265 268 ELSEIF ( PRESENT(localComm) .AND. .NOT. mpi_was_called ) THEN 269 WRITE(ldtxt(ii),cform_err) ; ii = ii + 1 266 270 WRITE(ldtxt(ii),*) ' lib_mpp: You cannot provide a local communicator ' ; ii = ii + 1 267 271 WRITE(ldtxt(ii),*) ' without calling MPI_Init before ! ' ; ii = ii + 1 … … 574 578 END SUBROUTINE mppscatter 575 579 580 581 SUBROUTINE mpp_delay_sum( cdname, cdelay, y_in, pout, ldlast, kcom ) 582 !!---------------------------------------------------------------------- 583 !! *** routine mpp_delay_sum *** 584 !! 585 !! ** Purpose : performed delayed mpp_sum, the result is received on next call 586 !! 587 !!---------------------------------------------------------------------- 588 CHARACTER(len=*), INTENT(in ) :: cdname ! name of the calling subroutine 589 CHARACTER(len=*), INTENT(in ) :: cdelay ! name (used as id) of the delayed operation 590 COMPLEX(wp), INTENT(in ), DIMENSION(:) :: y_in 591 REAL(wp), INTENT( out), DIMENSION(:) :: pout 592 LOGICAL, INTENT(in ) :: ldlast ! true if this is the last time we call this routine 593 INTEGER, INTENT(in ), OPTIONAL :: kcom 594 !! 595 INTEGER :: ji, isz 596 INTEGER :: idvar 597 INTEGER :: ierr, ilocalcomm 598 COMPLEX(wp), ALLOCATABLE, DIMENSION(:) :: ytmp 599 !!---------------------------------------------------------------------- 600 ilocalcomm = mpi_comm_oce 601 IF( PRESENT(kcom) ) ilocalcomm = kcom 602 603 isz = SIZE(y_in) 604 605 IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, ld_glb = .TRUE. ) 606 607 idvar = -1 608 DO ji = 1, nbdelay 609 IF( TRIM(cdelay) == TRIM(c_delaylist(ji)) ) idvar = ji 610 END DO 611 IF ( idvar == -1 ) CALL ctl_stop( 'STOP',' mpp_delay_sum : please add a new delayed exchange for '//TRIM(cdname) ) 612 613 IF ( ndelayid(idvar) == 0 ) THEN ! first call with restart: %z1d defined in iom_delay_rst 614 ! -------------------------- 615 IF ( SIZE(todelay(idvar)%z1d) /= isz ) THEN ! Check dimension coherence 616 IF(lwp) WRITE(numout,*) ' WARNING: the nb of delayed variables in restart file is not the model one' 617 DEALLOCATE(todelay(idvar)%z1d) 618 ndelayid(idvar) = -1 ! do as if we had no restart 619 ELSE 620 ALLOCATE(todelay(idvar)%y1d(isz)) 621 todelay(idvar)%y1d(:) = CMPLX(todelay(idvar)%z1d(:), 0., wp) ! create %y1d, complex variable needed by mpi_sumdd 622 END IF 623 ENDIF 624 625 IF( ndelayid(idvar) == -1 ) THEN ! first call without restart: define %y1d and %z1d from y_in with blocking allreduce 626 ! -------------------------- 627 ALLOCATE(todelay(idvar)%z1d(isz), todelay(idvar)%y1d(isz)) ! allocate also %z1d as used for the restart 628 CALL mpi_allreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ierr ) ! get %y1d 629 todelay(idvar)%z1d(:) = REAL(todelay(idvar)%y1d(:), wp) ! define %z1d from %y1d 630 ENDIF 631 632 IF( ndelayid(idvar) > 0 ) CALL mpp_delay_rcv( idvar ) ! make sure %z1d is received 633 634 ! send back pout from todelay(idvar)%z1d defined at previous call 635 pout(:) = todelay(idvar)%z1d(:) 636 637 ! send y_in into todelay(idvar)%y1d with a non-blocking communication 638 CALL mpi_iallreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ndelayid(idvar), ierr ) 639 640 END SUBROUTINE mpp_delay_sum 641 642 643 SUBROUTINE mpp_delay_max( cdname, cdelay, p_in, pout, ldlast, kcom ) 644 !!---------------------------------------------------------------------- 645 !! *** routine mpp_delay_max *** 646 !! 647 !! ** Purpose : performed delayed mpp_max, the result is received on next call 648 !! 649 !!---------------------------------------------------------------------- 650 CHARACTER(len=*), INTENT(in ) :: cdname ! name of the calling subroutine 651 CHARACTER(len=*), INTENT(in ) :: cdelay ! name (used as id) of the delayed operation 652 REAL(wp), INTENT(in ), DIMENSION(:) :: p_in ! 653 REAL(wp), INTENT( out), DIMENSION(:) :: pout ! 654 LOGICAL, INTENT(in ) :: ldlast ! true if this is the last time we call this routine 655 INTEGER, INTENT(in ), OPTIONAL :: kcom 656 !! 657 INTEGER :: ji, isz 658 INTEGER :: idvar 659 INTEGER :: ierr, ilocalcomm 660 !!---------------------------------------------------------------------- 661 ilocalcomm = mpi_comm_oce 662 IF( PRESENT(kcom) ) ilocalcomm = kcom 663 664 isz = SIZE(p_in) 665 666 IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, ld_glb = .TRUE. ) 667 668 idvar = -1 669 DO ji = 1, nbdelay 670 IF( TRIM(cdelay) == TRIM(c_delaylist(ji)) ) idvar = ji 671 END DO 672 IF ( idvar == -1 ) CALL ctl_stop( 'STOP',' mpp_delay_max : please add a new delayed exchange for '//TRIM(cdname) ) 673 674 IF ( ndelayid(idvar) == 0 ) THEN ! first call with restart: %z1d defined in iom_delay_rst 675 ! -------------------------- 676 IF ( SIZE(todelay(idvar)%z1d) /= isz ) THEN ! Check dimension coherence 677 IF(lwp) WRITE(numout,*) ' WARNING: the nb of delayed variables in restart file is not the model one' 678 DEALLOCATE(todelay(idvar)%z1d) 679 ndelayid(idvar) = -1 ! do as if we had no restart 680 END IF 681 ENDIF 682 683 IF( ndelayid(idvar) == -1 ) THEN ! first call without restart: define %z1d from p_in with a blocking allreduce 684 ! -------------------------- 685 ALLOCATE(todelay(idvar)%z1d(isz)) 686 CALL mpi_allreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ierr ) ! get %z1d 687 ENDIF 688 689 IF( ndelayid(idvar) > 0 ) CALL mpp_delay_rcv( idvar ) ! make sure %z1d is received 690 691 ! send back pout from todelay(idvar)%z1d defined at previous call 692 pout(:) = todelay(idvar)%z1d(:) 693 694 ! send p_in into todelay(idvar)%z1d with a non-blocking communication 695 CALL mpi_iallreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ndelayid(idvar), ierr ) 696 697 END SUBROUTINE mpp_delay_max 698 699 700 SUBROUTINE mpp_delay_rcv( kid ) 701 !!---------------------------------------------------------------------- 702 !! *** routine mpp_delay_rcv *** 703 !! 704 !! ** Purpose : force barrier for delayed mpp (needed for restart) 705 !! 706 !!---------------------------------------------------------------------- 707 INTEGER,INTENT(in ) :: kid 708 INTEGER :: ierr 709 !!---------------------------------------------------------------------- 710 IF( ndelayid(kid) /= -2 ) THEN 711 IF( ln_timing ) CALL tic_tac( .TRUE., ld_global = .TRUE.) 712 CALL mpi_wait( ndelayid(kid), MPI_STATUS_IGNORE, ierr ) ! make sure todelay(kid) is received 713 IF( ln_timing ) CALL tic_tac(.FALSE., ld_global = .TRUE.) 714 IF( ASSOCIATED(todelay(kid)%y1d) ) todelay(kid)%z1d(:) = REAL(todelay(kid)%y1d(:), wp) ! define %z1d from %y1d 715 ndelayid(kid) = -2 ! add flag to know that mpi_wait was already called on kid 716 ENDIF 717 END SUBROUTINE mpp_delay_rcv 718 719 576 720 !!---------------------------------------------------------------------- 577 721 !! *** mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real *** … … 579 723 !!---------------------------------------------------------------------- 580 724 !! 581 SUBROUTINE mppmax_a_int( ktab, kdim, kcom ) 582 !!---------------------------------------------------------------------- 583 INTEGER , INTENT(in ) :: kdim ! size of array 584 INTEGER , INTENT(inout), DIMENSION(kdim) :: ktab ! input array 585 INTEGER , INTENT(in ), OPTIONAL :: kcom ! 586 INTEGER :: ierror, ilocalcomm ! temporary integer 587 INTEGER, DIMENSION(kdim) :: iwork 588 !!---------------------------------------------------------------------- 589 ilocalcomm = mpi_comm_oce 590 IF( PRESENT(kcom) ) ilocalcomm = kcom 591 CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_max, ilocalcomm, ierror ) 592 ktab(:) = iwork(:) 593 END SUBROUTINE mppmax_a_int 594 !! 595 SUBROUTINE mppmax_int( ktab, kcom ) 596 !!---------------------------------------------------------------------- 597 INTEGER, INTENT(inout) :: ktab ! ??? 598 INTEGER, INTENT(in ), OPTIONAL :: kcom ! ??? 599 INTEGER :: ierror, iwork, ilocalcomm ! temporary integer 600 !!---------------------------------------------------------------------- 601 ilocalcomm = mpi_comm_oce 602 IF( PRESENT(kcom) ) ilocalcomm = kcom 603 CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_max, ilocalcomm, ierror ) 604 ktab = iwork 605 END SUBROUTINE mppmax_int 606 !! 607 SUBROUTINE mppmax_a_real( ptab, kdim, kcom ) 608 !!---------------------------------------------------------------------- 609 REAL(wp), DIMENSION(kdim), INTENT(inout) :: ptab 610 INTEGER , INTENT(in ) :: kdim 611 INTEGER , OPTIONAL , INTENT(in ) :: kcom 612 INTEGER :: ierror, ilocalcomm 613 REAL(wp), DIMENSION(kdim) :: zwork 614 !!---------------------------------------------------------------------- 615 ilocalcomm = mpi_comm_oce 616 IF( PRESENT(kcom) ) ilocalcomm = kcom 617 CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_max, ilocalcomm, ierror ) 618 ptab(:) = zwork(:) 619 END SUBROUTINE mppmax_a_real 620 !! 621 SUBROUTINE mppmax_real( ptab, kcom ) 622 !!---------------------------------------------------------------------- 623 REAL(wp), INTENT(inout) :: ptab ! ??? 624 INTEGER , INTENT(in ), OPTIONAL :: kcom ! ??? 625 INTEGER :: ierror, ilocalcomm 626 REAL(wp) :: zwork 627 !!---------------------------------------------------------------------- 628 ilocalcomm = mpi_comm_oce 629 IF( PRESENT(kcom) ) ilocalcomm = kcom! 630 CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_max, ilocalcomm, ierror ) 631 ptab = zwork 632 END SUBROUTINE mppmax_real 633 634 725 # define OPERATION_MAX 726 # define INTEGER_TYPE 727 # define DIM_0d 728 # define ROUTINE_ALLREDUCE mppmax_int 729 # include "mpp_allreduce_generic.h90" 730 # undef ROUTINE_ALLREDUCE 731 # undef DIM_0d 732 # define DIM_1d 733 # define ROUTINE_ALLREDUCE mppmax_a_int 734 # include "mpp_allreduce_generic.h90" 735 # undef ROUTINE_ALLREDUCE 736 # undef DIM_1d 737 # undef INTEGER_TYPE 738 ! 739 # define REAL_TYPE 740 # define DIM_0d 741 # define ROUTINE_ALLREDUCE mppmax_real 742 # include "mpp_allreduce_generic.h90" 743 # undef ROUTINE_ALLREDUCE 744 # undef DIM_0d 745 # define DIM_1d 746 # define ROUTINE_ALLREDUCE mppmax_a_real 747 # include "mpp_allreduce_generic.h90" 748 # undef ROUTINE_ALLREDUCE 749 # undef DIM_1d 750 # undef REAL_TYPE 751 # undef OPERATION_MAX 635 752 !!---------------------------------------------------------------------- 636 753 !! *** mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real *** … … 638 755 !!---------------------------------------------------------------------- 639 756 !! 640 SUBROUTINE mppmin_a_int( ktab, kdim, kcom ) 641 !!---------------------------------------------------------------------- 642 INTEGER , INTENT( in ) :: kdim ! size of array 643 INTEGER , INTENT(inout), DIMENSION(kdim) :: ktab ! input array 644 INTEGER , INTENT( in ), OPTIONAL :: kcom ! input array 645 !! 646 INTEGER :: ierror, ilocalcomm ! temporary integer 647 INTEGER, DIMENSION(kdim) :: iwork 648 !!---------------------------------------------------------------------- 649 ilocalcomm = mpi_comm_oce 650 IF( PRESENT(kcom) ) ilocalcomm = kcom 651 CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_min, ilocalcomm, ierror ) 652 ktab(:) = iwork(:) 653 END SUBROUTINE mppmin_a_int 654 !! 655 SUBROUTINE mppmin_int( ktab, kcom ) 656 !!---------------------------------------------------------------------- 657 INTEGER, INTENT(inout) :: ktab ! ??? 658 INTEGER , INTENT( in ), OPTIONAL :: kcom ! input array 659 !! 660 INTEGER :: ierror, iwork, ilocalcomm 661 !!---------------------------------------------------------------------- 662 ilocalcomm = mpi_comm_oce 663 IF( PRESENT(kcom) ) ilocalcomm = kcom 664 CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_min, ilocalcomm, ierror ) 665 ktab = iwork 666 END SUBROUTINE mppmin_int 667 !! 668 SUBROUTINE mppmin_a_real( ptab, kdim, kcom ) 669 !!---------------------------------------------------------------------- 670 INTEGER , INTENT(in ) :: kdim 671 REAL(wp), INTENT(inout), DIMENSION(kdim) :: ptab 672 INTEGER , INTENT(in ), OPTIONAL :: kcom 673 INTEGER :: ierror, ilocalcomm 674 REAL(wp), DIMENSION(kdim) :: zwork 675 !!----------------------------------------------------------------------- 676 ilocalcomm = mpi_comm_oce 677 IF( PRESENT(kcom) ) ilocalcomm = kcom 678 CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_min, ilocalcomm, ierror ) 679 ptab(:) = zwork(:) 680 END SUBROUTINE mppmin_a_real 681 !! 682 SUBROUTINE mppmin_real( ptab, kcom ) 683 !!----------------------------------------------------------------------- 684 REAL(wp), INTENT(inout) :: ptab ! 685 INTEGER , INTENT(in ), OPTIONAL :: kcom 686 INTEGER :: ierror, ilocalcomm 687 REAL(wp) :: zwork 688 !!----------------------------------------------------------------------- 689 ilocalcomm = mpi_comm_oce 690 IF( PRESENT(kcom) ) ilocalcomm = kcom 691 CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_min, ilocalcomm, ierror ) 692 ptab = zwork 693 END SUBROUTINE mppmin_real 694 757 # define OPERATION_MIN 758 # define INTEGER_TYPE 759 # define DIM_0d 760 # define ROUTINE_ALLREDUCE mppmin_int 761 # include "mpp_allreduce_generic.h90" 762 # undef ROUTINE_ALLREDUCE 763 # undef DIM_0d 764 # define DIM_1d 765 # define ROUTINE_ALLREDUCE mppmin_a_int 766 # include "mpp_allreduce_generic.h90" 767 # undef ROUTINE_ALLREDUCE 768 # undef DIM_1d 769 # undef INTEGER_TYPE 770 ! 771 # define REAL_TYPE 772 # define DIM_0d 773 # define ROUTINE_ALLREDUCE mppmin_real 774 # include "mpp_allreduce_generic.h90" 775 # undef ROUTINE_ALLREDUCE 776 # undef DIM_0d 777 # define DIM_1d 778 # define ROUTINE_ALLREDUCE mppmin_a_real 779 # include "mpp_allreduce_generic.h90" 780 # undef ROUTINE_ALLREDUCE 781 # undef DIM_1d 782 # undef REAL_TYPE 783 # undef OPERATION_MIN 695 784 696 785 !!---------------------------------------------------------------------- … … 700 789 !!---------------------------------------------------------------------- 701 790 !! 702 SUBROUTINE mppsum_a_int( ktab, kdim ) 703 !!---------------------------------------------------------------------- 704 INTEGER, INTENT(in ) :: kdim ! ??? 705 INTEGER, INTENT(inout), DIMENSION (kdim) :: ktab ! ??? 706 INTEGER :: ierror 707 INTEGER, DIMENSION (kdim) :: iwork 708 !!---------------------------------------------------------------------- 709 CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_sum, mpi_comm_oce, ierror ) 710 ktab(:) = iwork(:) 711 END SUBROUTINE mppsum_a_int 791 # define OPERATION_SUM 792 # define INTEGER_TYPE 793 # define DIM_0d 794 # define ROUTINE_ALLREDUCE mppsum_int 795 # include "mpp_allreduce_generic.h90" 796 # undef ROUTINE_ALLREDUCE 797 # undef DIM_0d 798 # define DIM_1d 799 # define ROUTINE_ALLREDUCE mppsum_a_int 800 # include "mpp_allreduce_generic.h90" 801 # undef ROUTINE_ALLREDUCE 802 # undef DIM_1d 803 # undef INTEGER_TYPE 804 ! 805 # define REAL_TYPE 806 # define DIM_0d 807 # define ROUTINE_ALLREDUCE mppsum_real 808 # include "mpp_allreduce_generic.h90" 809 # undef ROUTINE_ALLREDUCE 810 # undef DIM_0d 811 # define DIM_1d 812 # define ROUTINE_ALLREDUCE mppsum_a_real 813 # include "mpp_allreduce_generic.h90" 814 # undef ROUTINE_ALLREDUCE 815 # undef DIM_1d 816 # undef REAL_TYPE 817 # undef OPERATION_SUM 818 819 # define OPERATION_SUM_DD 820 # define COMPLEX_TYPE 821 # define DIM_0d 822 # define ROUTINE_ALLREDUCE mppsum_realdd 823 # include "mpp_allreduce_generic.h90" 824 # undef ROUTINE_ALLREDUCE 825 # undef DIM_0d 826 # define DIM_1d 827 # define ROUTINE_ALLREDUCE mppsum_a_realdd 828 # include "mpp_allreduce_generic.h90" 829 # undef ROUTINE_ALLREDUCE 830 # undef DIM_1d 831 # undef COMPLEX_TYPE 832 # undef OPERATION_SUM_DD 833 834 !!---------------------------------------------------------------------- 835 !! *** mpp_minloc2d, mpp_minloc3d, mpp_maxloc2d, mpp_maxloc3d 836 !! 837 !!---------------------------------------------------------------------- 712 838 !! 713 SUBROUTINE mppsum_int( ktab ) 714 !!---------------------------------------------------------------------- 715 INTEGER, INTENT(inout) :: ktab 716 INTEGER :: ierror, iwork 717 !!---------------------------------------------------------------------- 718 CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_sum, mpi_comm_oce, ierror ) 719 ktab = iwork 720 END SUBROUTINE mppsum_int 721 !! 722 SUBROUTINE mppsum_a_real( ptab, kdim, kcom ) 723 !!----------------------------------------------------------------------- 724 INTEGER , INTENT(in ) :: kdim ! size of ptab 725 REAL(wp), DIMENSION(kdim), INTENT(inout) :: ptab ! input array 726 INTEGER , OPTIONAL , INTENT(in ) :: kcom ! specific communicator 727 INTEGER :: ierror, ilocalcomm ! local integer 728 REAL(wp) :: zwork(kdim) ! local workspace 729 !!----------------------------------------------------------------------- 730 ilocalcomm = mpi_comm_oce 731 IF( PRESENT(kcom) ) ilocalcomm = kcom 732 CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_sum, ilocalcomm, ierror ) 733 ptab(:) = zwork(:) 734 END SUBROUTINE mppsum_a_real 735 !! 736 SUBROUTINE mppsum_real( ptab, kcom ) 737 !!----------------------------------------------------------------------- 738 REAL(wp) , INTENT(inout) :: ptab ! input scalar 739 INTEGER , OPTIONAL, INTENT(in ) :: kcom 740 INTEGER :: ierror, ilocalcomm 741 REAL(wp) :: zwork 742 !!----------------------------------------------------------------------- 743 ilocalcomm = mpi_comm_oce 744 IF( PRESENT(kcom) ) ilocalcomm = kcom 745 CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_sum, ilocalcomm, ierror ) 746 ptab = zwork 747 END SUBROUTINE mppsum_real 748 !! 749 SUBROUTINE mppsum_realdd( ytab, kcom ) 750 !!----------------------------------------------------------------------- 751 COMPLEX(wp) , INTENT(inout) :: ytab ! input scalar 752 INTEGER , OPTIONAL, INTENT(in ) :: kcom 753 INTEGER :: ierror, ilocalcomm 754 COMPLEX(wp) :: zwork 755 !!----------------------------------------------------------------------- 756 ilocalcomm = mpi_comm_oce 757 IF( PRESENT(kcom) ) ilocalcomm = kcom 758 CALL MPI_ALLREDUCE( ytab, zwork, 1, MPI_DOUBLE_COMPLEX, MPI_SUMDD, ilocalcomm, ierror ) 759 ytab = zwork 760 END SUBROUTINE mppsum_realdd 761 !! 762 SUBROUTINE mppsum_a_realdd( ytab, kdim, kcom ) 763 !!---------------------------------------------------------------------- 764 INTEGER , INTENT(in ) :: kdim ! size of ytab 765 COMPLEX(wp), DIMENSION(kdim), INTENT(inout) :: ytab ! input array 766 INTEGER , OPTIONAL , INTENT(in ) :: kcom 767 INTEGER:: ierror, ilocalcomm ! local integer 768 COMPLEX(wp), DIMENSION(kdim) :: zwork ! temporary workspace 769 !!----------------------------------------------------------------------- 770 ilocalcomm = mpi_comm_oce 771 IF( PRESENT(kcom) ) ilocalcomm = kcom 772 CALL MPI_ALLREDUCE( ytab, zwork, kdim, MPI_DOUBLE_COMPLEX, MPI_SUMDD, ilocalcomm, ierror ) 773 ytab(:) = zwork(:) 774 END SUBROUTINE mppsum_a_realdd 775 776 777 SUBROUTINE mppmax_real_multiple( pt1d, kdim, kcom ) 778 !!---------------------------------------------------------------------- 779 !! *** routine mppmax_real *** 780 !! 781 !! ** Purpose : Maximum across processor of each element of a 1D arrays 782 !! 783 !!---------------------------------------------------------------------- 784 REAL(wp), DIMENSION(kdim), INTENT(inout) :: pt1d ! 1D arrays 785 INTEGER , INTENT(in ) :: kdim 786 INTEGER , OPTIONAL , INTENT(in ) :: kcom ! local communicator 787 !! 788 INTEGER :: ierror, ilocalcomm 789 REAL(wp), DIMENSION(kdim) :: zwork 790 !!---------------------------------------------------------------------- 791 ilocalcomm = mpi_comm_oce 792 IF( PRESENT(kcom) ) ilocalcomm = kcom 793 ! 794 CALL mpi_allreduce( pt1d, zwork, kdim, mpi_double_precision, mpi_max, ilocalcomm, ierror ) 795 pt1d(:) = zwork(:) 796 ! 797 END SUBROUTINE mppmax_real_multiple 798 799 800 SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki,kj ) 801 !!------------------------------------------------------------------------ 802 !! *** routine mpp_minloc *** 803 !! 804 !! ** Purpose : Compute the global minimum of an array ptab 805 !! and also give its global position 806 !! 807 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 808 !! 809 !!-------------------------------------------------------------------------- 810 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: ptab ! Local 2D array 811 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: pmask ! Local mask 812 REAL(wp) , INTENT( out) :: pmin ! Global minimum of ptab 813 INTEGER , INTENT( out) :: ki, kj ! index of minimum in global frame 814 ! 815 INTEGER :: ierror 816 INTEGER , DIMENSION(2) :: ilocs 817 REAL(wp) :: zmin ! local minimum 818 REAL(wp), DIMENSION(2,1) :: zain, zaout 819 !!----------------------------------------------------------------------- 820 ! 821 zmin = MINVAL( ptab(:,:) , mask= pmask == 1._wp ) 822 ilocs = MINLOC( ptab(:,:) , mask= pmask == 1._wp ) 823 ! 824 ki = ilocs(1) + nimpp - 1 825 kj = ilocs(2) + njmpp - 1 826 ! 827 zain(1,:)=zmin 828 zain(2,:)=ki+10000.*kj 829 ! 830 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OCE,ierror) 831 ! 832 pmin = zaout(1,1) 833 kj = INT(zaout(2,1)/10000.) 834 ki = INT(zaout(2,1) - 10000.*kj ) 835 ! 836 END SUBROUTINE mpp_minloc2d 837 838 839 SUBROUTINE mpp_minloc3d( ptab, pmask, pmin, ki, kj ,kk) 840 !!------------------------------------------------------------------------ 841 !! *** routine mpp_minloc *** 842 !! 843 !! ** Purpose : Compute the global minimum of an array ptab 844 !! and also give its global position 845 !! 846 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 847 !! 848 !!-------------------------------------------------------------------------- 849 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: ptab ! Local 2D array 850 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pmask ! Local mask 851 REAL(wp) , INTENT( out) :: pmin ! Global minimum of ptab 852 INTEGER , INTENT( out) :: ki, kj, kk ! index of minimum in global frame 853 ! 854 INTEGER :: ierror 855 REAL(wp) :: zmin ! local minimum 856 INTEGER , DIMENSION(3) :: ilocs 857 REAL(wp), DIMENSION(2,1) :: zain, zaout 858 !!----------------------------------------------------------------------- 859 ! 860 zmin = MINVAL( ptab(:,:,:) , mask= pmask == 1._wp ) 861 ilocs = MINLOC( ptab(:,:,:) , mask= pmask == 1._wp ) 862 ! 863 ki = ilocs(1) + nimpp - 1 864 kj = ilocs(2) + njmpp - 1 865 kk = ilocs(3) 866 ! 867 zain(1,:) = zmin 868 zain(2,:) = ki + 10000.*kj + 100000000.*kk 869 ! 870 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OCE,ierror) 871 ! 872 pmin = zaout(1,1) 873 kk = INT( zaout(2,1) / 100000000. ) 874 kj = INT( zaout(2,1) - kk * 100000000. ) / 10000 875 ki = INT( zaout(2,1) - kk * 100000000. -kj * 10000. ) 876 ! 877 END SUBROUTINE mpp_minloc3d 878 879 880 SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj ) 881 !!------------------------------------------------------------------------ 882 !! *** routine mpp_maxloc *** 883 !! 884 !! ** Purpose : Compute the global maximum of an array ptab 885 !! and also give its global position 886 !! 887 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 888 !! 889 !!-------------------------------------------------------------------------- 890 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: ptab ! Local 2D array 891 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: pmask ! Local mask 892 REAL(wp) , INTENT( out) :: pmax ! Global maximum of ptab 893 INTEGER , INTENT( out) :: ki, kj ! index of maximum in global frame 894 !! 895 INTEGER :: ierror 896 INTEGER, DIMENSION (2) :: ilocs 897 REAL(wp) :: zmax ! local maximum 898 REAL(wp), DIMENSION(2,1) :: zain, zaout 899 !!----------------------------------------------------------------------- 900 ! 901 zmax = MAXVAL( ptab(:,:) , mask= pmask == 1._wp ) 902 ilocs = MAXLOC( ptab(:,:) , mask= pmask == 1._wp ) 903 ! 904 ki = ilocs(1) + nimpp - 1 905 kj = ilocs(2) + njmpp - 1 906 ! 907 zain(1,:) = zmax 908 zain(2,:) = ki + 10000. * kj 909 ! 910 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OCE,ierror) 911 ! 912 pmax = zaout(1,1) 913 kj = INT( zaout(2,1) / 10000. ) 914 ki = INT( zaout(2,1) - 10000.* kj ) 915 ! 916 END SUBROUTINE mpp_maxloc2d 917 918 919 SUBROUTINE mpp_maxloc3d( ptab, pmask, pmax, ki, kj, kk ) 920 !!------------------------------------------------------------------------ 921 !! *** routine mpp_maxloc *** 922 !! 923 !! ** Purpose : Compute the global maximum of an array ptab 924 !! and also give its global position 925 !! 926 !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC 927 !! 928 !!-------------------------------------------------------------------------- 929 REAL(wp), DIMENSION (:,:,:), INTENT(in ) :: ptab ! Local 2D array 930 REAL(wp), DIMENSION (:,:,:), INTENT(in ) :: pmask ! Local mask 931 REAL(wp) , INTENT( out) :: pmax ! Global maximum of ptab 932 INTEGER , INTENT( out) :: ki, kj, kk ! index of maximum in global frame 933 ! 934 INTEGER :: ierror ! local integer 935 REAL(wp) :: zmax ! local maximum 936 REAL(wp), DIMENSION(2,1) :: zain, zaout 937 INTEGER , DIMENSION(3) :: ilocs 938 !!----------------------------------------------------------------------- 939 ! 940 zmax = MAXVAL( ptab(:,:,:) , mask= pmask == 1._wp ) 941 ilocs = MAXLOC( ptab(:,:,:) , mask= pmask == 1._wp ) 942 ! 943 ki = ilocs(1) + nimpp - 1 944 kj = ilocs(2) + njmpp - 1 945 kk = ilocs(3) 946 ! 947 zain(1,:) = zmax 948 zain(2,:) = ki + 10000.*kj + 100000000.*kk 949 ! 950 CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OCE,ierror) 951 ! 952 pmax = zaout(1,1) 953 kk = INT( zaout(2,1) / 100000000. ) 954 kj = INT( zaout(2,1) - kk * 100000000. ) / 10000 955 ki = INT( zaout(2,1) - kk * 100000000. -kj * 10000. ) 956 ! 957 END SUBROUTINE mpp_maxloc3d 958 839 # define OPERATION_MINLOC 840 # define DIM_2d 841 # define ROUTINE_LOC mpp_minloc2d 842 # include "mpp_loc_generic.h90" 843 # undef ROUTINE_LOC 844 # undef DIM_2d 845 # define DIM_3d 846 # define ROUTINE_LOC mpp_minloc3d 847 # include "mpp_loc_generic.h90" 848 # undef ROUTINE_LOC 849 # undef DIM_3d 850 # undef OPERATION_MINLOC 851 852 # define OPERATION_MAXLOC 853 # define DIM_2d 854 # define ROUTINE_LOC mpp_maxloc2d 855 # include "mpp_loc_generic.h90" 856 # undef ROUTINE_LOC 857 # undef DIM_2d 858 # define DIM_3d 859 # define ROUTINE_LOC mpp_maxloc3d 860 # include "mpp_loc_generic.h90" 861 # undef ROUTINE_LOC 862 # undef DIM_3d 863 # undef OPERATION_MAXLOC 959 864 960 865 SUBROUTINE mppsync() … … 973 878 974 879 975 SUBROUTINE mppstop 880 SUBROUTINE mppstop( ldfinal, ld_force_abort ) 976 881 !!---------------------------------------------------------------------- 977 882 !! *** routine mppstop *** … … 980 885 !! 981 886 !!---------------------------------------------------------------------- 887 LOGICAL, OPTIONAL, INTENT(in) :: ldfinal ! source process number 888 LOGICAL, OPTIONAL, INTENT(in) :: ld_force_abort ! source process number 889 LOGICAL :: llfinal, ll_force_abort 982 890 INTEGER :: info 983 891 !!---------------------------------------------------------------------- 984 ! 985 CALL mppsync 986 CALL mpi_finalize( info ) 892 llfinal = .FALSE. 893 IF( PRESENT(ldfinal) ) llfinal = ldfinal 894 ll_force_abort = .FALSE. 895 IF( PRESENT(ld_force_abort) ) ll_force_abort = ld_force_abort 896 ! 897 IF(ll_force_abort) THEN 898 CALL mpi_abort( MPI_COMM_WORLD ) 899 ELSE 900 CALL mppsync 901 CALL mpi_finalize( info ) 902 ENDIF 903 IF( .NOT. llfinal ) STOP 123456 987 904 ! 988 905 END SUBROUTINE mppstop … … 999 916 ! 1000 917 END SUBROUTINE mpp_comm_free 1001 1002 1003 SUBROUTINE mpp_ini_ice( pindic, kumout )1004 !!----------------------------------------------------------------------1005 !! *** routine mpp_ini_ice ***1006 !!1007 !! ** Purpose : Initialize special communicator for ice areas1008 !! condition together with global variables needed in the ddmpp folding1009 !!1010 !! ** Method : - Look for ice processors in ice routines1011 !! - Put their number in nrank_ice1012 !! - Create groups for the world processors and the ice processors1013 !! - Create a communicator for ice processors1014 !!1015 !! ** output1016 !! njmppmax = njmpp for northern procs1017 !! ndim_rank_ice = number of processors with ice1018 !! nrank_ice (ndim_rank_ice) = ice processors1019 !! ngrp_iworld = group ID for the world processors1020 !! ngrp_ice = group ID for the ice processors1021 !! ncomm_ice = communicator for the ice procs.1022 !! n_ice_root = number (in the world) of proc 0 in the ice comm.1023 !!1024 !!----------------------------------------------------------------------1025 INTEGER, INTENT(in) :: pindic1026 INTEGER, INTENT(in) :: kumout ! ocean.output logical unit1027 !!1028 INTEGER :: jjproc1029 INTEGER :: ii, ierr1030 INTEGER, ALLOCATABLE, DIMENSION(:) :: kice1031 INTEGER, ALLOCATABLE, DIMENSION(:) :: zwork1032 !!----------------------------------------------------------------------1033 !1034 ALLOCATE( kice(jpnij), zwork(jpnij), STAT=ierr )1035 IF( ierr /= 0 ) THEN1036 WRITE(kumout, cform_err)1037 WRITE(kumout,*) 'mpp_ini_ice : failed to allocate 2, 1D arrays (jpnij in length)'1038 CALL mppstop1039 ENDIF1040 1041 ! Look for how many procs with sea-ice1042 !1043 kice = 01044 DO jjproc = 1, jpnij1045 IF( jjproc == narea .AND. pindic .GT. 0 ) kice(jjproc) = 11046 END DO1047 !1048 zwork = 01049 CALL MPI_ALLREDUCE( kice, zwork, jpnij, mpi_integer, mpi_sum, mpi_comm_oce, ierr )1050 ndim_rank_ice = SUM( zwork )1051 1052 ! Allocate the right size to nrank_north1053 IF( ALLOCATED ( nrank_ice ) ) DEALLOCATE( nrank_ice )1054 ALLOCATE( nrank_ice(ndim_rank_ice) )1055 !1056 ii = 01057 nrank_ice = 01058 DO jjproc = 1, jpnij1059 IF( zwork(jjproc) == 1) THEN1060 ii = ii + 11061 nrank_ice(ii) = jjproc -11062 ENDIF1063 END DO1064 1065 ! Create the world group1066 CALL MPI_COMM_GROUP( mpi_comm_oce, ngrp_iworld, ierr )1067 1068 ! Create the ice group from the world group1069 CALL MPI_GROUP_INCL( ngrp_iworld, ndim_rank_ice, nrank_ice, ngrp_ice, ierr )1070 1071 ! Create the ice communicator , ie the pool of procs with sea-ice1072 CALL MPI_COMM_CREATE( mpi_comm_oce, ngrp_ice, ncomm_ice, ierr )1073 1074 ! Find proc number in the world of proc 0 in the north1075 ! The following line seems to be useless, we just comment & keep it as reminder1076 ! CALL MPI_GROUP_TRANSLATE_RANKS(ngrp_ice,1,0,ngrp_iworld,n_ice_root,ierr)1077 !1078 CALL MPI_GROUP_FREE(ngrp_ice, ierr)1079 CALL MPI_GROUP_FREE(ngrp_iworld, ierr)1080 1081 DEALLOCATE(kice, zwork)1082 !1083 END SUBROUTINE mpp_ini_ice1084 918 1085 919 … … 1175 1009 l_znl_root = .FALSE. 1176 1010 kwork (1) = nimpp 1177 CALL mpp_min ( kwork(1), kcom = ncomm_znl)1011 CALL mpp_min ( 'lib_mpp', kwork(1), kcom = ncomm_znl) 1178 1012 IF ( nimpp == kwork(1)) l_znl_root = .TRUE. 1179 1013 END IF … … 1384 1218 ! 1385 1219 itaille = jpimax * ( ipj + 2*kextj ) 1220 ! 1221 IF( ln_timing ) CALL tic_tac(.TRUE.) 1386 1222 CALL MPI_ALLGATHER( znorthloc_e(1,1-kextj) , itaille, MPI_DOUBLE_PRECISION, & 1387 1223 & znorthgloio_e(1,1-kextj,1), itaille, MPI_DOUBLE_PRECISION, & 1388 1224 & ncomm_north, ierr ) 1225 ! 1226 IF( ln_timing ) CALL tic_tac(.FALSE.) 1389 1227 ! 1390 1228 DO jr = 1, ndim_rank_north ! recover the global north array … … 1418 1256 1419 1257 1420 SUBROUTINE mpp_lnk_2d_icb( pt2d, cd_type, psgn, kexti, kextj )1258 SUBROUTINE mpp_lnk_2d_icb( cdname, pt2d, cd_type, psgn, kexti, kextj ) 1421 1259 !!---------------------------------------------------------------------- 1422 1260 !! *** routine mpp_lnk_2d_icb *** … … 1440 1278 !! nono : number for local neighboring processors 1441 1279 !!---------------------------------------------------------------------- 1280 CHARACTER(len=*) , INTENT(in ) :: cdname ! name of the calling subroutine 1442 1281 REAL(wp), DIMENSION(1-kexti:jpi+kexti,1-kextj:jpj+kextj), INTENT(inout) :: pt2d ! 2D array with extra halo 1443 1282 CHARACTER(len=1) , INTENT(in ) :: cd_type ! nature of ptab array grid-points … … 1459 1298 iprecj = nn_hls + kextj 1460 1299 1300 IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, 1, 1, 1, ld_lbc = .TRUE. ) 1461 1301 1462 1302 ! 1. standard boundary treatment … … 1510 1350 ! ! Migrations 1511 1351 imigr = ipreci * ( jpj + 2*kextj ) 1352 ! 1353 IF( ln_timing ) CALL tic_tac(.TRUE.) 1512 1354 ! 1513 1355 SELECT CASE ( nbondi ) … … 1529 1371 END SELECT 1530 1372 ! 1373 IF( ln_timing ) CALL tic_tac(.FALSE.) 1374 ! 1531 1375 ! ! Write Dirichlet lateral conditions 1532 1376 iihom = jpi - nn_hls … … 1563 1407 ! ! Migrations 1564 1408 imigr = iprecj * ( jpi + 2*kexti ) 1409 ! 1410 IF( ln_timing ) CALL tic_tac(.TRUE.) 1565 1411 ! 1566 1412 SELECT CASE ( nbondj ) … … 1582 1428 END SELECT 1583 1429 ! 1430 IF( ln_timing ) CALL tic_tac(.FALSE.) 1431 ! 1584 1432 ! ! Write Dirichlet lateral conditions 1585 1433 ijhom = jpj - nn_hls … … 1602 1450 ! 1603 1451 END SUBROUTINE mpp_lnk_2d_icb 1452 1453 1454 SUBROUTINE mpp_report( cdname, kpk, kpl, kpf, ld_lbc, ld_glb ) 1455 !!---------------------------------------------------------------------- 1456 !! *** routine mpp_report *** 1457 !! 1458 !! ** Purpose : report use of mpp routines per time-setp 1459 !! 1460 !!---------------------------------------------------------------------- 1461 CHARACTER(len=*), INTENT(in ) :: cdname ! name of the calling subroutine 1462 INTEGER , OPTIONAL, INTENT(in ) :: kpk, kpl, kpf 1463 LOGICAL , OPTIONAL, INTENT(in ) :: ld_lbc, ld_glb 1464 !! 1465 LOGICAL :: ll_lbc, ll_glb 1466 INTEGER :: ji, jj, jk, jh, jf ! dummy loop indices 1467 !!---------------------------------------------------------------------- 1468 ! 1469 ll_lbc = .FALSE. 1470 IF( PRESENT(ld_lbc) ) ll_lbc = ld_lbc 1471 ll_glb = .FALSE. 1472 IF( PRESENT(ld_glb) ) ll_glb = ld_glb 1473 ! 1474 ! find the smallest common frequency: default = frequency product, if multiple, choose the larger of the 2 frequency 1475 IF( ncom_dttrc /= 1 ) CALL ctl_stop( 'STOP', 'mpp_report, ncom_dttrc /= 1 not coded...' ) 1476 ncom_freq = ncom_fsbc 1477 ! 1478 IF ( ncom_stp == nit000+ncom_freq ) THEN ! avoid to count extra communications in potential initializations at nit000 1479 IF( ll_lbc ) THEN 1480 IF( .NOT. ALLOCATED(ncomm_sequence) ) ALLOCATE( ncomm_sequence(ncom_rec_max,2) ) 1481 IF( .NOT. ALLOCATED( crname_lbc) ) ALLOCATE( crname_lbc(ncom_rec_max ) ) 1482 n_sequence_lbc = n_sequence_lbc + 1 1483 IF( n_sequence_lbc > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' ) ! deadlock 1484 crname_lbc(n_sequence_lbc) = cdname ! keep the name of the calling routine 1485 ncomm_sequence(n_sequence_lbc,1) = kpk*kpl ! size of 3rd and 4th dimensions 1486 ncomm_sequence(n_sequence_lbc,2) = kpf ! number of arrays to be treated (multi) 1487 ENDIF 1488 IF( ll_glb ) THEN 1489 IF( .NOT. ALLOCATED(crname_glb) ) ALLOCATE( crname_glb(ncom_rec_max) ) 1490 n_sequence_glb = n_sequence_glb + 1 1491 IF( n_sequence_glb > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' ) ! deadlock 1492 crname_glb(n_sequence_glb) = cdname ! keep the name of the calling routine 1493 ENDIF 1494 ELSE IF ( ncom_stp == nit000+2*ncom_freq ) THEN 1495 CALL ctl_opn( numcom, 'communication_report.txt', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea ) 1496 WRITE(numcom,*) ' ' 1497 WRITE(numcom,*) ' ------------------------------------------------------------' 1498 WRITE(numcom,*) ' Communication pattern report (second oce+sbc+top time step):' 1499 WRITE(numcom,*) ' ------------------------------------------------------------' 1500 WRITE(numcom,*) ' ' 1501 WRITE(numcom,'(A,I4)') ' Exchanged halos : ', n_sequence_lbc 1502 jj = 0; jk = 0; jf = 0; jh = 0 1503 DO ji = 1, n_sequence_lbc 1504 IF ( ncomm_sequence(ji,1) .GT. 1 ) jk = jk + 1 1505 IF ( ncomm_sequence(ji,2) .GT. 1 ) jf = jf + 1 1506 IF ( ncomm_sequence(ji,1) .GT. 1 .AND. ncomm_sequence(ji,2) .GT. 1 ) jj = jj + 1 1507 jh = MAX (jh, ncomm_sequence(ji,1)*ncomm_sequence(ji,2)) 1508 END DO 1509 WRITE(numcom,'(A,I3)') ' 3D Exchanged halos : ', jk 1510 WRITE(numcom,'(A,I3)') ' Multi arrays exchanged halos : ', jf 1511 WRITE(numcom,'(A,I3)') ' from which 3D : ', jj 1512 WRITE(numcom,'(A,I10)') ' Array max size : ', jh*jpi*jpj 1513 WRITE(numcom,*) ' ' 1514 WRITE(numcom,*) ' lbc_lnk called' 1515 jj = 1 1516 DO ji = 2, n_sequence_lbc 1517 IF( crname_lbc(ji-1) /= crname_lbc(ji) ) THEN 1518 WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_lbc(ji-1)) 1519 jj = 0 1520 END IF 1521 jj = jj + 1 1522 END DO 1523 WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_lbc(n_sequence_lbc)) 1524 WRITE(numcom,*) ' ' 1525 IF ( n_sequence_glb > 0 ) THEN 1526 WRITE(numcom,'(A,I4)') ' Global communications : ', n_sequence_glb 1527 jj = 1 1528 DO ji = 2, n_sequence_glb 1529 IF( crname_glb(ji-1) /= crname_glb(ji) ) THEN 1530 WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_glb(ji-1)) 1531 jj = 0 1532 END IF 1533 jj = jj + 1 1534 END DO 1535 WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_glb(n_sequence_glb)) 1536 DEALLOCATE(crname_glb) 1537 ELSE 1538 WRITE(numcom,*) ' No MPI global communication ' 1539 ENDIF 1540 WRITE(numcom,*) ' ' 1541 WRITE(numcom,*) ' -----------------------------------------------' 1542 WRITE(numcom,*) ' ' 1543 DEALLOCATE(ncomm_sequence) 1544 DEALLOCATE(crname_lbc) 1545 ENDIF 1546 END SUBROUTINE mpp_report 1547 1548 1549 SUBROUTINE tic_tac (ld_tic, ld_global) 1550 1551 LOGICAL, INTENT(IN) :: ld_tic 1552 LOGICAL, OPTIONAL, INTENT(IN) :: ld_global 1553 REAL(wp), DIMENSION(2), SAVE :: tic_wt 1554 REAL(wp), SAVE :: tic_ct = 0._wp 1555 INTEGER :: ii 1556 1557 IF( ncom_stp <= nit000 ) RETURN 1558 IF( ncom_stp == nitend ) RETURN 1559 ii = 1 1560 IF( PRESENT( ld_global ) ) THEN 1561 IF( ld_global ) ii = 2 1562 END IF 1563 1564 IF ( ld_tic ) THEN 1565 tic_wt(ii) = MPI_Wtime() ! start count tic->tac (waiting time) 1566 IF ( tic_ct > 0.0_wp ) compute_time = compute_time + MPI_Wtime() - tic_ct ! cumulate count tac->tic 1567 ELSE 1568 waiting_time(ii) = waiting_time(ii) + MPI_Wtime() - tic_wt(ii) ! cumulate count tic->tac 1569 tic_ct = MPI_Wtime() ! start count tac->tic (waiting time) 1570 ENDIF 1571 1572 END SUBROUTINE tic_tac 1573 1604 1574 1605 1575 #else … … 1610 1580 1611 1581 INTERFACE mpp_sum 1612 MODULE PROCEDURE mpp _sum_a2s, mpp_sum_as, mpp_sum_ai, mpp_sum_s, mpp_sum_i, mppsum_realdd, mppsum_a_realdd1582 MODULE PROCEDURE mppsum_int, mppsum_a_int, mppsum_real, mppsum_a_real, mppsum_realdd, mppsum_a_realdd 1613 1583 END INTERFACE 1614 1584 INTERFACE mpp_max … … 1624 1594 MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d 1625 1595 END INTERFACE 1626 INTERFACE mpp_max_multiple1627 MODULE PROCEDURE mppmax_real_multiple1628 END INTERFACE1629 1596 1630 1597 LOGICAL, PUBLIC, PARAMETER :: lk_mpp = .FALSE. !: mpp flag 1631 1598 LOGICAL, PUBLIC :: ln_nnogather !: namelist control of northfold comms (needed here in case "key_mpp_mpi" is not used) 1632 INTEGER :: ncomm_ice1633 1599 INTEGER, PUBLIC :: mpi_comm_oce ! opa local communicator 1600 1601 INTEGER, PARAMETER, PUBLIC :: nbdelay = 0 ! make sure we don't enter loops: DO ji = 1, nbdelay 1602 CHARACTER(len=32), DIMENSION(1), PUBLIC :: c_delaylist = 'empty' 1603 CHARACTER(len=32), DIMENSION(1), PUBLIC :: c_delaycpnt = 'empty' 1604 TYPE :: DELAYARR 1605 REAL( wp), POINTER, DIMENSION(:) :: z1d => NULL() 1606 COMPLEX(wp), POINTER, DIMENSION(:) :: y1d => NULL() 1607 END TYPE DELAYARR 1608 TYPE( DELAYARR ), DIMENSION(1), PUBLIC :: todelay 1609 INTEGER, PUBLIC, DIMENSION(1) :: ndelayid = -1 1634 1610 !!---------------------------------------------------------------------- 1635 1611 CONTAINS … … 1654 1630 END SUBROUTINE mppsync 1655 1631 1656 SUBROUTINE mpp_sum_as( parr, kdim, kcom ) ! Dummy routine 1657 REAL , DIMENSION(:) :: parr 1658 INTEGER :: kdim 1659 INTEGER, OPTIONAL :: kcom 1660 WRITE(*,*) 'mpp_sum_as: You should not have seen this print! error?', kdim, parr(1), kcom 1661 END SUBROUTINE mpp_sum_as 1662 1663 SUBROUTINE mpp_sum_a2s( parr, kdim, kcom ) ! Dummy routine 1664 REAL , DIMENSION(:,:) :: parr 1665 INTEGER :: kdim 1666 INTEGER, OPTIONAL :: kcom 1667 WRITE(*,*) 'mpp_sum_a2s: You should not have seen this print! error?', kdim, parr(1,1), kcom 1668 END SUBROUTINE mpp_sum_a2s 1669 1670 SUBROUTINE mpp_sum_ai( karr, kdim, kcom ) ! Dummy routine 1671 INTEGER, DIMENSION(:) :: karr 1672 INTEGER :: kdim 1673 INTEGER, OPTIONAL :: kcom 1674 WRITE(*,*) 'mpp_sum_ai: You should not have seen this print! error?', kdim, karr(1), kcom 1675 END SUBROUTINE mpp_sum_ai 1676 1677 SUBROUTINE mpp_sum_s( psca, kcom ) ! Dummy routine 1678 REAL :: psca 1679 INTEGER, OPTIONAL :: kcom 1680 WRITE(*,*) 'mpp_sum_s: You should not have seen this print! error?', psca, kcom 1681 END SUBROUTINE mpp_sum_s 1682 1683 SUBROUTINE mpp_sum_i( kint, kcom ) ! Dummy routine 1684 integer :: kint 1685 INTEGER, OPTIONAL :: kcom 1686 WRITE(*,*) 'mpp_sum_i: You should not have seen this print! error?', kint, kcom 1687 END SUBROUTINE mpp_sum_i 1688 1689 SUBROUTINE mppsum_realdd( ytab, kcom ) 1690 COMPLEX(wp), INTENT(inout) :: ytab ! input scalar 1691 INTEGER , INTENT( in ), OPTIONAL :: kcom 1692 WRITE(*,*) 'mppsum_realdd: You should not have seen this print! error?', ytab 1693 END SUBROUTINE mppsum_realdd 1694 1695 SUBROUTINE mppsum_a_realdd( ytab, kdim, kcom ) 1696 INTEGER , INTENT( in ) :: kdim ! size of ytab 1697 COMPLEX(wp), DIMENSION(kdim), INTENT( inout ) :: ytab ! input array 1698 INTEGER , INTENT( in ), OPTIONAL :: kcom 1699 WRITE(*,*) 'mppsum_a_realdd: You should not have seen this print! error?', kdim, ytab(1), kcom 1700 END SUBROUTINE mppsum_a_realdd 1701 1702 SUBROUTINE mppmax_a_real( parr, kdim, kcom ) 1703 REAL , DIMENSION(:) :: parr 1704 INTEGER :: kdim 1705 INTEGER, OPTIONAL :: kcom 1706 WRITE(*,*) 'mppmax_a_real: You should not have seen this print! error?', kdim, parr(1), kcom 1707 END SUBROUTINE mppmax_a_real 1708 1709 SUBROUTINE mppmax_real( psca, kcom ) 1710 REAL :: psca 1711 INTEGER, OPTIONAL :: kcom 1712 WRITE(*,*) 'mppmax_real: You should not have seen this print! error?', psca, kcom 1713 END SUBROUTINE mppmax_real 1714 1715 SUBROUTINE mppmin_a_real( parr, kdim, kcom ) 1716 REAL , DIMENSION(:) :: parr 1717 INTEGER :: kdim 1718 INTEGER, OPTIONAL :: kcom 1719 WRITE(*,*) 'mppmin_a_real: You should not have seen this print! error?', kdim, parr(1), kcom 1720 END SUBROUTINE mppmin_a_real 1721 1722 SUBROUTINE mppmin_real( psca, kcom ) 1723 REAL :: psca 1724 INTEGER, OPTIONAL :: kcom 1725 WRITE(*,*) 'mppmin_real: You should not have seen this print! error?', psca, kcom 1726 END SUBROUTINE mppmin_real 1727 1728 SUBROUTINE mppmax_a_int( karr, kdim ,kcom) 1729 INTEGER, DIMENSION(:) :: karr 1730 INTEGER :: kdim 1731 INTEGER, OPTIONAL :: kcom 1732 WRITE(*,*) 'mppmax_a_int: You should not have seen this print! error?', kdim, karr(1), kcom 1733 END SUBROUTINE mppmax_a_int 1734 1735 SUBROUTINE mppmax_int( kint, kcom) 1736 INTEGER :: kint 1737 INTEGER, OPTIONAL :: kcom 1738 WRITE(*,*) 'mppmax_int: You should not have seen this print! error?', kint, kcom 1739 END SUBROUTINE mppmax_int 1740 1741 SUBROUTINE mppmin_a_int( karr, kdim, kcom ) 1742 INTEGER, DIMENSION(:) :: karr 1743 INTEGER :: kdim 1744 INTEGER, OPTIONAL :: kcom 1745 WRITE(*,*) 'mppmin_a_int: You should not have seen this print! error?', kdim, karr(1), kcom 1746 END SUBROUTINE mppmin_a_int 1747 1748 SUBROUTINE mppmin_int( kint, kcom ) 1749 INTEGER :: kint 1750 INTEGER, OPTIONAL :: kcom 1751 WRITE(*,*) 'mppmin_int: You should not have seen this print! error?', kint, kcom 1752 END SUBROUTINE mppmin_int 1753 1754 SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki, kj ) 1755 REAL :: pmin 1756 REAL , DIMENSION (:,:) :: ptab, pmask 1757 INTEGER :: ki, kj 1758 WRITE(*,*) 'mpp_minloc2d: You should not have seen this print! error?', pmin, ki, kj, ptab(1,1), pmask(1,1) 1759 END SUBROUTINE mpp_minloc2d 1760 1761 SUBROUTINE mpp_minloc3d( ptab, pmask, pmin, ki, kj, kk ) 1762 REAL :: pmin 1763 REAL , DIMENSION (:,:,:) :: ptab, pmask 1764 INTEGER :: ki, kj, kk 1765 WRITE(*,*) 'mpp_minloc3d: You should not have seen this print! error?', pmin, ki, kj, kk, ptab(1,1,1), pmask(1,1,1) 1766 END SUBROUTINE mpp_minloc3d 1767 1768 SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj ) 1769 REAL :: pmax 1770 REAL , DIMENSION (:,:) :: ptab, pmask 1771 INTEGER :: ki, kj 1772 WRITE(*,*) 'mpp_maxloc2d: You should not have seen this print! error?', pmax, ki, kj, ptab(1,1), pmask(1,1) 1773 END SUBROUTINE mpp_maxloc2d 1774 1775 SUBROUTINE mpp_maxloc3d( ptab, pmask, pmax, ki, kj, kk ) 1776 REAL :: pmax 1777 REAL , DIMENSION (:,:,:) :: ptab, pmask 1778 INTEGER :: ki, kj, kk 1779 WRITE(*,*) 'mpp_maxloc3d: You should not have seen this print! error?', pmax, ki, kj, kk, ptab(1,1,1), pmask(1,1,1) 1780 END SUBROUTINE mpp_maxloc3d 1781 1782 SUBROUTINE mppstop 1632 !!---------------------------------------------------------------------- 1633 !! *** mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real *** 1634 !! 1635 !!---------------------------------------------------------------------- 1636 !! 1637 # define OPERATION_MAX 1638 # define INTEGER_TYPE 1639 # define DIM_0d 1640 # define ROUTINE_ALLREDUCE mppmax_int 1641 # include "mpp_allreduce_generic.h90" 1642 # undef ROUTINE_ALLREDUCE 1643 # undef DIM_0d 1644 # define DIM_1d 1645 # define ROUTINE_ALLREDUCE mppmax_a_int 1646 # include "mpp_allreduce_generic.h90" 1647 # undef ROUTINE_ALLREDUCE 1648 # undef DIM_1d 1649 # undef INTEGER_TYPE 1650 ! 1651 # define REAL_TYPE 1652 # define DIM_0d 1653 # define ROUTINE_ALLREDUCE mppmax_real 1654 # include "mpp_allreduce_generic.h90" 1655 # undef ROUTINE_ALLREDUCE 1656 # undef DIM_0d 1657 # define DIM_1d 1658 # define ROUTINE_ALLREDUCE mppmax_a_real 1659 # include "mpp_allreduce_generic.h90" 1660 # undef ROUTINE_ALLREDUCE 1661 # undef DIM_1d 1662 # undef REAL_TYPE 1663 # undef OPERATION_MAX 1664 !!---------------------------------------------------------------------- 1665 !! *** mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real *** 1666 !! 1667 !!---------------------------------------------------------------------- 1668 !! 1669 # define OPERATION_MIN 1670 # define INTEGER_TYPE 1671 # define DIM_0d 1672 # define ROUTINE_ALLREDUCE mppmin_int 1673 # include "mpp_allreduce_generic.h90" 1674 # undef ROUTINE_ALLREDUCE 1675 # undef DIM_0d 1676 # define DIM_1d 1677 # define ROUTINE_ALLREDUCE mppmin_a_int 1678 # include "mpp_allreduce_generic.h90" 1679 # undef ROUTINE_ALLREDUCE 1680 # undef DIM_1d 1681 # undef INTEGER_TYPE 1682 ! 1683 # define REAL_TYPE 1684 # define DIM_0d 1685 # define ROUTINE_ALLREDUCE mppmin_real 1686 # include "mpp_allreduce_generic.h90" 1687 # undef ROUTINE_ALLREDUCE 1688 # undef DIM_0d 1689 # define DIM_1d 1690 # define ROUTINE_ALLREDUCE mppmin_a_real 1691 # include "mpp_allreduce_generic.h90" 1692 # undef ROUTINE_ALLREDUCE 1693 # undef DIM_1d 1694 # undef REAL_TYPE 1695 # undef OPERATION_MIN 1696 1697 !!---------------------------------------------------------------------- 1698 !! *** mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real *** 1699 !! 1700 !! Global sum of 1D array or a variable (integer, real or complex) 1701 !!---------------------------------------------------------------------- 1702 !! 1703 # define OPERATION_SUM 1704 # define INTEGER_TYPE 1705 # define DIM_0d 1706 # define ROUTINE_ALLREDUCE mppsum_int 1707 # include "mpp_allreduce_generic.h90" 1708 # undef ROUTINE_ALLREDUCE 1709 # undef DIM_0d 1710 # define DIM_1d 1711 # define ROUTINE_ALLREDUCE mppsum_a_int 1712 # include "mpp_allreduce_generic.h90" 1713 # undef ROUTINE_ALLREDUCE 1714 # undef DIM_1d 1715 # undef INTEGER_TYPE 1716 ! 1717 # define REAL_TYPE 1718 # define DIM_0d 1719 # define ROUTINE_ALLREDUCE mppsum_real 1720 # include "mpp_allreduce_generic.h90" 1721 # undef ROUTINE_ALLREDUCE 1722 # undef DIM_0d 1723 # define DIM_1d 1724 # define ROUTINE_ALLREDUCE mppsum_a_real 1725 # include "mpp_allreduce_generic.h90" 1726 # undef ROUTINE_ALLREDUCE 1727 # undef DIM_1d 1728 # undef REAL_TYPE 1729 # undef OPERATION_SUM 1730 1731 # define OPERATION_SUM_DD 1732 # define COMPLEX_TYPE 1733 # define DIM_0d 1734 # define ROUTINE_ALLREDUCE mppsum_realdd 1735 # include "mpp_allreduce_generic.h90" 1736 # undef ROUTINE_ALLREDUCE 1737 # undef DIM_0d 1738 # define DIM_1d 1739 # define ROUTINE_ALLREDUCE mppsum_a_realdd 1740 # include "mpp_allreduce_generic.h90" 1741 # undef ROUTINE_ALLREDUCE 1742 # undef DIM_1d 1743 # undef COMPLEX_TYPE 1744 # undef OPERATION_SUM_DD 1745 1746 !!---------------------------------------------------------------------- 1747 !! *** mpp_minloc2d, mpp_minloc3d, mpp_maxloc2d, mpp_maxloc3d 1748 !! 1749 !!---------------------------------------------------------------------- 1750 !! 1751 # define OPERATION_MINLOC 1752 # define DIM_2d 1753 # define ROUTINE_LOC mpp_minloc2d 1754 # include "mpp_loc_generic.h90" 1755 # undef ROUTINE_LOC 1756 # undef DIM_2d 1757 # define DIM_3d 1758 # define ROUTINE_LOC mpp_minloc3d 1759 # include "mpp_loc_generic.h90" 1760 # undef ROUTINE_LOC 1761 # undef DIM_3d 1762 # undef OPERATION_MINLOC 1763 1764 # define OPERATION_MAXLOC 1765 # define DIM_2d 1766 # define ROUTINE_LOC mpp_maxloc2d 1767 # include "mpp_loc_generic.h90" 1768 # undef ROUTINE_LOC 1769 # undef DIM_2d 1770 # define DIM_3d 1771 # define ROUTINE_LOC mpp_maxloc3d 1772 # include "mpp_loc_generic.h90" 1773 # undef ROUTINE_LOC 1774 # undef DIM_3d 1775 # undef OPERATION_MAXLOC 1776 1777 SUBROUTINE mpp_delay_sum( cdname, cdelay, y_in, pout, ldlast, kcom ) 1778 CHARACTER(len=*), INTENT(in ) :: cdname ! name of the calling subroutine 1779 CHARACTER(len=*), INTENT(in ) :: cdelay ! name (used as id) of the delayed operation 1780 COMPLEX(wp), INTENT(in ), DIMENSION(:) :: y_in 1781 REAL(wp), INTENT( out), DIMENSION(:) :: pout 1782 LOGICAL, INTENT(in ) :: ldlast ! true if this is the last time we call this routine 1783 INTEGER, INTENT(in ), OPTIONAL :: kcom 1784 ! 1785 pout(:) = REAL(y_in(:), wp) 1786 END SUBROUTINE mpp_delay_sum 1787 1788 SUBROUTINE mpp_delay_max( cdname, cdelay, p_in, pout, ldlast, kcom ) 1789 CHARACTER(len=*), INTENT(in ) :: cdname ! name of the calling subroutine 1790 CHARACTER(len=*), INTENT(in ) :: cdelay ! name (used as id) of the delayed operation 1791 REAL(wp), INTENT(in ), DIMENSION(:) :: p_in 1792 REAL(wp), INTENT( out), DIMENSION(:) :: pout 1793 LOGICAL, INTENT(in ) :: ldlast ! true if this is the last time we call this routine 1794 INTEGER, INTENT(in ), OPTIONAL :: kcom 1795 ! 1796 pout(:) = p_in(:) 1797 END SUBROUTINE mpp_delay_max 1798 1799 SUBROUTINE mpp_delay_rcv( kid ) 1800 INTEGER,INTENT(in ) :: kid 1801 WRITE(*,*) 'mpp_delay_rcv: You should not have seen this print! error?', kid 1802 END SUBROUTINE mpp_delay_rcv 1803 1804 SUBROUTINE mppstop( ldfinal, ld_force_abort ) 1805 LOGICAL, OPTIONAL, INTENT(in) :: ldfinal ! source process number 1806 LOGICAL, OPTIONAL, INTENT(in) :: ld_force_abort ! source process number 1783 1807 STOP ! non MPP case, just stop the run 1784 1808 END SUBROUTINE mppstop 1785 1786 SUBROUTINE mpp_ini_ice( kcom, knum )1787 INTEGER :: kcom, knum1788 WRITE(*,*) 'mpp_ini_ice: You should not have seen this print! error?', kcom, knum1789 END SUBROUTINE mpp_ini_ice1790 1809 1791 1810 SUBROUTINE mpp_ini_znl( knum ) … … 1799 1818 END SUBROUTINE mpp_comm_free 1800 1819 1801 SUBROUTINE mppmax_real_multiple( ptab, kdim , kcom )1802 REAL, DIMENSION(:) :: ptab !1803 INTEGER :: kdim !1804 INTEGER, OPTIONAL :: kcom !1805 WRITE(*,*) 'mppmax_real_multiple: You should not have seen this print! error?', ptab(1), kdim1806 END SUBROUTINE mppmax_real_multiple1807 1808 1820 #endif 1809 1821 … … 1825 1837 ! 1826 1838 nstop = nstop + 1 1827 IF(lwp) THEN 1828 WRITE(numout,cform_err) 1829 IF( PRESENT(cd1 ) ) WRITE(numout,*) cd1 1830 IF( PRESENT(cd2 ) ) WRITE(numout,*) cd2 1831 IF( PRESENT(cd3 ) ) WRITE(numout,*) cd3 1832 IF( PRESENT(cd4 ) ) WRITE(numout,*) cd4 1833 IF( PRESENT(cd5 ) ) WRITE(numout,*) cd5 1834 IF( PRESENT(cd6 ) ) WRITE(numout,*) cd6 1835 IF( PRESENT(cd7 ) ) WRITE(numout,*) cd7 1836 IF( PRESENT(cd8 ) ) WRITE(numout,*) cd8 1837 IF( PRESENT(cd9 ) ) WRITE(numout,*) cd9 1838 IF( PRESENT(cd10) ) WRITE(numout,*) cd10 1839 ENDIF 1839 1840 ! force to open ocean.output file 1841 IF( numout == 6 ) CALL ctl_opn( numout, 'ocean.output', 'APPEND', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. ) 1842 1843 WRITE(numout,cform_err) 1844 IF( PRESENT(cd1 ) ) WRITE(numout,*) TRIM(cd1) 1845 IF( PRESENT(cd2 ) ) WRITE(numout,*) TRIM(cd2) 1846 IF( PRESENT(cd3 ) ) WRITE(numout,*) TRIM(cd3) 1847 IF( PRESENT(cd4 ) ) WRITE(numout,*) TRIM(cd4) 1848 IF( PRESENT(cd5 ) ) WRITE(numout,*) TRIM(cd5) 1849 IF( PRESENT(cd6 ) ) WRITE(numout,*) TRIM(cd6) 1850 IF( PRESENT(cd7 ) ) WRITE(numout,*) TRIM(cd7) 1851 IF( PRESENT(cd8 ) ) WRITE(numout,*) TRIM(cd8) 1852 IF( PRESENT(cd9 ) ) WRITE(numout,*) TRIM(cd9) 1853 IF( PRESENT(cd10) ) WRITE(numout,*) TRIM(cd10) 1854 1840 1855 CALL FLUSH(numout ) 1841 1856 IF( numstp /= -1 ) CALL FLUSH(numstp ) … … 1844 1859 ! 1845 1860 IF( cd1 == 'STOP' ) THEN 1846 IF(lwp)WRITE(numout,*) 'huge E-R-R-O-R : immediate stop'1847 CALL mppstop( )1861 WRITE(numout,*) 'huge E-R-R-O-R : immediate stop' 1862 CALL mppstop(ld_force_abort = .true.) 1848 1863 ENDIF 1849 1864 ! … … 1866 1881 IF(lwp) THEN 1867 1882 WRITE(numout,cform_war) 1868 IF( PRESENT(cd1 ) ) WRITE(numout,*) cd11869 IF( PRESENT(cd2 ) ) WRITE(numout,*) cd21870 IF( PRESENT(cd3 ) ) WRITE(numout,*) cd31871 IF( PRESENT(cd4 ) ) WRITE(numout,*) cd41872 IF( PRESENT(cd5 ) ) WRITE(numout,*) cd51873 IF( PRESENT(cd6 ) ) WRITE(numout,*) cd61874 IF( PRESENT(cd7 ) ) WRITE(numout,*) cd71875 IF( PRESENT(cd8 ) ) WRITE(numout,*) cd81876 IF( PRESENT(cd9 ) ) WRITE(numout,*) cd91877 IF( PRESENT(cd10) ) WRITE(numout,*) cd101883 IF( PRESENT(cd1 ) ) WRITE(numout,*) TRIM(cd1) 1884 IF( PRESENT(cd2 ) ) WRITE(numout,*) TRIM(cd2) 1885 IF( PRESENT(cd3 ) ) WRITE(numout,*) TRIM(cd3) 1886 IF( PRESENT(cd4 ) ) WRITE(numout,*) TRIM(cd4) 1887 IF( PRESENT(cd5 ) ) WRITE(numout,*) TRIM(cd5) 1888 IF( PRESENT(cd6 ) ) WRITE(numout,*) TRIM(cd6) 1889 IF( PRESENT(cd7 ) ) WRITE(numout,*) TRIM(cd7) 1890 IF( PRESENT(cd8 ) ) WRITE(numout,*) TRIM(cd8) 1891 IF( PRESENT(cd9 ) ) WRITE(numout,*) TRIM(cd9) 1892 IF( PRESENT(cd10) ) WRITE(numout,*) TRIM(cd10) 1878 1893 ENDIF 1879 1894 CALL FLUSH(numout) … … 1916 1931 knum=get_unit() 1917 1932 #endif 1933 IF( TRIM(cdfile) == '/dev/null' ) clfile = TRIM(cdfile) ! force the use of /dev/null 1918 1934 ! 1919 1935 iost=0 1920 IF( cdacce(1:6) == 'DIRECT' ) THEN 1921 OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh, ERR=100, IOSTAT=iost ) 1936 IF( cdacce(1:6) == 'DIRECT' ) THEN ! cdacce has always more than 6 characters 1937 OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh , ERR=100, IOSTAT=iost ) 1938 ELSE IF( TRIM(cdstat) == 'APPEND' ) THEN ! cdstat can have less than 6 characters 1939 OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS='UNKNOWN', POSITION='APPEND', ERR=100, IOSTAT=iost ) 1922 1940 ELSE 1923 OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat , ERR=100, IOSTAT=iost ) 1924 ENDIF 1941 OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat , ERR=100, IOSTAT=iost ) 1942 ENDIF 1943 IF( iost /= 0 .AND. TRIM(clfile) == '/dev/null' ) & ! for windows 1944 & OPEN(UNIT=knum,FILE='NUL', FORM=cdform, ACCESS=cdacce, STATUS=cdstat , ERR=100, IOSTAT=iost ) 1925 1945 IF( iost == 0 ) THEN 1926 1946 IF(ldwp) THEN 1927 WRITE(kout,*) ' file : ', clfile,' open ok'1947 WRITE(kout,*) ' file : ', TRIM(clfile),' open ok' 1928 1948 WRITE(kout,*) ' unit = ', knum 1929 1949 WRITE(kout,*) ' status = ', cdstat … … 1937 1957 IF(ldwp) THEN 1938 1958 WRITE(kout,*) 1939 WRITE(kout,*) ' ===>>>> : bad opening file: ', clfile1959 WRITE(kout,*) ' ===>>>> : bad opening file: ', TRIM(clfile) 1940 1960 WRITE(kout,*) ' ======= === ' 1941 1961 WRITE(kout,*) ' unit = ', knum … … 1948 1968 ELSE !!! Force writing to make sure we get the information - at least once - in this violent STOP!! 1949 1969 WRITE(*,*) 1950 WRITE(*,*) ' ===>>>> : bad opening file: ', clfile1970 WRITE(*,*) ' ===>>>> : bad opening file: ', TRIM(clfile) 1951 1971 WRITE(*,*) ' ======= === ' 1952 1972 WRITE(*,*) ' unit = ', knum
Note: See TracChangeset
for help on using the changeset viewer.