From 2799bb01914414b6f44a143cf01e373d5bda556a Mon Sep 17 00:00:00 2001 From: Antoine Cyril David Hoffmann <ahoffman@spcpc606.epfl.ch> Date: Mon, 15 Mar 2021 16:12:28 +0100 Subject: [PATCH] Adaptation for 2D parallel runs with mpi along p and kr --- src/basic_mod.F90 | 16 +++---- src/ppinit.F90 | 103 ++++++++++++++++++---------------------------- 2 files changed, 50 insertions(+), 69 deletions(-) diff --git a/src/basic_mod.F90 b/src/basic_mod.F90 index 5a7e4a80..15ca6353 100644 --- a/src/basic_mod.F90 +++ b/src/basic_mod.F90 @@ -10,10 +10,10 @@ MODULE basic real(dp) :: tmax = 100000.0 ! Maximum simulation time real(dp) :: dt = 1.0 ! Time step real(dp) :: time = 0 ! Current simulation time (Init from restart file) - + INTEGER :: comm0 ! Default communicator with a topology INTEGER :: commp, commr ! Communicators for 1-dim cartesian subgrids of comm0 - + INTEGER :: jobnum = 0 ! Job number INTEGER :: step = 0 ! Calculation step of this run INTEGER :: cstep = 0 ! Current step number (Init from restart file) @@ -24,8 +24,10 @@ MODULE basic INTEGER :: ierr ! flag for MPI error INTEGER :: my_id ! Rank in COMM_WORLD INTEGER :: num_procs ! number of MPI processes - INTEGER :: ncp, ncr ! Number of processes in p and r - INTEGER :: me_0, me_p, me_z ! Ranks in comm0, commp, commz + INTEGER :: num_procs_p, num_procs_kr ! Number of processes in p and r + INTEGER :: rank_0, rank_p, rank_r! Ranks in comm0, commp, commr + INTEGER :: nbr_L, nbr_R ! Left and right neighbours (along p) + INTEGER :: nbr_T, nbr_B ! Top and bottom neighbours (along kr) INTEGER :: iframe1d ! counting the number of times 1d datasets are outputed (for diagnose) INTEGER :: iframe2d ! counting the number of times 2d datasets are outputed (for diagnose) @@ -38,9 +40,9 @@ MODULE basic ! To measure computation time real :: start, finish - real(dp) :: t0_rhs, t0_adv_field, t0_poisson, t0_Sapj, t0_diag, t0_checkfield, t0_step - real(dp) :: t1_rhs, t1_adv_field, t1_poisson, t1_Sapj, t1_diag, t1_checkfield, t1_step - real(dp) :: tc_rhs, tc_adv_field, tc_poisson, tc_Sapj, tc_diag, tc_checkfield, tc_step + real(dp) :: t0_rhs, t0_adv_field, t0_poisson, t0_Sapj, t0_diag, t0_checkfield, t0_step, t0_comm + real(dp) :: t1_rhs, t1_adv_field, t1_poisson, t1_Sapj, t1_diag, t1_checkfield, t1_step, t1_comm + real(dp) :: tc_rhs, tc_adv_field, tc_poisson, tc_Sapj, tc_diag, tc_checkfield, tc_step, tc_comm real(dp):: maxruntime = 1e9 ! Maximum simulation CPU time INTERFACE allocate_array diff --git a/src/ppinit.F90 b/src/ppinit.F90 index 9e8b8632..68dc30bc 100644 --- a/src/ppinit.F90 +++ b/src/ppinit.F90 @@ -12,79 +12,58 @@ SUBROUTINE ppinit LOGICAL :: periods(ndims) = .FALSE., reorder=.FALSE. CHARACTER(len=32) :: str INTEGER :: nargs, i, l - INTEGER :: nghb_L, nghb_R ! Left and right neighbors along p - INTEGER :: source CALL MPI_INIT(ierr) - ! CALL MPI_INIT_THREAD(MPI_THREAD_SINGLE,version_prov,ierr) - ! CALL MPI_INIT_THREAD(MPI_THREAD_FUNNELED,version_prov,ierr) CALL MPI_COMM_RANK (MPI_COMM_WORLD, my_id, ierr) CALL MPI_COMM_SIZE (MPI_COMM_WORLD, num_procs, ierr) - - ! nargs = COMMAND_ARGUMENT_COUNT() - ! IF( nargs .NE. 0 .AND. nargs .NE. ndims ) THEN - ! IF(my_id .EQ. 0) WRITE(*, '(a,i4,a)') 'Number of arguments not equal to NDIMS =', ndims, '!' - ! CALL MPI_ABORT(MPI_COMM_WORLD, -1, ierr) - ! END IF - ! ! - ! IF( nargs .NE. 0 ) THEN - ! DO i=1,nargs - ! CALL GET_COMMAND_ARGUMENT(i, str, l, ierr) - ! READ(str(1:l),'(i3)') dims(i) - ! ncp = dims(1) ! Number of processes along p - ! ncr = dims(2) ! Number of processes along kr - ! END DO - ! IF( PRODUCT(dims) .NE. num_procs ) THEN - ! IF(my_id .EQ. 0) WRITE(*, '(a,i4,a,i4)') 'Product of dims: ', PRODUCT(dims), " is not consistent WITH NPROCS=",num_procs - ! CALL MPI_ABORT(MPI_COMM_WORLD, -2, ierr) - ! END IF - ! ELSE - ! CALL MPI_DIMS_CREATE(num_procs, ndims, dims, ierr) - ! END IF - ! ! - ! !periodicity in p - ! periods(1)=.FALSE. - ! !periodicity in kr - ! periods(2)=.FALSE. - ! CALL MPI_CART_CREATE(MPI_COMM_WORLD, ndims, dims, periods, reorder, comm0, ierr) + nargs = COMMAND_ARGUMENT_COUNT() + IF( nargs .NE. 0 .AND. nargs .NE. ndims ) THEN + IF(my_id .EQ. 0) WRITE(*, '(a,i4,a)') 'Number of arguments not equal to NDIMS =', ndims, '!' + CALL MPI_ABORT(MPI_COMM_WORLD, -1, ierr) + END IF + ! + IF( nargs .NE. 0 ) THEN + DO i=1,nargs + CALL GET_COMMAND_ARGUMENT(i, str, l, ierr) + READ(str(1:l),'(i3)') dims(i) + END DO + IF( PRODUCT(dims) .NE. num_procs ) THEN + IF(my_id .EQ. 0) WRITE(*, '(a,i4,a,i4)') 'Product of dims: ', PRODUCT(dims), " is not consistent WITH NPROCS=",num_procs + CALL MPI_ABORT(MPI_COMM_WORLD, -2, ierr) + END IF + ELSE + ! CALL MPI_DIMS_CREATE(num_procs, ndims, dims, ierr) + dims(1) = 1 + dims(2) = num_procs + END IF - ! CALL MPI_COMM_RANK(comm0, me_0, ierr) + num_procs_p = dims(1) ! Number of processes along p + num_procs_kr = dims(2) ! Number of processes along kr - ! CALL MPI_CART_COORDS(comm0,me_0,ndims,coords,ierr) + ! + !periodicity in p + periods(1)=.FALSE. + !periodicity in kr + periods(2)=.FALSE. - ! DO i=0,num_procs-1 - ! CALL mpi_barrier(MPI_COMM_WORLD, ierr) - ! IF (my_id .EQ. i) THEN - ! WRITE(*,*) 'The coords of process ',me_0,' are: ', coords + CALL MPI_CART_CREATE(MPI_COMM_WORLD, ndims, dims, periods, reorder, comm0, ierr) - ! IF( coords(1) .GT. 0 ) THEN - ! CALL MPI_CART_SHIFT(comm0, 0, -1, source , nghb_L , ierr) - ! CALL MPI_CART_COORDS(comm0,nghb_L,ndims,coords_L,ierr) - ! WRITE(*,*) coords_L,' are L from ', coords - ! ELSE - ! nghb_L = MPI_PROC_NULL - ! ENDIF + CALL MPI_COMM_RANK(comm0, rank_0, ierr) + CALL MPI_CART_COORDS(comm0,rank_0,ndims,coords,ierr) - ! IF( coords(1) .LT. dims(1)-1 ) THEN - ! CALL MPI_CART_SHIFT(comm0, 0, +1, source , nghb_R , ierr) - ! CALL MPI_CART_COORDS(comm0,nghb_R,ndims,coords_R,ierr) - ! WRITE(*,*) coords_R,' are R from ', coords - ! ELSE - ! nghb_R = MPI_PROC_NULL - ! ENDIF - ! ENDIF - ! ENDDO - - ! ! - ! ! Partitions 2-dim cartesian topology of comm0 into 1-dim cartesian subgrids - ! ! - ! CALL MPI_CART_SUB (comm0, (/.TRUE.,.FALSE./), commp, ierr) - ! CALL MPI_CART_SUB (comm0, (/.FALSE.,.TRUE./), commr, ierr) - - ! CALL MPI_COMM_RANK(commp, me_p, ierr) - ! CALL MPI_COMM_RANK(commr, me_r, ierr) + ! + ! Partitions 2-dim cartesian topology of comm0 into 1-dim cartesian subgrids + ! + CALL MPI_CART_SUB (comm0, (/.TRUE.,.FALSE./), commp, ierr) + CALL MPI_CART_SUB (comm0, (/.FALSE.,.TRUE./), commr, ierr) + ! Find id inside the sub communicators + CALL MPI_COMM_RANK(commp, rank_p, ierr) + CALL MPI_COMM_RANK(commr, rank_r, ierr) + ! Find neighbours + CALL MPI_CART_SHIFT(comm0, 0, 1, nbr_L, nbr_R, ierr) + CALL MPI_CART_SHIFT(comm0, 1, 1, nbr_B, nbr_T, ierr) END SUBROUTINE ppinit -- GitLab