KMR
testf.f90
1 ! testf.f90 (2014-02-04)
2 
3 module testfn
4  implicit none
5 
6  integer, target :: foobar
7 
8  type :: tuple2
9  integer :: a0, a1
10  end type tuple2
11 
12 contains
13 
14  function upcase(s) result(ss)
15  use iso_c_binding
16  implicit none
17  character(kind=c_char,len=*), intent(in) :: s
18  character(kind=c_char,len=len(s)) :: ss
19  integer :: n, i
20  character(*), parameter :: UC = "ABCDEFGHIJKLMNOPQRSTUVWXYZ"
21  character(*), parameter :: LC = "abcdefghijklmnopqrstuvwxyz"
22  ss = s
23  do i = 1, len_trim(s)
24  n = index(lc, s(i:i))
25  if (n /= 0) ss(i:i) = uc(n:n)
26  end do
27  end function upcase
28 
29  integer(c_int) function addfivekeysfn(kv, kvi, kvo, p, i) bind(c) result(zz)
30  use iso_c_binding
31  use kmrf
32  implicit none
33  type(kmr_kv_box), value, intent(in) :: kv
34  type(c_ptr), value, intent(in) :: kvi, kvo
35  type(c_ptr), value, intent(in) :: p
36  integer(c_long), value, intent(in) :: i
37  type(kmr_kv_box) :: nkv
38  integer :: j
39  character(kind=c_char,len=1) :: strj
40  character(kind=c_char,len=5), target :: k
41  character(kind=c_char,len=7), target :: v
42  !print *, "mapfn:kvo=", kmr_ptrint(kvo)
43  do j = 0, 4
44  nkv%klen = 5
45  nkv%vlen = 7
46  write(strj, "(I1)") j
47  k = (c_char_"key" // strj // c_null_char)
48  v = (c_char_"value" // strj // c_null_char)
49  nkv%k = kmr_strint(k)
50  nkv%v = kmr_strint(v)
51  zz = kmr_add_kv(kvo, nkv)
52  end do
53  zz = 0
54  end function addfivekeysfn
55 
56  integer(c_int) function replacevaluefn(kv, kvi, kvo, p, i) bind(c) result(zz)
57  use iso_c_binding
58  use kmrf
59  implicit none
60 
61  type(kmr_kv_box), value, intent(in) :: kv
62  type(c_ptr), value, intent(in) :: kvi, kvo
63  type(c_ptr), value, intent(in) :: p
64  integer(c_long), value, intent(in) :: i
65  type(kmr_kv_box) :: nkv
66  character(kind=c_char,len=kv%vlen), target :: v0
67  character(kind=c_char,len=kv%vlen), target :: v1
68  integer(c_int) :: n
69  nkv%klen = kv%klen
70  nkv%vlen = kv%vlen
71  nkv%k = kv%k
72  n = kmr_intstr(kv%v, v0, kv%vlen)
73  v1 = upcase(v0)
74  nkv%v = kmr_strint(v1)
75  zz = kmr_add_kv(kvo, nkv)
76  end function replacevaluefn
77 
78  integer(c_int) function starttaskfn(kv, kvi, kvo, p, i) bind(c) result(zz)
79  use iso_c_binding
80  use kmrf
81  implicit none
82  type(kmr_kv_box), value, intent(in) :: kv
83  type(c_ptr), value, intent(in) :: kvi, kvo
84  type(c_ptr), value, intent(in) :: p
85  integer(c_long), value, intent(in) :: i
86  character(kind=c_char,len=kv%klen), target :: k
87  character(kind=c_char,len=kv%vlen), target :: v
88  integer(c_int) :: n, rank
89  n = kmr_intstr(kv%k, k, kv%klen)
90  n = kmr_intstr(kv%v, v, kv%vlen)
91  rank = kmr_get_rank(kvi)
92  print "(A,A,A,A,A,A,I0)", "Start master/worker task...", &
93  " key=", k(1:kv%klen-1), &
94  " value=", v(1:kv%vlen-1), " rank=", rank
95  zz = kmr_add_kv(kvo, kv)
96  end function starttaskfn
97 
98  integer(c_int) function printpairsfn(kv, n, kvi, kvo, p) bind(c) result(zz)
99  use iso_c_binding
100  use kmrf
101  implicit none
102  type(kmr_kv_box), intent(in) :: kv(*)
103  integer(c_long), value, intent(in) :: n
104  type(c_ptr), value, intent(in) :: kvi, kvo
105  type(c_ptr), value, intent(in) :: p
106  character(kind=c_char,len=5), target :: k
107  character(kind=c_char,len=7), target :: v
108  integer(c_long) :: i
109  integer(c_int) :: nn, rank
110  rank = kmr_get_rank(kvi);
111  print "(A,I0,A,I0)", "Reducing count=", n, " rank=", rank
112  do i = 1, n
113  nn = kmr_intstr(kv(i)%k, k, kv(i)%klen)
114  nn = kmr_intstr(kv(i)%v, v, kv(i)%vlen)
115  print "(A,I0,A,A,A,A)", "index=", i, " key=", k(1:kv(i)%klen-1), &
116  " value=", v(1:kv(i)%vlen-1)
117  end do
118  zz = 0
119  end function printpairsfn
120 
121  integer(c_int) function addstructfn(kv, kvi, kvo, p, i) bind(c) result(zz)
122  use iso_c_binding
123  use kmrf
124  implicit none
125  type(kmr_kv_box), value, intent(in) :: kv
126  type(c_ptr), value, intent(in) :: kvi, kvo
127  type(c_ptr), value, intent(in) :: p
128  integer(c_long), value, intent(in) :: i
129  type(kmr_kv_box) :: nkv
130  character(kind=c_char,len=3) :: stra0
131  character(kind=c_char,len=3) :: stra1
132  character(kind=c_char,len=7), target :: k
133  character(kind=c_char,len=9), target :: v
134  type(tuple2), pointer :: ptr
135  call c_f_pointer(p, ptr)
136  write(stra0, "(I3)") ptr%a0
137  write(stra1, "(I3)") ptr%a1
138  k = (c_char_"key" // stra0 // c_null_char)
139  v = (c_char_"value" // stra1 // c_null_char)
140  nkv%k = kmr_strint(k)
141  nkv%v = kmr_strint(v)
142  nkv%klen = 7
143  nkv%vlen = 9
144  zz = kmr_add_kv(kvo, nkv)
145  end function addstructfn
146 
147  integer(c_int) function addcommandsfn(kv, kvi, kvo, p, ii) bind(c) result(zz)
148  use iso_c_binding
149  use kmrf
150  implicit none
151  type(kmr_kv_box), value, intent(in) :: kv
152  type(c_ptr), value, intent(in) :: kvi, kvo
153  type(c_ptr), value, intent(in) :: p
154  integer(c_long), value, intent(in) :: ii
155  type(kmr_kv_box) :: nkv
156  character(kind=c_char,len=4), target :: k
157  character(kind=c_char,len=100), target :: v
158  character(kind=c_char,len=1), target :: reply
159  integer :: i
160  if (c_associated(p)) then
161  reply = c_char_"1"
162  else
163  reply = c_char_"0"
164  end if
165  nkv%klen = 4
166  nkv%vlen = (11 + 8 + 2 + 3 + 3 + 3)
167  k = (c_char_"key" // c_null_char)
168  v = (c_char_"maxprocs=2" // c_char_" " &
169  // c_char_"./a.out" // c_char_" " &
170  // reply // c_char_" " &
171  // c_char_"a0" // c_char_" " &
172  // c_char_"a1" // c_char_" " &
173  // c_char_"a2" // c_null_char)
174  nkv%k = kmr_strint(k)
175  nkv%v = kmr_strint(v)
176  do i = 1, 4
177  zz = kmr_add_kv(kvo, nkv)
178  end do
179  zz = 0
180  end function addcommandsfn
181 
182  integer(c_int) function waitprocfn(kv, kvi, kvo, p, ii) bind(c) result(zz)
183  use iso_c_binding
184  use kmrf
185  implicit none
186  include "mpif.h"
187  type(kmr_kv_box), value, intent(in) :: kv
188  type(c_ptr), value, intent(in) :: kvi, kvo
189  type(c_ptr), value, intent(in) :: p
190  integer(c_long), value, intent(in) :: ii
191  type(c_ptr) :: mr
192  integer :: ic
193  integer :: ierr
194  mr = kmr_get_context_of_kvs(kvi)
195  ierr = kmr_get_spawner_communicator(mr, ii, ic)
196  print "(A,(I0),A)", "(waitprocfn sleeping(3)... icomm=", ic, ")"
197  !call mpi_comm_free(ic, ierr)
198  call sleep(3)
199  zz = 0
200  end function waitprocfn
201 
202 end module testfn
203 
204 program main
205  use iso_c_binding
206  use kmrf
207  use testfn
208  implicit none
209  include "mpif.h"
210 
211  integer :: argc, ierr, rank, nprocs, thlv
212  character(len=128) argi
213  type(c_ptr) :: mr
214  type(c_ptr) :: kvs0, kvs1, kvs2, kvs3, kvs4, kvs5
215  type(c_ptr) :: kvs10
216  type(c_ptr) :: kvs20, kvs21, kvs22, kvs23, kvs24, kvs25
217  integer(c_long) :: cnt0, cnt1, cnt2, cnt3, cnt5
218  integer(c_int) :: opt
219  character(len=128) keepstack
220  integer :: parent, sz, peernprocs
221  logical :: maybespawned, needreply
222  integer :: iargc
223  character(len=MPI_MAX_PROCESSOR_NAME) :: name
224  integer :: namelen
225  type(tuple2), target :: tuple
226  integer :: sig
227 
228  ! Disable backtracing (which unwinds stack) on ILL/ABRT/SEGV
229 
230  call getenv("KEEPSTACK", keepstack)
231  if (keepstack(1:2) /= "") then
232  call signal(4, 1, sig)
233  call signal(6, 1, sig)
234  call signal(11, 1, sig)
235  end if
236 
237  ! PARENT BE WITHOUT ARGUMENTS; CHILD BE WITH ARGUMENTS.
238 
239  argc = iargc()
240  if (argc >= 1) then
241  maybespawned = .true.
242  else
243  maybespawned = .false.
244  end if
245 
246  call mpi_init_thread(mpi_thread_serialized, thlv, ierr);
247  call mpi_comm_rank(mpi_comm_world, rank, ierr)
248  call mpi_comm_size(mpi_comm_world, nprocs, ierr)
249  call mpi_get_processor_name(name, namelen, ierr)
250 
251  ! CHILD PART:
252 
253  if (maybespawned) then
254  call mpi_comm_get_parent(parent, ierr)
255  if (parent == mpi_comm_null) then
256  print *, "NO PARENTS"
257  call mpi_abort(1);
258  end if
259  call mpi_comm_remote_size(parent, peernprocs, ierr)
260  call kmr_assert(peernprocs == 1, "peernprocs == 1")
261  print *, "Spawned process runs (", trim(name), ")..."
262  call getarg(1, argi)
263  if (argi(1:1) == "0") then
264  needreply = .false.
265  else
266  needreply = .true.
267  end if
268  call sleep(1)
269  if (needreply) then
270  sz = 0
271  call mpi_send(sz, 0, mpi_byte, 0, 500, parent, ierr)
272  print *, "(spawned process send reply (", trim(name), "))"
273  end if
274  call mpi_comm_free(parent, ierr)
275  call mpi_finalize(ierr)
276  print *, "Spawned process runs (", trim(name), ") DONE"
277  call exit(0)
278  end if
279 
280  ! NORMAL (PARENT) PART:
281 
282  ierr = kmr_init()
283 
284  mr = kmr_create_context(mpi_comm_world, mpi_info_null)
285 
286  kvs0 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
287  ierr = kmr_map_on_rank_zero(kvs0, c_null_ptr, 0, addfivekeysfn)
288 
289  ierr = kmr_get_element_count(kvs0, cnt0)
290  !print *, "cnt0=", cnt0
291  !ierr = kmr_dump_kvs(kvs0, 0)
292  call kmr_assert(cnt0 == 5, "cnt0 == 5")
293 
294  kvs1 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
295  ierr = kmr_replicate(kvs0, kvs1, 0);
296 
297  ierr = kmr_get_element_count(kvs1, cnt1)
298  !print *, "cnt1=", cnt1
299  !ierr = kmr_dump_kvs(kvs1, 0)
300  call kmr_assert(cnt1 == (5 * nprocs), "cnt1 == (5 * nprocs)")
301 
302  kvs2 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
303  ierr = kmr_map(kvs1, kvs2, c_null_ptr, 0, replacevaluefn)
304 
305  ierr = kmr_get_element_count(kvs2, cnt2)
306  if (rank == 0) print "(A,I0)", "cnt2=", cnt2
307  !ierr = kmr_dump_kvs(kvs2, 0)
308  call kmr_assert(cnt2 == (5 * nprocs), "cnt2 == (5 * nprocs)")
309 
310  kvs3 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
311  opt = kmr_rank_zero
312  ierr = kmr_replicate(kvs2, kvs3, opt);
313  !ierr = kmr_dump_kvs(kvs3, 0)
314 
315  ierr = kmr_local_element_count(kvs3, cnt3)
316  if (rank == 0) then
317  call kmr_assert(cnt3 == (5 * nprocs), "cnt3 == (5 * nprocs)")
318  else
319  call kmr_assert(cnt3 == 0, "cnt3 == 0")
320  end if
321 
322  kvs4 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
323  opt = kmr_nothreading
324  ierr = 7
325  do while (ierr /= 0)
326  ierr = kmr_map_ms(kvs3, kvs4, c_null_ptr, opt, starttaskfn);
327  end do
328 
329  kvs5 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
330  ierr = kmr_shuffle(kvs4, kvs5, 0);
331 
332  !ierr = kmr_dump_kvs(kvs5, 0)
333  ierr = kmr_get_element_count(kvs5, cnt5)
334  if (rank == 0) print "(A,I0)", "cnt5=", cnt5
335  call kmr_assert(cnt5 == (5 * nprocs), "cnt5 == (5 * nprocs)")
336 
337  call sleep(1)
338 
339  ierr = kmr_reduce(kvs5, c_null_ptr, c_null_ptr, 0, printpairsfn);
340 
341  call sleep(1)
342 
343  ! PASS STRUCTURE POINTERS
344 
345  tuple%a0 = 333
346  tuple%a1 = 555
347 
348  kvs10 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
349  ierr = kmr_map_on_rank_zero(kvs10, c_loc(tuple), 0, addstructfn)
350  ierr = kmr_dump_kvs(kvs10, 0)
351  ierr = kmr_free_kvs(kvs10)
352 
353  call sleep(5)
354  if (rank == 0) print *, "Run spawn test with at least 4 dynamic processes"
355 
356  ! SPAWN WAITING IN MAP-FN
357 
358  call sleep(5)
359  if (rank == 0) print *, "SPAWN WAITING IN MAP-FN"
360 
361  kvs20 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
362  ierr = kmr_map_on_rank_zero(kvs20, c_null_ptr, 0, addcommandsfn)
363 
364  kvs21 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
365  opt = ior(kmr_one_by_one, kmr_separator_space)
366  ierr = kmr_map_via_spawn(kvs20, kvs21, c_loc(foobar), &
367  mpi_info_null, opt, waitprocfn);
368  ierr = kmr_free_kvs(kvs21)
369 
370  call sleep(5)
371  call mpi_barrier(mpi_comm_world, ierr);
372 
373  ! SPAWN WAITING WITH REPLY_EACH
374 
375  call sleep(5)
376  if (rank == 0) print *, "SPAWN WAITING WITH REPLY_EACH"
377 
378  kvs22 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
379  ierr = kmr_map_on_rank_zero(kvs22, mr, 0, addcommandsfn)
380 
381  kvs23 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
382  opt = ior(kmr_reply_each, kmr_separator_space)
383  ierr = kmr_map_via_spawn(kvs22, kvs23, c_null_ptr, &
384  mpi_info_null, opt, kmr_nullmapfn);
385  ierr = kmr_free_kvs(kvs23)
386 
387  call sleep(5)
388  call mpi_barrier(mpi_comm_world, ierr);
389 
390  ! SPAWN NO WAITING
391 
392  if (rank == 0) print *, "SPAWN NO WAITING"
393  if (rank == 0) print *, "THIS PART FAILS WHEN LESS THAN 8 DYNAMIC PROCESSES"
394 
395  kvs24 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
396  ierr = kmr_map_on_rank_zero(kvs24, c_null_ptr, 0, addcommandsfn)
397 
398  kvs25 = kmr_create_kvs(mr, kmr_kv_opaque, kmr_kv_opaque)
399  opt = kmr_separator_space
400  ierr = kmr_map_via_spawn(kvs24, kvs25, c_null_ptr, &
401  mpi_info_null, opt, kmr_nullmapfn);
402  ierr = kmr_free_kvs(kvs25)
403 
404  call sleep(5)
405  call mpi_barrier(mpi_comm_world, ierr);
406 
407  if (rank == 0) print *, "TEST DONE"
408 
409  ierr = kmr_free_context(mr)
410  ierr = kmr_fin()
411 
412  call mpi_finalize(ierr)
413  if (rank == 0) print *, "OK"
414 end program main
#define kmr_reduce(KVI, KVO, ARG, OPT, R)
Reduces key-value pairs.
Definition: kmr.h:88
Converts in reverse of kmr_strint().
Definition: kmrf.F90:186
#define kmr_create_kvs(MR, KF, VF)
Makes a new key-value stream (of type KMR_KVS) with the specified field datatypes.
Definition: kmr.h:71
(See kmr_get_context_of_kvs() in C).
Definition: kmrf.F90:304
MPI_Comm * kmr_get_spawner_communicator(KMR *mr, long index)
Obtains (a reference to) a parent inter-communicator of a spawned process.
Definition: kmrmapms.c:1916
int kmr_shuffle(KMR_KVS *kvi, KMR_KVS *kvo, struct kmr_option opt)
Shuffles key-value pairs to the appropriate destination ranks.
Definition: kmrbase.c:2094
(See kmr_free_context() in C).
Definition: kmrf.F90:293
(See kmr_add_kv() in C).
Definition: kmrf.F90:339
(See kmr_fin() in C).
Definition: kmrf.F90:273
Returns the element count local on each node.
Definition: kmrf.F90:250
Gets MPI rank from a key-value stream.
Definition: kmrf.F90:210
#define kmr_map(KVI, KVO, ARG, OPT, M)
Maps simply.
Definition: kmr.h:82
(See kmr_get_element_count() in C).
Definition: kmrf.F90:499
#define kmr_init()
Sets up the environment.
Definition: kmr.h:794
(See kmr_free_kvs() in C).
Definition: kmrf.F90:329
int kmr_map_via_spawn(KMR_KVS *kvi, KMR_KVS *kvo, void *arg, MPI_Info info, struct kmr_spawn_option opt, kmr_mapfn_t mapfn)
Maps on processes started by MPI_Comm_spawn().
Definition: kmrmapms.c:1992
int kmr_map_ms(KMR_KVS *kvi, KMR_KVS *kvo, void *arg, struct kmr_option opt, kmr_mapfn_t m)
Maps in master-worker mode.
Definition: kmrmapms.c:344
int kmr_replicate(KMR_KVS *kvi, KMR_KVS *kvo, struct kmr_option opt)
Replicates key-value pairs to be visible on all ranks, that is, it has the effect of bcast or all-gat...
Definition: kmrbase.c:2240
Converts a character array to a (pointer value) integer for key/value (it is casting in C)...
Definition: kmrf.F90:174
(See kmr_dump_kvs() in C).
Definition: kmrf.F90:487
int kmr_map_on_rank_zero(KMR_KVS *kvo, void *arg, struct kmr_option opt, kmr_mapfn_t m)
Maps on rank0 only.
Definition: kmrbase.c:1514
KMR * kmr_create_context(const MPI_Comm comm, const MPI_Info conf, const char *name)
Makes a new KMR context (a context has type KMR).
Definition: kmrbase.c:168