[963] | 1 | MODULE transfert_mpi_collectives_mod |
---|
[1055] | 2 | IMPLICIT NONE |
---|
[963] | 3 | |
---|
| 4 | INTERFACE bcast_mpi |
---|
[1055] | 5 | MODULE PROCEDURE bcast_mpi_c |
---|
| 6 | MODULE PROCEDURE bcast_mpi_i, bcast_mpi_i1, bcast_mpi_i2, bcast_mpi_i3, bcast_mpi_i4 |
---|
| 7 | MODULE PROCEDURE bcast_mpi_r, bcast_mpi_r1, bcast_mpi_r2, bcast_mpi_r3, bcast_mpi_r4 |
---|
| 8 | MODULE PROCEDURE bcast_mpi_l, bcast_mpi_l1, bcast_mpi_l2, bcast_mpi_l3, bcast_mpi_l4 |
---|
[963] | 9 | END INTERFACE |
---|
| 10 | |
---|
| 11 | CONTAINS |
---|
| 12 | |
---|
| 13 | SUBROUTINE gather_field(field_loc,field_glo) |
---|
[1055] | 14 | USE field_mod |
---|
| 15 | USE domain_mod |
---|
| 16 | USE mpi_mod |
---|
| 17 | USE mpipara |
---|
[963] | 18 | TYPE(t_field),POINTER :: field_loc(:) |
---|
| 19 | TYPE(t_field),POINTER :: field_glo(:) |
---|
| 20 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 21 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 22 | INTEGER :: ireq,nreq |
---|
| 23 | INTEGER :: ind_glo,ind_loc |
---|
| 24 | |
---|
| 25 | IF (.NOT. using_mpi) THEN |
---|
[1055] | 26 | DO ind_loc=1,ndomain |
---|
| 27 | field_glo(ind_loc)%rval4d=field_loc(ind_loc)%rval4d |
---|
| 28 | ENDDO |
---|
[963] | 29 | |
---|
| 30 | ELSE |
---|
[1055] | 31 | nreq=ndomain |
---|
| 32 | IF (mpi_rank==0) nreq=nreq+ndomain_glo |
---|
| 33 | ALLOCATE(mpi_req(nreq)) |
---|
| 34 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
[963] | 35 | |
---|
[1055] | 36 | ireq=0 |
---|
| 37 | IF (mpi_rank==0) THEN |
---|
| 38 | DO ind_glo=1,ndomain_glo |
---|
| 39 | ireq=ireq+1 |
---|
| 40 | CALL MPI_IRECV(field_glo(ind_glo)%rval4d,size(field_glo(ind_glo)%rval4d) , MPI_REAL8 , & |
---|
| 41 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 42 | ENDDO |
---|
| 43 | ENDIF |
---|
[963] | 44 | |
---|
[1055] | 45 | DO ind_loc=1,ndomain |
---|
[963] | 46 | ireq=ireq+1 |
---|
| 47 | CALL MPI_ISEND(field_loc(ind_loc)%rval4d,size(field_loc(ind_loc)%rval4d) , MPI_REAL8 , & |
---|
[1055] | 48 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 49 | ENDDO |
---|
[963] | 50 | |
---|
[1055] | 51 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
[963] | 52 | |
---|
| 53 | ENDIF |
---|
| 54 | |
---|
| 55 | END SUBROUTINE gather_field |
---|
| 56 | |
---|
| 57 | SUBROUTINE bcast_field(field_glo) |
---|
[1055] | 58 | USE field_mod |
---|
| 59 | USE domain_mod |
---|
| 60 | USE mpi_mod |
---|
| 61 | USE mpipara |
---|
[963] | 62 | TYPE(t_field),POINTER :: field_glo(:) |
---|
| 63 | INTEGER :: ind_glo |
---|
| 64 | |
---|
[1055] | 65 | IF (using_mpi) THEN |
---|
| 66 | DO ind_glo=1,ndomain_glo |
---|
| 67 | CALL MPI_BCAST(field_glo(ind_glo)%rval4d, size(field_glo(ind_glo)%rval4d) , MPI_REAL8, 0, comm_icosa, ierr) |
---|
| 68 | ENDDO |
---|
| 69 | ENDIF |
---|
[963] | 70 | |
---|
| 71 | END SUBROUTINE bcast_field |
---|
| 72 | |
---|
| 73 | SUBROUTINE scatter_field(field_glo,field_loc) |
---|
[1055] | 74 | USE field_mod |
---|
| 75 | USE domain_mod |
---|
| 76 | USE mpi_mod |
---|
| 77 | USE mpipara |
---|
[963] | 78 | TYPE(t_field),POINTER :: field_glo(:) |
---|
| 79 | TYPE(t_field),POINTER :: field_loc(:) |
---|
| 80 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 81 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 82 | INTEGER :: ireq,nreq |
---|
| 83 | INTEGER :: ind_glo,ind_loc |
---|
| 84 | |
---|
| 85 | IF (.NOT. using_mpi) THEN |
---|
[1055] | 86 | DO ind_loc=1,ndomain |
---|
| 87 | field_loc(ind_loc)%rval4d=field_glo(ind_loc)%rval4d |
---|
| 88 | ENDDO |
---|
[963] | 89 | |
---|
| 90 | ELSE |
---|
[1055] | 91 | nreq=ndomain |
---|
| 92 | IF (mpi_rank==0) nreq=nreq+ndomain_glo |
---|
| 93 | ALLOCATE(mpi_req(nreq)) |
---|
| 94 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
[963] | 95 | |
---|
[1055] | 96 | ireq=0 |
---|
| 97 | IF (mpi_rank==0) THEN |
---|
| 98 | DO ind_glo=1,ndomain_glo |
---|
| 99 | ireq=ireq+1 |
---|
| 100 | CALL MPI_ISEND(field_glo(ind_glo)%rval4d,size(field_glo(ind_glo)%rval4d) , MPI_REAL8 , & |
---|
| 101 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 102 | ENDDO |
---|
| 103 | ENDIF |
---|
[963] | 104 | |
---|
[1055] | 105 | DO ind_loc=1,ndomain |
---|
[963] | 106 | ireq=ireq+1 |
---|
| 107 | CALL MPI_IRECV(field_loc(ind_loc)%rval4d,size(field_loc(ind_loc)%rval4d) , MPI_REAL8 , & |
---|
[1055] | 108 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 109 | ENDDO |
---|
[963] | 110 | |
---|
[1055] | 111 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
[963] | 112 | |
---|
| 113 | ENDIF |
---|
| 114 | |
---|
| 115 | END SUBROUTINE scatter_field |
---|
| 116 | |
---|
[1055] | 117 | !=================== Broadcast routines for strings ==================! |
---|
[963] | 118 | |
---|
[1055] | 119 | SUBROUTINE bcast_mpi_cgen(var,nb) |
---|
| 120 | USE mpi_mod |
---|
| 121 | USE mpipara |
---|
| 122 | CHARACTER(LEN=*),INTENT(INOUT) :: Var |
---|
| 123 | INTEGER,INTENT(IN) :: nb |
---|
| 124 | IF (.NOT. using_mpi) RETURN |
---|
| 125 | CALL MPI_BCAST(Var,nb,MPI_CHARACTER,mpi_master,comm_icosa,ierr) |
---|
| 126 | END SUBROUTINE bcast_mpi_cgen |
---|
[963] | 127 | |
---|
| 128 | SUBROUTINE bcast_mpi_c(var1) |
---|
| 129 | CHARACTER(LEN=*),INTENT(INOUT) :: Var1 |
---|
| 130 | CALL bcast_mpi_cgen(Var1,len(Var1)) |
---|
| 131 | END SUBROUTINE bcast_mpi_c |
---|
| 132 | |
---|
[1055] | 133 | !============= Broadcast routines for INTEGER scalars and arrays ============! |
---|
[963] | 134 | |
---|
[1055] | 135 | SUBROUTINE bcast_mpi_igen(var,nb) |
---|
| 136 | USE mpi_mod |
---|
| 137 | USE mpipara |
---|
| 138 | INTEGER, INTENT(IN) :: nb |
---|
| 139 | INTEGER, DIMENSION(nb), INTENT(INOUT) :: var |
---|
| 140 | IF (using_mpi) CALL MPI_BCAST(Var, nb, MPI_INTEGER, mpi_master, comm_icosa, ierr) |
---|
| 141 | END SUBROUTINE bcast_mpi_igen |
---|
| 142 | |
---|
[963] | 143 | SUBROUTINE bcast_mpi_i(var) |
---|
[1055] | 144 | USE mpipara |
---|
| 145 | INTEGER, INTENT(INOUT) :: var |
---|
| 146 | INTEGER :: var_tmp(1) |
---|
[963] | 147 | IF (is_mpi_master) var_tmp(1)=var |
---|
[1055] | 148 | CALL bcast_mpi_igen(var_tmp,1) |
---|
[963] | 149 | var=var_tmp(1) |
---|
| 150 | END SUBROUTINE bcast_mpi_i |
---|
| 151 | |
---|
| 152 | SUBROUTINE bcast_mpi_i1(var) |
---|
[1055] | 153 | INTEGER, INTENT(INOUT) :: var(:) |
---|
| 154 | CALL bcast_mpi_igen(var,size(var)) |
---|
[963] | 155 | END SUBROUTINE bcast_mpi_i1 |
---|
| 156 | |
---|
| 157 | SUBROUTINE bcast_mpi_i2(var) |
---|
[1055] | 158 | INTEGER, INTENT(INOUT) :: var(:,:) |
---|
| 159 | CALL bcast_mpi_igen(var,size(var)) |
---|
[963] | 160 | END SUBROUTINE bcast_mpi_i2 |
---|
| 161 | |
---|
| 162 | SUBROUTINE bcast_mpi_i3(var) |
---|
[1055] | 163 | INTEGER, INTENT(INOUT) :: var(:,:,:) |
---|
| 164 | CALL bcast_mpi_igen(var,size(var)) |
---|
[963] | 165 | END SUBROUTINE bcast_mpi_i3 |
---|
| 166 | |
---|
| 167 | SUBROUTINE bcast_mpi_i4(var) |
---|
[1055] | 168 | INTEGER, INTENT(INOUT) :: var(:,:,:,:) |
---|
| 169 | CALL bcast_mpi_igen(var,size(var)) |
---|
[963] | 170 | END SUBROUTINE bcast_mpi_i4 |
---|
| 171 | |
---|
[1055] | 172 | !============= Broadcast routines for REAL scalars and arrays ============! |
---|
[963] | 173 | |
---|
[1055] | 174 | SUBROUTINE bcast_mpi_rgen(var,nb) |
---|
| 175 | USE mpi_mod |
---|
| 176 | USE mpipara |
---|
| 177 | INTEGER, INTENT(IN) :: nb |
---|
| 178 | REAL, DIMENSION(nb), INTENT(INOUT) :: var |
---|
| 179 | IF (using_mpi) CALL MPI_BCAST(Var, nb, MPI_REAL8, mpi_master, comm_icosa, ierr) |
---|
| 180 | END SUBROUTINE bcast_mpi_rgen |
---|
[963] | 181 | |
---|
| 182 | SUBROUTINE bcast_mpi_r(var) |
---|
[1055] | 183 | USE mpipara |
---|
| 184 | REAL, INTENT(INOUT) :: var |
---|
| 185 | REAL :: var_tmp(1) |
---|
[963] | 186 | IF (is_mpi_master) var_tmp(1)=var |
---|
[1055] | 187 | CALL bcast_mpi_rgen(var_tmp,1) |
---|
[963] | 188 | var=var_tmp(1) |
---|
| 189 | END SUBROUTINE bcast_mpi_r |
---|
| 190 | |
---|
| 191 | SUBROUTINE bcast_mpi_r1(var) |
---|
[1055] | 192 | REAL, INTENT(INOUT) :: var(:) |
---|
| 193 | CALL bcast_mpi_rgen(var,size(var)) |
---|
[963] | 194 | END SUBROUTINE bcast_mpi_r1 |
---|
| 195 | |
---|
| 196 | SUBROUTINE bcast_mpi_r2(var) |
---|
[1055] | 197 | REAL, INTENT(INOUT) :: var(:,:) |
---|
| 198 | CALL bcast_mpi_rgen(var,size(var)) |
---|
[963] | 199 | END SUBROUTINE bcast_mpi_r2 |
---|
| 200 | |
---|
| 201 | SUBROUTINE bcast_mpi_r3(var) |
---|
[1055] | 202 | REAL, INTENT(INOUT) :: var(:,:,:) |
---|
| 203 | CALL bcast_mpi_rgen(var,size(var)) |
---|
[963] | 204 | END SUBROUTINE bcast_mpi_r3 |
---|
| 205 | |
---|
| 206 | SUBROUTINE bcast_mpi_r4(var) |
---|
[1055] | 207 | REAL, INTENT(INOUT) :: var(:,:,:,:) |
---|
| 208 | CALL bcast_mpi_rgen(var,size(var)) |
---|
| 209 | END SUBROUTINE bcast_mpi_r4 |
---|
[963] | 210 | |
---|
[1055] | 211 | !============= Broadcast routines for LOGICAL scalars and arrays ============! |
---|
[963] | 212 | |
---|
[1055] | 213 | SUBROUTINE bcast_mpi_lgen(var,nb) |
---|
| 214 | USE mpi_mod |
---|
| 215 | USE mpipara |
---|
| 216 | INTEGER, INTENT(IN) :: nb |
---|
| 217 | LOGICAL, DIMENSION(nb), INTENT(INOUT) :: var |
---|
| 218 | IF (using_mpi) CALL MPI_BCAST(Var, nb, MPI_LOGICAL, mpi_master, comm_icosa, ierr) |
---|
| 219 | END SUBROUTINE bcast_mpi_lgen |
---|
[963] | 220 | |
---|
| 221 | SUBROUTINE bcast_mpi_l(var) |
---|
[1055] | 222 | USE mpipara |
---|
| 223 | LOGICAL, INTENT(INOUT) :: var |
---|
| 224 | LOGICAL :: var_tmp(1) |
---|
[963] | 225 | IF (is_mpi_master) var_tmp(1)=var |
---|
[1055] | 226 | CALL bcast_mpi_lgen(var_tmp,1) |
---|
[963] | 227 | var=var_tmp(1) |
---|
| 228 | END SUBROUTINE bcast_mpi_l |
---|
| 229 | |
---|
| 230 | SUBROUTINE bcast_mpi_l1(var) |
---|
[1055] | 231 | LOGICAL, INTENT(INOUT) :: var(:) |
---|
| 232 | CALL bcast_mpi_lgen(var,size(var)) |
---|
[963] | 233 | END SUBROUTINE bcast_mpi_l1 |
---|
| 234 | |
---|
| 235 | SUBROUTINE bcast_mpi_l2(var) |
---|
[1055] | 236 | LOGICAL, INTENT(INOUT) :: var(:,:) |
---|
| 237 | CALL bcast_mpi_lgen(var,size(var)) |
---|
[963] | 238 | END SUBROUTINE bcast_mpi_l2 |
---|
| 239 | |
---|
| 240 | SUBROUTINE bcast_mpi_l3(var) |
---|
[1055] | 241 | LOGICAL, INTENT(INOUT) :: var(:,:,:) |
---|
| 242 | CALL bcast_mpi_lgen(var,size(var)) |
---|
[963] | 243 | END SUBROUTINE bcast_mpi_l3 |
---|
| 244 | |
---|
| 245 | SUBROUTINE bcast_mpi_l4(var) |
---|
[1055] | 246 | LOGICAL, INTENT(INOUT) :: var(:,:,:,:) |
---|
| 247 | CALL bcast_mpi_lgen(var,size(var)) |
---|
[963] | 248 | END SUBROUTINE bcast_mpi_l4 |
---|
| 249 | |
---|
| 250 | END MODULE transfert_mpi_collectives_mod |
---|