FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_comm_f.F90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
5
7contains
8
9 !C
10 !C***
11 !C*** hecmw_barrier
12 !C***
13 !C
14 subroutine hecmw_barrier (hecMESH)
15 use hecmw_util
16 implicit none
17 type (hecmwST_local_mesh) :: hecMESH
18#ifndef HECMW_SERIAL
19 integer(kind=kint):: ierr
20 call mpi_barrier (hecmesh%MPI_COMM, ierr)
21#endif
22 end subroutine hecmw_barrier
23
24 subroutine hecmw_scatterv_dp(sbuf, sc, disp, rbuf, rc, root, comm)
25 use hecmw_util
26 implicit none
27 integer(kind=kint) :: sc !send counts
28 double precision :: sbuf(sc) !send buffer
29 integer(kind=kint) :: disp !displacement
30 integer(kind=kint) :: rc !receive counts
31 double precision :: rbuf(rc) !receive buffer
32 integer(kind=kint) :: root
33 integer(kind=kint) :: comm
34#ifndef HECMW_SERIAL
35 integer(kind=kint) :: ierr
36 call mpi_scatterv( sbuf, sc, disp, mpi_double_precision, &
37 rbuf, rc, mpi_double_precision, &
38 root, comm, ierr )
39#else
40 rbuf(1) = sbuf(1)
41#endif
42 end subroutine hecmw_scatterv_dp
43
44#ifndef HECMW_SERIAL
45 function hecmw_operation_hec2mpi(operation)
46 use hecmw_util
47 implicit none
48 integer(kind=kint) :: hecmw_operation_hec2mpi
49 integer(kind=kint) :: operation
51 if (operation == hecmw_sum) then
53 elseif (operation == hecmw_prod) then
55 elseif (operation == hecmw_max) then
57 elseif (operation == hecmw_min) then
59 endif
60 end function hecmw_operation_hec2mpi
61#endif
62
63 subroutine hecmw_scatter_int_1(sbuf, rval, root, comm)
64 use hecmw_util
65 implicit none
66 integer(kind=kint) :: sbuf(*) !send buffer
67 integer(kind=kint) :: rval !receive value
68 integer(kind=kint) :: root
69 integer(kind=kint) :: comm
70#ifndef HECMW_SERIAL
71 integer(kind=kint) :: ierr
72 call mpi_scatter( sbuf, 1, mpi_integer, &
73 & rval, 1, mpi_integer, root, comm, ierr )
74#else
75 rval=sbuf(1)
76#endif
77 end subroutine hecmw_scatter_int_1
78
79 subroutine hecmw_gather_int_1(sval, rbuf, root, comm)
80 use hecmw_util
81 implicit none
82 integer(kind=kint) :: sval !send buffer
83 integer(kind=kint) :: rbuf(*) !receive buffer
84 integer(kind=kint) :: root
85 integer(kind=kint) :: comm
86#ifndef HECMW_SERIAL
87 integer(kind=kint) :: ierr
88 call mpi_gather( sval, 1, mpi_integer, &
89 & rbuf, 1, mpi_integer, root, comm, ierr )
90#else
91 rbuf(1)=sval
92#endif
93 end subroutine hecmw_gather_int_1
94
95 subroutine hecmw_allgather_int_1(sval, rbuf, comm)
96 use hecmw_util
97 implicit none
98 integer(kind=kint) :: sval !send buffer
99 integer(kind=kint) :: rbuf(*) !receive buffer
100 integer(kind=kint) :: comm
101#ifndef HECMW_SERIAL
102 integer(kind=kint) :: ierr
103 call mpi_allgather( sval, 1, mpi_integer, &
104 & rbuf, 1, mpi_integer, comm, ierr )
105#else
106 rbuf(1)=sval
107#endif
108 end subroutine hecmw_allgather_int_1
109
110 subroutine hecmw_scatterv_real(sbuf, scs, disp, &
111 & rbuf, rc, root, comm)
112 use hecmw_util
113 implicit none
114 real(kind=kreal) :: sbuf(*) !send buffer
115 integer(kind=kint) :: scs(*) !send counts
116 integer(kind=kint) :: disp(*) !displacement
117 real(kind=kreal) :: rbuf(*) !receive buffer
118 integer(kind=kint) :: rc !receive counts
119 integer(kind=kint) :: root
120 integer(kind=kint) :: comm
121#ifndef HECMW_SERIAL
122 integer(kind=kint) :: ierr
123 call mpi_scatterv( sbuf, scs, disp, mpi_real8, &
124 & rbuf, rc, mpi_real8, root, comm, ierr )
125#else
126 rbuf(1:rc)=sbuf(1:rc)
127#endif
128 end subroutine hecmw_scatterv_real
129
130 subroutine hecmw_gatherv_real(sbuf, sc, &
131 & rbuf, rcs, disp, root, comm)
132 use hecmw_util
133 implicit none
134 real(kind=kreal) :: sbuf(*) !send buffer
135 integer(kind=kint) :: sc !send counts
136 real(kind=kreal) :: rbuf(*) !receive buffer
137 integer(kind=kint) :: rcs(*) !receive counts
138 integer(kind=kint) :: disp(*) !displacement
139 integer(kind=kint) :: root
140 integer(kind=kint) :: comm
141#ifndef HECMW_SERIAL
142 integer(kind=kint) :: ierr
143 call mpi_gatherv( sbuf, sc, mpi_real8, &
144 & rbuf, rcs, disp, mpi_real8, root, comm, ierr )
145#else
146 rbuf(1:sc)=sbuf(1:sc)
147#endif
148 end subroutine hecmw_gatherv_real
149
150 subroutine hecmw_gatherv_int(sbuf, sc, &
151 & rbuf, rcs, disp, root, comm)
152 use hecmw_util
153 implicit none
154 integer(kind=kint) :: sbuf(*) !send buffer
155 integer(kind=kint) :: sc !send counts
156 integer(kind=kint) :: rbuf(*) !receive buffer
157 integer(kind=kint) :: rcs(*) !receive counts
158 integer(kind=kint) :: disp(*) !displacement
159 integer(kind=kint) :: root
160 integer(kind=kint) :: comm
161#ifndef HECMW_SERIAL
162 integer(kind=kint) :: ierr
163 call mpi_gatherv( sbuf, sc, mpi_integer, &
164 & rbuf, rcs, disp, mpi_integer, root, comm, ierr )
165#else
166 rbuf(1:sc)=sbuf(1:sc)
167#endif
168 end subroutine hecmw_gatherv_int
169
170 subroutine hecmw_allreduce_int_1(sval, rval, op, comm)
171 use hecmw_util
172 implicit none
173 integer(kind=kint) :: sval !send val
174 integer(kind=kint) :: rval !receive val
175 integer(kind=kint):: op, comm, ierr
176#ifndef HECMW_SERIAL
177 call mpi_allreduce(sval, rval, 1, mpi_integer, &
178 & hecmw_operation_hec2mpi(op), comm, ierr)
179#else
180 rval=sval
181#endif
182 end subroutine hecmw_allreduce_int_1
183
184 subroutine hecmw_alltoall_int(sbuf, sc, rbuf, rc, comm)
185 use hecmw_util
186 implicit none
187 integer(kind=kint) :: sbuf(*)
188 integer(kind=kint) :: sc
189 integer(kind=kint) :: rbuf(*)
190 integer(kind=kint) :: rc
191 integer(kind=kint) :: comm
192#ifndef HECMW_SERIAL
193 integer(kind=kint) :: ierr
194 call mpi_alltoall( sbuf, sc, mpi_integer, &
195 & rbuf, rc, mpi_integer, comm, ierr )
196#else
197 rbuf(1:sc)=sbuf(1:sc)
198#endif
199 end subroutine hecmw_alltoall_int
200
201 subroutine hecmw_isend_int(sbuf, sc, dest, &
202 & tag, comm, req)
203 use hecmw_util
204 implicit none
205 integer(kind=kint) :: sbuf(*)
206 integer(kind=kint) :: sc
207 integer(kind=kint) :: dest
208 integer(kind=kint) :: tag
209 integer(kind=kint) :: comm
210 integer(kind=kint) :: req
211#ifndef HECMW_SERIAL
212 integer(kind=kint) :: ierr
213 call mpi_isend(sbuf, sc, mpi_integer, &
214 & dest, tag, comm, req, ierr)
215#endif
216 end subroutine hecmw_isend_int
217
218 subroutine hecmw_isend_r(sbuf, sc, dest, &
219 & tag, comm, req)
220 use hecmw_util
221 implicit none
222 integer(kind=kint) :: sc
223 double precision, dimension(sc) :: sbuf
224 integer(kind=kint) :: dest
225 integer(kind=kint) :: tag
226 integer(kind=kint) :: comm
227 integer(kind=kint) :: req
228#ifndef HECMW_SERIAL
229 integer(kind=kint) :: ierr
230 call mpi_isend(sbuf, sc, mpi_double_precision, &
231 & dest, tag, comm, req, ierr)
232#endif
233 end subroutine hecmw_isend_r
234
235 subroutine hecmw_irecv_int(rbuf, rc, source, &
236 & tag, comm, req)
237 use hecmw_util
238 implicit none
239 integer(kind=kint) :: rbuf(*)
240 integer(kind=kint) :: rc
241 integer(kind=kint) :: source
242 integer(kind=kint) :: tag
243 integer(kind=kint) :: comm
244 integer(kind=kint) :: req
245#ifndef HECMW_SERIAL
246 integer(kind=kint) :: ierr
247 call mpi_irecv(rbuf, rc, mpi_integer, &
248 & source, tag, comm, req, ierr)
249#endif
250 end subroutine hecmw_irecv_int
251
252 subroutine hecmw_irecv_r(rbuf, rc, source, &
253 & tag, comm, req)
254 use hecmw_util
255 implicit none
256 integer(kind=kint) :: rc
257 double precision, dimension(rc) :: rbuf
258 integer(kind=kint) :: source
259 integer(kind=kint) :: tag
260 integer(kind=kint) :: comm
261 integer(kind=kint) :: req
262#ifndef HECMW_SERIAL
263 integer(kind=kint) :: ierr
264 call mpi_irecv(rbuf, rc, mpi_double_precision, &
265 & source, tag, comm, req, ierr)
266#endif
267 end subroutine hecmw_irecv_r
268
269 subroutine hecmw_waitall(cnt, reqs, stats)
270 use hecmw_util
271 implicit none
272 integer(kind=kint) :: cnt
273 integer(kind=kint) :: reqs(*)
274 integer(kind=kint) :: stats(HECMW_STATUS_SIZE,*)
275#ifndef HECMW_SERIAL
276 integer(kind=kint) :: ierr
277 call mpi_waitall(cnt, reqs, stats, ierr)
278#endif
279 end subroutine hecmw_waitall
280
281 subroutine hecmw_recv_int(rbuf, rc, source, &
282 & tag, comm, stat)
283 use hecmw_util
284 implicit none
285 integer(kind=kint) :: rbuf(*)
286 integer(kind=kint) :: rc
287 integer(kind=kint) :: source
288 integer(kind=kint) :: tag
289 integer(kind=kint) :: comm
290 integer(kind=kint) :: stat(HECMW_STATUS_SIZE)
291#ifndef HECMW_SERIAL
292 integer(kind=kint) :: ierr
293 call mpi_recv(rbuf, rc, mpi_integer, &
294 & source, tag, comm, stat, ierr)
295#endif
296 end subroutine hecmw_recv_int
297
298 subroutine hecmw_recv_r(rbuf, rc, source, &
299 & tag, comm, stat)
300 use hecmw_util
301 implicit none
302 integer(kind=kint) :: rc
303 double precision, dimension(rc) :: rbuf
304 integer(kind=kint) :: source
305 integer(kind=kint) :: tag
306 integer(kind=kint) :: comm
307 integer(kind=kint) :: stat(HECMW_STATUS_SIZE)
308#ifndef HECMW_SERIAL
309 integer(kind=kint) :: ierr
310 call mpi_recv(rbuf, rc, mpi_double_precision, &
311 & source, tag, comm, stat, ierr)
312#endif
313 end subroutine hecmw_recv_r
314 !C
315 !C***
316 !C*** hecmw_allREDUCE
317 !C***
318 !C
319 subroutine hecmw_allreduce_dp(val,VALM,n,hec_op,comm )
320 use hecmw_util
321 implicit none
322 integer(kind=kint) :: n, hec_op,op, comm, ierr
323 double precision, dimension(n) :: val
324 double precision, dimension(n) :: VALM
325#ifndef HECMW_SERIAL
326 select case( hec_op )
327 case ( hecmw_sum )
328 op = mpi_sum
329 case ( hecmw_prod )
330 op = mpi_prod
331 case ( hecmw_max )
332 op = mpi_max
333 case ( hecmw_min )
334 op = mpi_min
335 end select
336 call mpi_allreduce(val,valm,n,mpi_double_precision,op,comm,ierr)
337#else
338 integer(kind=kint) :: i
339 do i=1,n
340 valm(i) = val(i)
341 end do
342#endif
343 end subroutine hecmw_allreduce_dp
344
345 subroutine hecmw_allreduce_dp1(s1,s2,hec_op,comm )
346 use hecmw_util
347 implicit none
348 integer(kind=kint) :: hec_op, comm
349 double precision :: s1, s2
350 double precision, dimension(1) :: val
351 double precision, dimension(1) :: VALM
352#ifndef HECMW_SERIAL
353 val(1) = s1
354 valm(1) = s2
355 call hecmw_allreduce_dp(val,valm,1,hec_op,comm )
356 s1 = val(1)
357 s2 = valm(1)
358#else
359 s2 = s1
360#endif
361 end subroutine hecmw_allreduce_dp1
362 !C
363 !C***
364 !C*** hecmw_allREDUCE_R
365 !C***
366 !C
367 subroutine hecmw_allreduce_r (hecMESH, val, n, ntag)
368 use hecmw_util
369 implicit none
370 integer(kind=kint):: n, ntag
371 real(kind=kreal), dimension(n) :: val
372 type (hecmwST_local_mesh) :: hecMESH
373#ifndef HECMW_SERIAL
374 integer(kind=kint):: ierr
375 real(kind=kreal), dimension(:), allocatable :: valm
376
377 allocate (valm(n))
378 valm= 0.d0
379 if (ntag .eq. hecmw_sum) then
380 call mpi_allreduce &
381 & (val, valm, n, mpi_double_precision, mpi_sum, &
382 & hecmesh%MPI_COMM, ierr)
383 endif
384
385 if (ntag .eq. hecmw_max) then
386 call mpi_allreduce &
387 & (val, valm, n, mpi_double_precision, mpi_max, &
388 & hecmesh%MPI_COMM, ierr)
389 endif
390
391 if (ntag .eq. hecmw_min) then
392 call mpi_allreduce &
393 & (val, valm, n, mpi_double_precision, mpi_min, &
394 & hecmesh%MPI_COMM, ierr)
395 endif
396
397 val= valm
398 deallocate (valm)
399#endif
400 end subroutine hecmw_allreduce_r
401
402 subroutine hecmw_allreduce_r1 (hecMESH, s, ntag)
403 use hecmw_util
404 implicit none
405 integer(kind=kint):: ntag
406 real(kind=kreal) :: s
407 type (hecmwST_local_mesh) :: hecMESH
408#ifndef HECMW_SERIAL
409 real(kind=kreal), dimension(1) :: val
410 val(1) = s
411 call hecmw_allreduce_r(hecmesh, val, 1, ntag )
412 s = val(1)
413#endif
414 end subroutine hecmw_allreduce_r1
415
416 !C
417 !C***
418 !C*** hecmw_allREDUCE_I
419 !C***
420 !C
421 subroutine hecmw_allreduce_i(hecMESH, val, n, ntag)
422 use hecmw_util
423 implicit none
424 integer(kind=kint):: n, ntag
425 integer(kind=kint), dimension(n) :: val
426 type (hecmwST_local_mesh) :: hecMESH
427#ifndef HECMW_SERIAL
428 integer(kind=kint):: ierr
429 integer(kind=kint), dimension(:), allocatable :: VALM
430
431 allocate (valm(n))
432 valm= 0
433 if (ntag .eq. hecmw_sum) then
434 call mpi_allreduce &
435 & (val, valm, n, mpi_integer, mpi_sum, &
436 & hecmesh%MPI_COMM, ierr)
437 endif
438
439 if (ntag .eq. hecmw_max) then
440 call mpi_allreduce &
441 & (val, valm, n, mpi_integer, mpi_max, &
442 & hecmesh%MPI_COMM, ierr)
443 endif
444
445 if (ntag .eq. hecmw_min) then
446 call mpi_allreduce &
447 & (val, valm, n, mpi_integer, mpi_min, &
448 & hecmesh%MPI_COMM, ierr)
449 endif
450
451
452 val= valm
453 deallocate (valm)
454#endif
455 end subroutine hecmw_allreduce_i
456
457 subroutine hecmw_allreduce_i1 (hecMESH, s, ntag)
458 use hecmw_util
459 implicit none
460 integer(kind=kint):: ntag, s
461 type (hecmwST_local_mesh) :: hecMESH
462#ifndef HECMW_SERIAL
463 integer(kind=kint), dimension(1) :: val
464
465 val(1) = s
466 call hecmw_allreduce_i(hecmesh, val, 1, ntag )
467 s = val(1)
468#endif
469 end subroutine hecmw_allreduce_i1
470
471 !C
472 !C***
473 !C*** hecmw_bcast_R
474 !C***
475 !C
476 subroutine hecmw_bcast_r (hecMESH, val, n, nbase)
477 use hecmw_util
478 implicit none
479 integer(kind=kint):: n, nbase
480 real(kind=kreal), dimension(n) :: val
481 type (hecmwST_local_mesh) :: hecMESH
482#ifndef HECMW_SERIAL
483 integer(kind=kint):: ierr
484 call mpi_bcast (val, n, mpi_double_precision, nbase, hecmesh%MPI_COMM, ierr)
485#endif
486 end subroutine hecmw_bcast_r
487
488 subroutine hecmw_bcast_r_comm (val, n, nbase, comm)
489 use hecmw_util
490 implicit none
491 integer(kind=kint):: n, nbase
492 real(kind=kreal), dimension(n) :: val
493 integer(kind=kint):: comm
494#ifndef HECMW_SERIAL
495 integer(kind=kint):: ierr
496 call mpi_bcast (val, n, mpi_double_precision, nbase, comm, ierr)
497#endif
498 end subroutine hecmw_bcast_r_comm
499
500 subroutine hecmw_bcast_r1 (hecMESH, s, nbase)
501 use hecmw_util
502 implicit none
503 integer(kind=kint):: nbase, ierr
504 real(kind=kreal) :: s
505 type (hecmwST_local_mesh) :: hecMESH
506#ifndef HECMW_SERIAL
507 real(kind=kreal), dimension(1) :: val
508 val(1)=s
509 call mpi_bcast (val, 1, mpi_double_precision, nbase, hecmesh%MPI_COMM, ierr)
510 s = val(1)
511#endif
512 end subroutine hecmw_bcast_r1
513
514 subroutine hecmw_bcast_r1_comm (s, nbase, comm)
515 use hecmw_util
516 implicit none
517 integer(kind=kint):: nbase
518 real(kind=kreal) :: s
519 integer(kind=kint):: comm
520#ifndef HECMW_SERIAL
521 integer(kind=kint):: ierr
522 real(kind=kreal), dimension(1) :: val
523 val(1)=s
524 call mpi_bcast (val, 1, mpi_double_precision, nbase, comm, ierr)
525 s = val(1)
526#endif
527 end subroutine hecmw_bcast_r1_comm
528 !C
529 !C***
530 !C*** hecmw_bcast_I
531 !C***
532 !C
533 subroutine hecmw_bcast_i (hecMESH, val, n, nbase)
534 use hecmw_util
535 implicit none
536 integer(kind=kint):: n, nbase
537 integer(kind=kint), dimension(n) :: val
538 type (hecmwST_local_mesh) :: hecMESH
539#ifndef HECMW_SERIAL
540 integer(kind=kint):: ierr
541 call mpi_bcast (val, n, mpi_integer, nbase, hecmesh%MPI_COMM, ierr)
542#endif
543 end subroutine hecmw_bcast_i
544
545 subroutine hecmw_bcast_i_comm (val, n, nbase, comm)
546 use hecmw_util
547 implicit none
548 integer(kind=kint):: n, nbase
549 integer(kind=kint), dimension(n) :: val
550 integer(kind=kint):: comm
551#ifndef HECMW_SERIAL
552 integer(kind=kint):: ierr
553 call mpi_bcast (val, n, mpi_integer, nbase, comm, ierr)
554#endif
555 end subroutine hecmw_bcast_i_comm
556
557 subroutine hecmw_bcast_i1 (hecMESH, s, nbase)
558 use hecmw_util
559 implicit none
560 integer(kind=kint):: nbase, s
561 type (hecmwST_local_mesh) :: hecMESH
562#ifndef HECMW_SERIAL
563 integer(kind=kint):: ierr
564 integer(kind=kint), dimension(1) :: val
565 val(1) = s
566 call mpi_bcast (val, 1, mpi_integer, nbase, hecmesh%MPI_COMM, ierr)
567 s = val(1)
568#endif
569 end subroutine hecmw_bcast_i1
570
571 subroutine hecmw_bcast_i1_comm (s, nbase, comm)
572 use hecmw_util
573 implicit none
574 integer(kind=kint):: nbase, s
575 integer(kind=kint):: comm
576#ifndef HECMW_SERIAL
577 integer(kind=kint):: ierr
578 integer(kind=kint), dimension(1) :: val
579 val(1) = s
580 call mpi_bcast (val, 1, mpi_integer, nbase, comm, ierr)
581 s = val(1)
582#endif
583 end subroutine hecmw_bcast_i1_comm
584 !C
585 !C***
586 !C*** hecmw_bcast_C
587 !C***
588 !C
589 subroutine hecmw_bcast_c (hecMESH, val, n, nn, nbase)
590 use hecmw_util
591 implicit none
592 integer(kind=kint):: n, nn, nbase
593 character(len=n) :: val(nn)
594 type (hecmwST_local_mesh) :: hecMESH
595#ifndef HECMW_SERIAL
596 integer(kind=kint):: ierr
597 call mpi_bcast (val, n*nn, mpi_character, nbase, hecmesh%MPI_COMM,&
598 & ierr)
599#endif
600 end subroutine hecmw_bcast_c
601
602 subroutine hecmw_bcast_c_comm (val, n, nn, nbase, comm)
603 use hecmw_util
604 implicit none
605 integer(kind=kint):: n, nn, nbase
606 character(len=n) :: val(nn)
607 integer(kind=kint):: comm
608#ifndef HECMW_SERIAL
609 integer(kind=kint):: ierr
610 call mpi_bcast (val, n*nn, mpi_character, nbase, comm,&
611 & ierr)
612#endif
613 end subroutine hecmw_bcast_c_comm
614
615 !C
616 !C***
617 !C*** hecmw_update_1_R
618 !C***
619 !C
620 !C 1-DOF, REAL
621 !C
622 subroutine hecmw_update_1_r (hecMESH, val, n)
623 use hecmw_util
625
626 implicit none
627 integer(kind=kint):: n
628 real(kind=kreal), dimension(n) :: val
629 type (hecmwST_local_mesh) :: hecMESH
630#ifndef HECMW_SERIAL
631 integer(kind=kint):: ns, nr
632 real(kind=kreal), dimension(:), allocatable :: ws, wr
633
634 if( hecmesh%n_neighbor_pe == 0 ) return
635
636 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
637 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
638
639 allocate (ws(ns), wr(nr))
641 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
642 & hecmesh%import_index, hecmesh%import_item, &
643 & hecmesh%export_index, hecmesh%export_item, &
644 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
645 deallocate (ws, wr)
646#endif
647 end subroutine hecmw_update_1_r
648
649 !C
650 !C***
651 !C*** hecmw_update_2_R
652 !C***
653 !C
654 !C 2-DOF, REAL
655 !C
656 subroutine hecmw_update_2_r (hecMESH, val, n)
657 use hecmw_util
659
660 implicit none
661 integer(kind=kint):: n
662 real(kind=kreal), dimension(2*n) :: val
663 type (hecmwST_local_mesh) :: hecMESH
664#ifndef HECMW_SERIAL
665 integer(kind=kint):: ns, nr
666 real(kind=kreal), dimension(:), allocatable :: ws, wr
667
668 if( hecmesh%n_neighbor_pe == 0 ) return
669
670 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
671 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
672
673 allocate (ws(2*ns), wr(2*nr))
675 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
676 & hecmesh%import_index, hecmesh%import_item, &
677 & hecmesh%export_index, hecmesh%export_item, &
678 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
679 deallocate (ws, wr)
680#endif
681 end subroutine hecmw_update_2_r
682
683 !C
684 !C***
685 !C*** hecmw_update_3_R
686 !C***
687 !C
688 !C 3-DOF, REAL
689 !C
690 subroutine hecmw_update_3_r (hecMESH, val, n)
691 use hecmw_util
693
694 implicit none
695 integer(kind=kint):: n
696 real(kind=kreal), dimension(3*n) :: val
697 type (hecmwST_local_mesh) :: hecMESH
698#ifndef HECMW_SERIAL
699 integer(kind=kint):: ns, nr
700 real(kind=kreal), dimension(:), allocatable :: ws, wr
701
702 if( hecmesh%n_neighbor_pe == 0 ) return
703
704 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
705 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
706
707 allocate (ws(3*ns), wr(3*nr))
709 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
710 & hecmesh%import_index, hecmesh%import_item, &
711 & hecmesh%export_index, hecmesh%export_item, &
712 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
713 deallocate (ws, wr)
714#endif
715 end subroutine hecmw_update_3_r
716
717 !C
718 !C***
719 !C*** hecmw_update_3_R_async
720 !C***
721 !C
722 !C 3-DOF, REAL
723 !C
724 subroutine hecmw_update_3_r_async (hecMESH, val, n, ireq)
725 use hecmw_util
727 implicit none
728 integer(kind=kint):: n, ireq
729 real(kind=kreal), dimension(3*n) :: val
730 type (hecmwST_local_mesh) :: hecMESH
731#ifndef HECMW_SERIAL
732 if( hecmesh%n_neighbor_pe == 0 ) return
733
735 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
736 & hecmesh%import_index, hecmesh%import_item, &
737 & hecmesh%export_index, hecmesh%export_item, &
738 & val , hecmesh%MPI_COMM, hecmesh%my_rank, ireq)
739#endif
740 end subroutine hecmw_update_3_r_async
741
742 !C
743 !C***
744 !C*** hecmw_update_3_R_wait
745 !C***
746 !C
747 !C 3-DOF, REAL
748 !C
749 subroutine hecmw_update_3_r_wait (hecMESH, ireq)
750 use hecmw_util
752 implicit none
753 integer(kind=kint):: ireq
754 type (hecmwST_local_mesh) :: hecMESH
755#ifndef HECMW_SERIAL
756 if( hecmesh%n_neighbor_pe == 0 ) return
757
759#endif
760 end subroutine hecmw_update_3_r_wait
761
762 !C
763 !C***
764 !C*** hecmw_update_4_R
765 !C***
766 !C
767 !C 4-DOF, REAL
768 !C
769 subroutine hecmw_update_4_r (hecMESH, val, n)
770 use hecmw_util
772
773 implicit none
774 integer(kind=kint):: n
775 real(kind=kreal), dimension(4*n) :: val
776 type (hecmwST_local_mesh) :: hecMESH
777#ifndef HECMW_SERIAL
778 integer(kind=kint):: ns, nr
779 real(kind=kreal), dimension(:), allocatable :: ws, wr
780
781 if( hecmesh%n_neighbor_pe == 0 ) return
782
783 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
784 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
785
786 allocate (ws(4*ns), wr(4*nr))
788 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
789 & hecmesh%import_index, hecmesh%import_item, &
790 & hecmesh%export_index, hecmesh%export_item, &
791 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
792 deallocate (ws, wr)
793#endif
794 end subroutine hecmw_update_4_r
795
796 !C
797 !C***
798 !C*** hecmw_update_6_R
799 !C***
800 !C
801 !C 6-DOF, REAL
802 !C
803 subroutine hecmw_update_6_r (hecMESH, val, n)
804 use hecmw_util
806
807 implicit none
808 integer(kind=kint):: n
809 real(kind=kreal), dimension(6*n) :: val
810 type (hecmwST_local_mesh) :: hecMESH
811#ifndef HECMW_SERIAL
812 integer(kind=kint):: ns, nr
813 real(kind=kreal), dimension(:), allocatable :: ws, wr
814
815 if( hecmesh%n_neighbor_pe == 0 ) return
816
817 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
818 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
819
820 allocate (ws(6*ns), wr(6*nr))
822 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
823 & hecmesh%import_index, hecmesh%import_item, &
824 & hecmesh%export_index, hecmesh%export_item, &
825 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
826 deallocate (ws, wr)
827#endif
828 end subroutine hecmw_update_6_r
829
830 !C
831 !C***
832 !C*** hecmw_update_m_R
833 !C***
834 !C
835 !C m-DOF, REAL
836 !C
837 subroutine hecmw_update_m_r (hecMESH, val, n, m)
838 use hecmw_util
840
841 implicit none
842 integer(kind=kint):: n, m
843 real(kind=kreal), dimension(m*n) :: val
844 type (hecmwST_local_mesh) :: hecMESH
845#ifndef HECMW_SERIAL
846 integer(kind=kint):: ns, nr
847 real(kind=kreal), dimension(:), allocatable :: ws, wr
848
849 if( hecmesh%n_neighbor_pe == 0 ) return
850
851 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
852 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
853
854 allocate (ws(m*ns), wr(m*nr))
856 & ( n, m, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
857 & hecmesh%import_index, hecmesh%import_item, &
858 & hecmesh%export_index, hecmesh%export_item, &
859 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
860 deallocate (ws, wr)
861#endif
862 end subroutine hecmw_update_m_r
863
864 !C
865 !C***
866 !C*** hecmw_update_1_I
867 !C***
868 !C
869 !C 1-DOF, INTEGER
870 !C
871 subroutine hecmw_update_1_i (hecMESH, val, n)
872 use hecmw_util
874
875 implicit none
876 integer(kind=kint):: n
877 integer(kind=kint), dimension(n) :: val
878 type (hecmwST_local_mesh) :: hecMESH
879#ifndef HECMW_SERIAL
880 integer(kind=kint):: ns, nr
881 integer(kind=kint), dimension(:), allocatable :: WS, WR
882
883 if( hecmesh%n_neighbor_pe == 0 ) return
884
885 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
886 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
887
888 allocate (ws(ns), wr(nr))
890 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
891 & hecmesh%import_index, hecmesh%import_item, &
892 & hecmesh%export_index, hecmesh%export_item, &
893 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
894 deallocate (ws, wr)
895#endif
896 end subroutine hecmw_update_1_i
897
898 !C
899 !C***
900 !C*** hecmw_update_2_I
901 !C***
902 !C
903 !C 2-DOF, INTEGER
904 !C
905 subroutine hecmw_update_2_i (hecMESH, val, n)
906 use hecmw_util
908
909 implicit none
910 integer(kind=kint):: n
911 integer(kind=kint), dimension(2*n) :: val
912 type (hecmwST_local_mesh) :: hecMESH
913#ifndef HECMW_SERIAL
914 integer(kind=kint):: ns, nr
915 integer(kind=kint), dimension(:), allocatable :: WS, WR
916
917 if( hecmesh%n_neighbor_pe == 0 ) return
918
919 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
920 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
921
922 allocate (ws(2*ns), wr(2*nr))
924 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
925 & hecmesh%import_index, hecmesh%import_item, &
926 & hecmesh%export_index, hecmesh%export_item, &
927 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
928 deallocate (ws, wr)
929#endif
930 end subroutine hecmw_update_2_i
931
932 !C
933 !C***
934 !C*** hecmw_update_3_I
935 !C***
936 !C
937 !C 3-DOF, INTEGER
938 !C
939 subroutine hecmw_update_3_i (hecMESH, val, n)
940 use hecmw_util
942
943 implicit none
944 integer(kind=kint):: n
945 integer(kind=kint), dimension(3*n) :: val
946 type (hecmwST_local_mesh) :: hecMESH
947#ifndef HECMW_SERIAL
948 integer(kind=kint):: ns, nr
949 integer(kind=kint), dimension(:), allocatable :: WS, WR
950
951 if( hecmesh%n_neighbor_pe == 0 ) return
952
953 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
954 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
955
956 allocate (ws(3*ns), wr(3*nr))
958 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
959 & hecmesh%import_index, hecmesh%import_item, &
960 & hecmesh%export_index, hecmesh%export_item, &
961 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
962 deallocate (ws, wr)
963#endif
964 end subroutine hecmw_update_3_i
965
966 !C
967 !C***
968 !C*** hecmw_update_4_I
969 !C***
970 !C
971 !C 4-DOF, INTEGER
972 !C
973 subroutine hecmw_update_4_i (hecMESH, val, n)
974 use hecmw_util
976
977 implicit none
978 integer(kind=kint):: n
979 integer(kind=kint), dimension(4*n) :: val
980 type (hecmwST_local_mesh) :: hecMESH
981#ifndef HECMW_SERIAL
982 integer(kind=kint):: ns, nr
983 integer(kind=kint), dimension(:), allocatable :: WS, WR
984
985 if( hecmesh%n_neighbor_pe == 0 ) return
986
987 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
988 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
989
990 allocate (ws(4*ns), wr(4*nr))
992 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
993 & hecmesh%import_index, hecmesh%import_item, &
994 & hecmesh%export_index, hecmesh%export_item, &
995 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
996 deallocate (ws, wr)
997#endif
998 end subroutine hecmw_update_4_i
999
1000 !C
1001 !C***
1002 !C*** hecmw_update_6_I
1003 !C***
1004 !C
1005 !C 6-DOF, INTEGER
1006 !C
1007 subroutine hecmw_update_6_i (hecMESH, val, n)
1008 use hecmw_util
1010
1011 implicit none
1012 integer(kind=kint):: n
1013 integer(kind=kint), dimension(6*n) :: val
1014 type (hecmwST_local_mesh) :: hecMESH
1015#ifndef HECMW_SERIAL
1016 integer(kind=kint):: ns, nr
1017 integer(kind=kint), dimension(:), allocatable :: WS, WR
1018
1019 if( hecmesh%n_neighbor_pe == 0 ) return
1020
1021 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
1022 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
1023
1024 allocate (ws(6*ns), wr(6*nr))
1026 & ( n, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
1027 & hecmesh%import_index, hecmesh%import_item, &
1028 & hecmesh%export_index, hecmesh%export_item, &
1029 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
1030 deallocate (ws, wr)
1031#endif
1032 end subroutine hecmw_update_6_i
1033
1034 !C
1035 !C***
1036 !C*** hecmw_update_m_I
1037 !C***
1038 !C
1039 !C m-DOF, REAL
1040 !C
1041 subroutine hecmw_update_m_i (hecMESH, val, n, m)
1042 use hecmw_util
1044
1045 implicit none
1046 integer(kind=kint):: n, m
1047 integer(kind=kint), dimension(m*n) :: val
1048 type (hecmwST_local_mesh) :: hecMESH
1049#ifndef HECMW_SERIAL
1050 integer(kind=kint):: ns, nr
1051 integer(kind=kint), dimension(:), allocatable :: WS, WR
1052
1053 if( hecmesh%n_neighbor_pe == 0 ) return
1054
1055 ns = hecmesh%export_index(hecmesh%n_neighbor_pe)
1056 nr = hecmesh%import_index(hecmesh%n_neighbor_pe)
1057
1058 allocate (ws(m*ns), wr(m*nr))
1060 & ( n, m, hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
1061 & hecmesh%import_index, hecmesh%import_item, &
1062 & hecmesh%export_index, hecmesh%export_item, &
1063 & ws, wr, val , hecmesh%MPI_COMM, hecmesh%my_rank)
1064 deallocate (ws, wr)
1065#endif
1066 end subroutine hecmw_update_m_i
1067
1068end module m_hecmw_comm_f
subroutine hecmw_solve_send_recv_11(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_11i(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_22(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_22i(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine, public hecmw_solve_isend_irecv_33_wait(ireq)
subroutine, public hecmw_solve_isend_irecv_33(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, x, solver_comm, my_rank, ireq)
subroutine, public hecmw_solve_send_recv_33(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_33i(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_44(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_44i(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_66(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_66i(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_mm(n, m, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
subroutine hecmw_solve_send_recv_mmi(n, m, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=kint), parameter hecmw_sum
integer(kind=kint), parameter hecmw_prod
integer(kind=kint), parameter hecmw_max
integer(kind=4), parameter kreal
integer(kind=kint), parameter hecmw_min
subroutine hecmw_bcast_c_comm(val, n, nn, nbase, comm)
subroutine hecmw_gatherv_real(sbuf, sc, rbuf, rcs, disp, root, comm)
subroutine hecmw_isend_int(sbuf, sc, dest, tag, comm, req)
subroutine hecmw_bcast_i1(hecmesh, s, nbase)
subroutine hecmw_scatterv_real(sbuf, scs, disp, rbuf, rc, root, comm)
subroutine hecmw_allreduce_i(hecmesh, val, n, ntag)
subroutine hecmw_bcast_i(hecmesh, val, n, nbase)
subroutine hecmw_allreduce_dp(val, valm, n, hec_op, comm)
subroutine hecmw_recv_int(rbuf, rc, source, tag, comm, stat)
subroutine hecmw_bcast_r_comm(val, n, nbase, comm)
subroutine hecmw_update_2_r(hecmesh, val, n)
subroutine hecmw_isend_r(sbuf, sc, dest, tag, comm, req)
subroutine hecmw_update_3_r_async(hecmesh, val, n, ireq)
subroutine hecmw_update_4_i(hecmesh, val, n)
subroutine hecmw_gather_int_1(sval, rbuf, root, comm)
subroutine hecmw_bcast_i_comm(val, n, nbase, comm)
subroutine hecmw_barrier(hecmesh)
subroutine hecmw_allgather_int_1(sval, rbuf, comm)
subroutine hecmw_allreduce_dp1(s1, s2, hec_op, comm)
subroutine hecmw_waitall(cnt, reqs, stats)
subroutine hecmw_update_3_i(hecmesh, val, n)
subroutine hecmw_update_6_r(hecmesh, val, n)
subroutine hecmw_bcast_c(hecmesh, val, n, nn, nbase)
subroutine hecmw_bcast_r1(hecmesh, s, nbase)
subroutine hecmw_irecv_int(rbuf, rc, source, tag, comm, req)
subroutine hecmw_update_m_r(hecmesh, val, n, m)
subroutine hecmw_bcast_r1_comm(s, nbase, comm)
subroutine hecmw_allreduce_r(hecmesh, val, n, ntag)
subroutine hecmw_scatter_int_1(sbuf, rval, root, comm)
subroutine hecmw_irecv_r(rbuf, rc, source, tag, comm, req)
subroutine hecmw_update_4_r(hecmesh, val, n)
subroutine hecmw_update_m_i(hecmesh, val, n, m)
subroutine hecmw_update_3_r(hecmesh, val, n)
subroutine hecmw_allreduce_i1(hecmesh, s, ntag)
subroutine hecmw_recv_r(rbuf, rc, source, tag, comm, stat)
subroutine hecmw_alltoall_int(sbuf, sc, rbuf, rc, comm)
subroutine hecmw_update_1_i(hecmesh, val, n)
subroutine hecmw_update_1_r(hecmesh, val, n)
subroutine hecmw_update_3_r_wait(hecmesh, ireq)
subroutine hecmw_scatterv_dp(sbuf, sc, disp, rbuf, rc, root, comm)
subroutine hecmw_allreduce_int_1(sval, rval, op, comm)
subroutine hecmw_update_2_i(hecmesh, val, n)
subroutine hecmw_bcast_r(hecmesh, val, n, nbase)
subroutine hecmw_gatherv_int(sbuf, sc, rbuf, rcs, disp, root, comm)
subroutine hecmw_allreduce_r1(hecmesh, s, ntag)
subroutine hecmw_update_6_i(hecmesh, val, n)
integer(kind=kint) function hecmw_operation_hec2mpi(operation)
subroutine hecmw_bcast_i1_comm(s, nbase, comm)