1! 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt. 6! See Docs/Contributors.txt for a list of contributors. 7! 8 module m_mpi_utils 9 10 use precision, only: dp, sp 11 use sys, only: die 12#ifdef MPI 13 use mpi_siesta 14 implicit none 15 integer, private :: MPIerror 16#else 17 implicit none 18#endif 19 public :: globalize_max, globalize_sum, broadcast 20 public :: globalize_or, globalize_min 21 public :: barrier 22 private 23 24 interface globalize_max 25 module procedure globalize_max_dp, globalize_max_int 26 end interface 27 interface globalize_min 28 module procedure globalize_min_dp, globalize_min_int 29 end interface 30 interface globalize_sum 31 module procedure globalize_sum_dp, globalize_sum_int 32 module procedure globalize_sum_v_dp 33 module procedure globalize_sum_vv_dp 34 end interface 35 interface broadcast 36 module procedure broadcast_dp, broadcast_int, broadcast_logical 37 module procedure broadcast_sp, broadcast_char 38 module procedure broadcast_v_dp, broadcast_v_int 39 module procedure broadcast_vv_dp, broadcast_vv_int 40 module procedure broadcast_vvv_dp, broadcast_vvv_int 41 end interface 42 43 CONTAINS 44 45 subroutine barrier() 46#ifdef MPI 47 call MPI_Barrier(MPI_Comm_World,MPIerror) 48#endif 49 end subroutine barrier 50 51 subroutine broadcast_dp(scalar,comm) 52 real(dp), intent(inout) :: scalar 53 integer, intent(in), optional :: comm 54 55#ifdef MPI 56 integer :: mpi_comm 57 mpi_comm = mpi_comm_world 58 if (present(comm)) then 59 mpi_comm = comm 60 endif 61 call MPI_Bcast(scalar,1,MPI_double_precision,0, 62 $ MPI_Comm,MPIerror) 63#endif 64 end subroutine broadcast_dp 65 66 subroutine broadcast_v_dp(a,comm) 67 real(dp), dimension(:), intent(inout) :: a 68 integer, intent(in), optional :: comm 69 70#ifdef MPI 71 integer :: mpi_comm 72 mpi_comm = mpi_comm_world 73 if (present(comm)) then 74 mpi_comm = comm 75 endif 76 call MPI_Bcast(a,size(a),MPI_double_precision,0, 77 $ MPI_Comm,MPIerror) 78#endif 79 end subroutine broadcast_v_dp 80 81 subroutine broadcast_vv_dp(a,comm) 82 real(dp), dimension(:,:), intent(inout) :: a 83!! Only for contiguous array sections !! 84 integer, intent(in), optional :: comm 85 86#ifdef MPI 87 integer :: mpi_comm 88 mpi_comm = mpi_comm_world 89 if (present(comm)) then 90 mpi_comm = comm 91 endif 92 call MPI_Bcast(a(1,1),size(a),MPI_double_precision,0, 93 $ MPI_Comm,MPIerror) 94#endif 95 end subroutine broadcast_vv_dp 96 97 subroutine broadcast_vvv_dp(a,comm) 98 real(dp), dimension(:,:,:), intent(inout) :: a 99!! Only for contiguous array sections !! 100 integer, intent(in), optional :: comm 101 102#ifdef MPI 103 integer :: mpi_comm 104 mpi_comm = mpi_comm_world 105 if (present(comm)) then 106 mpi_comm = comm 107 endif 108 call MPI_Bcast(a(1,1,1),size(a),MPI_double_precision,0, 109 $ MPI_Comm,MPIerror) 110#endif 111 end subroutine broadcast_vvv_dp 112 113 subroutine broadcast_sp(scalar,comm) 114 real(sp), intent(inout) :: scalar 115 integer, intent(in), optional :: comm 116 117#ifdef MPI 118 integer :: mpi_comm 119 mpi_comm = mpi_comm_world 120 if (present(comm)) then 121 mpi_comm = comm 122 endif 123 call MPI_Bcast(scalar,1,MPI_real,0, 124 $ MPI_Comm,MPIerror) 125#endif 126 end subroutine broadcast_sp 127 128 subroutine broadcast_int(scalar,comm) 129 integer, intent(inout) :: scalar 130 integer, intent(in), optional :: comm 131 132#ifdef MPI 133 integer :: mpi_comm 134 mpi_comm = mpi_comm_world 135 if (present(comm)) then 136 mpi_comm = comm 137 endif 138 call MPI_Bcast(scalar,1,MPI_Integer,0, 139 $ MPI_Comm,MPIerror) 140#endif 141 end subroutine broadcast_int 142 143 subroutine broadcast_v_int(a,comm) 144 integer, dimension(:), intent(inout) :: a 145 integer, intent(in), optional :: comm 146 147#ifdef MPI 148 integer :: mpi_comm 149 mpi_comm = mpi_comm_world 150 if (present(comm)) then 151 mpi_comm = comm 152 endif 153 call MPI_Bcast(a,size(a),MPI_integer,0, 154 $ MPI_Comm,MPIerror) 155#endif 156 end subroutine broadcast_v_int 157 158 subroutine broadcast_vv_int(a,comm) 159 integer, dimension(:,:), intent(inout) :: a 160 integer, intent(in), optional :: comm 161 162#ifdef MPI 163 integer :: mpi_comm 164 mpi_comm = mpi_comm_world 165 if (present(comm)) then 166 mpi_comm = comm 167 endif 168 call MPI_Bcast(a(1,1),size(a),MPI_integer,0, 169 $ MPI_Comm,MPIerror) 170#endif 171 end subroutine broadcast_vv_int 172 173 subroutine broadcast_vvv_int(a,comm) 174 integer, dimension(:,:,:), intent(inout) :: a 175 integer, intent(in), optional :: comm 176 177#ifdef MPI 178 integer :: mpi_comm 179 mpi_comm = mpi_comm_world 180 if (present(comm)) then 181 mpi_comm = comm 182 endif 183 call MPI_Bcast(a(1,1,1),size(a),MPI_integer,0, 184 $ MPI_Comm,MPIerror) 185#endif 186 end subroutine broadcast_vvv_int 187 188 subroutine broadcast_char(str,comm) 189 character(len=*), intent(inout) :: str 190 integer, intent(in), optional :: comm 191 192#ifdef MPI 193 integer :: mpi_comm 194 mpi_comm = mpi_comm_world 195 if (present(comm)) then 196 mpi_comm = comm 197 endif 198 call MPI_Bcast(str,len(str),MPI_Character,0, 199 $ MPI_Comm,MPIerror) 200#endif 201 end subroutine broadcast_char 202 203 subroutine broadcast_logical(scalar,comm) 204 logical, intent(inout) :: scalar 205 integer, intent(in), optional :: comm 206 207#ifdef MPI 208 integer :: mpi_comm 209 mpi_comm = mpi_comm_world 210 if (present(comm)) then 211 mpi_comm = comm 212 endif 213 call MPI_Bcast(scalar,1,MPI_Logical,0, 214 $ MPI_Comm,MPIerror) 215#endif 216 end subroutine broadcast_logical 217 218 219 subroutine Globalize_or(local,global,comm) 220 logical, intent(in) :: local 221 logical, intent(out) :: global 222 integer, intent(in), optional :: comm 223 224#ifdef MPI 225 integer :: mpi_comm 226 mpi_comm = mpi_comm_world 227 if (present(comm)) then 228 mpi_comm = comm 229 endif 230 call MPI_AllReduce(local,global,1,MPI_Logical, 231 $ MPI_LOR,MPI_Comm,MPIerror) 232#else 233 global = local 234#endif 235 end subroutine Globalize_or 236 237 subroutine Globalize_sum_dp(local,global,comm) 238 real(dp), intent(in) :: local 239 real(dp), intent(out) :: global 240 integer, intent(in), optional :: comm 241 242#ifdef MPI 243 integer :: mpi_comm 244 mpi_comm = mpi_comm_world 245 if (present(comm)) then 246 mpi_comm = comm 247 endif 248 call MPI_AllReduce(local,global,1,MPI_double_precision, 249 $ MPI_sum,MPI_Comm,MPIerror) 250#else 251 global = local 252#endif 253 end subroutine Globalize_sum_dp 254 255 subroutine Globalize_sum_v_dp(local,global,comm) 256 real(dp), intent(in), dimension(:) :: local 257 real(dp), intent(out), dimension(:) :: global 258 integer, intent(in), optional :: comm 259 260 integer :: n 261 integer :: mpi_comm 262 263 n = size(local) 264 if ( n /= size(global)) 265 $ call die("Globalize_sum_v_dp error") 266 267#ifdef MPI 268 mpi_comm = mpi_comm_world 269 if (present(comm)) then 270 mpi_comm = comm 271 endif 272 call MPI_AllReduce(local,global,n,MPI_double_precision, 273 $ MPI_sum,MPI_Comm,MPIerror) 274#else 275 global = local 276#endif 277 end subroutine Globalize_sum_v_dp 278 279 subroutine Globalize_sum_vv_dp(local,global,comm) 280 real(dp), intent(in), dimension(:,:) :: local 281 real(dp), intent(out), dimension(:,:) :: global 282 integer, intent(in), optional :: comm 283 284 integer :: n 285 integer :: mpi_comm 286 287 n = size(local) 288 if ( n /= size(global)) 289 $ call die("Globalize_sum_v_dp error") 290 291#ifdef MPI 292 mpi_comm = mpi_comm_world 293 if (present(comm)) then 294 mpi_comm = comm 295 endif 296 call MPI_AllReduce(local(1,1),global(1,1),n, 297 $ MPI_double_precision, 298 $ MPI_sum,MPI_Comm,MPIerror) 299 300#else 301 global = local 302#endif 303 end subroutine Globalize_sum_vv_dp 304 305 subroutine Globalize_sum_int(local,global,comm) 306 integer, intent(in) :: local 307 integer, intent(out) :: global 308 integer, intent(in), optional :: comm 309 310 integer :: mpi_comm 311#ifdef MPI 312 mpi_comm = mpi_comm_world 313 if (present(comm)) then 314 mpi_comm = comm 315 endif 316 call MPI_AllReduce(local,global,1,MPI_integer, 317 $ MPI_sum,MPI_Comm,MPIerror) 318#else 319 global = local 320#endif 321 end subroutine Globalize_sum_int 322 323 subroutine Globalize_max_dp(local,global,comm) 324 real(dp), intent(in) :: local 325 real(dp), intent(out) :: global 326 integer, intent(in), optional :: comm 327 328#ifdef MPI 329 integer :: mpi_comm 330 mpi_comm = mpi_comm_world 331 if (present(comm)) then 332 mpi_comm = comm 333 endif 334 call MPI_AllReduce(local,global,1,MPI_double_precision, 335 $ MPI_max,MPI_Comm,MPIerror) 336#else 337 global = local 338#endif 339 end subroutine Globalize_max_dp 340 341 subroutine Globalize_max_int(local,global,comm) 342 integer, intent(in) :: local 343 integer, intent(out) :: global 344 integer, intent(in), optional :: comm 345 346#ifdef MPI 347 integer :: mpi_comm 348 mpi_comm = mpi_comm_world 349 if (present(comm)) then 350 mpi_comm = comm 351 endif 352 call MPI_AllReduce(local,global,1,MPI_integer, 353 $ MPI_max,MPI_Comm,MPIerror) 354#else 355 global = local 356#endif 357 end subroutine Globalize_max_int 358 359 subroutine Globalize_min_dp(local,global,comm) 360 real(dp), intent(in) :: local 361 real(dp), intent(out) :: global 362 integer, intent(in), optional :: comm 363 364#ifdef MPI 365 integer :: mpi_comm 366 mpi_comm = mpi_comm_world 367 if (present(comm)) then 368 mpi_comm = comm 369 endif 370 call MPI_AllReduce(local,global,1,MPI_double_precision, 371 $ MPI_min,MPI_Comm,MPIerror) 372#else 373 global = local 374#endif 375 end subroutine Globalize_min_dp 376 377 subroutine Globalize_min_int(local,global,comm) 378 integer, intent(in) :: local 379 integer, intent(out) :: global 380 integer, intent(in), optional :: comm 381 382#ifdef MPI 383 integer :: mpi_comm 384 mpi_comm = mpi_comm_world 385 if (present(comm)) then 386 mpi_comm = comm 387 endif 388 call MPI_AllReduce(local,global,1,MPI_integer, 389 $ MPI_min,MPI_Comm,MPIerror) 390#else 391 global = local 392#endif 393 end subroutine Globalize_min_int 394 395 end module m_mpi_utils 396 397