1#if HAVE_CONFIG_H
2#   include "config.fh"
3#endif
4c $Id: test.F,v 1.64.2.11 2007-04-06 22:37:35 d3g293 Exp $
5c vector boxes lack arithmetic precision
6#if defined(FUJITSU)
7# define THRESH 1d-12
8# define THRESHF 1e-5
9#else
10# define THRESH 1d-13
11# define THRESHF 2e-5
12#endif
13
14#define MISMATCH(x,y) abs(x-y)/max(1d0,abs(x)).gt.THRESH
15#define MISMATCHF(x,y) abs(x-y)/max(1.0,abs(x)).gt.THRESHF
16
17c#define NEW_API
18c#define MIRROR
19#define GA3
20#define NGA_GATSCAT
21c#define BLOCK_CYCLIC
22c#define USE_SCALAPACK_DISTR
23c#define USE_RESTRICTED
24
25#ifdef USE_RESTRICTED
26#  define NEW_API
27#endif
28
29#define MEM_INC 1000
30
31#ifdef BLOCK_CYCLIC
32#  define NEW_API
33#  undef MIRROR
34#else
35#  undef USE_SCALAPAC_DISTR
36#endif
37
38      program main
39      implicit none
40#include "mafdecls.fh"
41#include "global.fh"
42#include "testutil.fh"
43      integer heap, stack, fudge, ma_heap, me, nproc, map(4096), block
44      integer g_s, ndim, dim1, i
45      logical status
46      parameter (heap=200*200*4, fudge=100, stack=200*200)
47c
48c***  Intitialize a message passing library
49c
50#include "mp3.fh"
51c
52c***  Initialize GA
53c
54c     There are 2 choices: ga_initialize or ga_initialize_ltd.
55c     In the first case, there is no explicit limit on memory usage.
56c     In the second, user can set limit (per processor) in bytes.
57c
58      call ga_initialize()
59      nproc = ga_nnodes()
60      me = ga_nodeid()
61c     we can also use GA_set_memory_limit BEFORE first ga_create call
62c
63      ma_heap = heap/nproc + fudge
64      ma_heap = 2*ma_heap
65#ifdef USE_RESTRICTED
66      ma_heap = 2*ma_heap
67#endif
68      call GA_set_memory_limit(util_mdtob(ma_heap))
69c
70      if(ga_nodeid().eq.0)then
71#ifdef MIRROR
72         print *,' Performing tests on Mirrored Arrays'
73#endif
74         print *,' GA initialized '
75         call ffflush(6)
76      endif
77c
78c***  Initialize the MA package
79c     MA must be initialized before any global array is allocated
80c
81      status = ma_init(MT_DCPL, stack, ma_heap)
82      if (.not. status) call ga_error('ma_init failed',-1)
83c
84c     Uncomment the below line to register external memory allocator
85c     for dynamic arrays inside GA routines.
86c      call register_ext_memory()
87c
88      if(me.eq.(nproc-1))then
89        print *, 'using ', nproc,' process(es) ', ga_cluster_nnodes(),
90     $           ' cluster nodes'
91        print *,'process ', me, ' is on node ',ga_cluster_nodeid(),
92     $          ' with ',  ga_cluster_nprocs(-1), ' processes'
93        call ffflush(6)
94      endif
95c
96c   create array to force staggering of memory and uneven distribution
97c   of pointers
98c
99      dim1 = MEM_INC
100      map(1) = 1
101      do i = 2, nproc
102        map(i) = MEM_INC*(i-1)+1
103        dim1 = dim1 + MEM_INC*i
104      end do
105      g_s = ga_create_handle()
106      ndim = 1
107      call ga_set_data(g_s,ndim,dim1,MT_INT)
108      call ga_set_array_name(g_s,'s')
109      call ga_set_irreg_distr(g_s,map,nproc)
110
111c
112c***  Check support for single precision complex arrays
113c
114      if (me.eq.0) then
115         write(6,*)
116         write(6,*) ' CHECKING SINGLE COMPLEX  '
117         write(6,*)
118         call ffflush(6)
119      endif
120
121      call check_complex_float()
122#if 1
123c
124c***  Check support for double precision complex arrays
125c
126      if (me.eq.0) then
127         write(6,*)
128         write(6,*) ' CHECKING DOUBLE COMPLEX  '
129         write(6,*)
130         call ffflush(6)
131      endif
132
133      call check_complex()
134#endif
135      if(me.eq.0) call ga_print_stats()
136      if(me.eq.0) print *,' '
137      if(me.eq.0) print *,'All tests successful '
138      status = ga_destroy(g_s)
139c
140c***  Tidy up the GA package
141c
142      call ga_terminate()
143c
144c***  Tidy up after message-passing library
145c
146      call MP_FINALIZE()
147c
148      end
149
150
151
152
153      subroutine check_complex_float()
154      implicit none
155#include "mafdecls.fh"
156#include "global.fh"
157#include "testutil.fh"
158c
159      integer n,m
160      parameter (n = 60)
161      parameter (m = 2*n)
162      real areal(n,n),aimg(n,n)
163      real breal(n,n),bimg(n,n)
164#ifdef MIRROR
165      integer ndim, dims(2), chunk(2), p_mirror
166#else
167#  ifdef NEW_API
168      integer ndim, dims(2), chunk(2), p_mirror
169#  endif
170#endif
171      integer iv(m), jv(m)
172      logical status
173      integer g_a, g_b
174      integer iran, i,j, loop,nloop,maxloop, ilo, ihi, jlo, jhi, itmp
175      integer nproc, me, int, ij, inc, ii, jj, nnodes
176      parameter (maxloop = 100)
177      integer maxproc
178      parameter (maxproc = 4096)
179      double precision crap, real
180      double precision nwords
181      complex   x, sum1, sum2, factor
182      integer lprocs, inode, iproc, lproc
183      integer lo(2), hi(2)
184#ifdef USE_RESTRICTED
185      integer num_rstrctd
186      integer rstrctd_list(maxproc/2)
187#endif
188      intrinsic int
189      iran(i) = int(drand(0)*real(i)) + 1
190c
191      nproc = ga_nnodes()
192      me    = ga_nodeid()
193      inode = ga_cluster_nodeid()
194      lprocs = ga_cluster_nprocs(inode)
195      nnodes = ga_cluster_nnodes()
196      iproc = mod(me,lprocs)
197      nloop = Min(maxloop,n)
198#ifdef USE_RESTRICTED
199      num_rstrctd = nproc/2
200      if (num_rstrctd.eq.0) num_rstrctd = 1
201      do i = 1, num_rstrctd
202        rstrctd_list(i) = (num_rstrctd/2) + i-1
203      end do
204#endif
205c
206c     a() is a local copy of what the global array should start as
207c
208      do j = 1, n
209         do i = 1, n
210#ifndef MIRROR
211            areal(i,j) = real(i-1)
212            aimg(i,j)  = real((j-1)*n)
213#else
214            areal(i,j) = real(inode) + real(i-1)
215            aimg(i,j)  = 0.0d00 + real((j-1)*n)
216#endif
217            breal(i,j) = -1d0
218            bimg(i,j)  = 1d0
219         enddo
220      enddo
221c
222c     Create a global array
223c
224c     print *,ga_nodeid(), ' creating array'
225      call ffflush(6)
226c     call setdbg(1)
227#ifdef NEW_API
228      ndim = 2
229      dims(1) = n
230      dims(2) = n
231      g_a = ga_create_handle()
232      call ga_set_data(g_a,ndim,dims,MT_SCPL)
233      call ga_set_array_name(g_a,'a')
234#ifdef USE_RESTRICTED
235      call ga_set_restricted(g_a, rstrctd_list, num_rstrctd)
236#endif
237#  ifdef MIRROR
238      p_mirror = ga_pgroup_get_mirror()
239      call ga_set_pgroup(g_a,p_mirror)
240#  endif
241      status = ga_allocate(g_a)
242#else
243#  ifndef MIRROR
244      status = ga_create(MT_SCPL, n, n, 'a', 0, 0, g_a)
245#  else
246      ndim = 2
247      dims(1) = n
248      dims(2) = n
249      chunk(1) = 0
250      chunk(2) = 0
251      p_mirror = ga_pgroup_get_mirror()
252      status = nga_create_config(MT_SCPL, ndim, dims, 'a', chunk,
253     +                           p_mirror, g_a)
254#  endif
255#endif
256      if (.not. status) then
257         write(6,*) ' ga_create failed'
258         call ga_error('... exiting ',0)
259      endif
260#ifdef NEW_API
261      g_b = ga_create_handle()
262      call ga_set_data(g_b,ndim,dims,MT_SCPL)
263      call ga_set_array_name(g_b,'b')
264#  ifdef MIRROR
265      call ga_set_pgroup(g_b,p_mirror)
266#  endif
267      if (.not.ga_allocate(g_b)) then
268#else
269#  ifndef MIRROR
270      if (.not. ga_create(MT_SCPL, n, n, 'b', 0, 0, g_b)) then
271#  else
272      if (.not. nga_create_config(MT_SCPL, ndim, dims, 'b', chunk,
273     _                            p_mirror, g_b)) then
274#  endif
275#endif
276         call ga_error('ga_create failed for second array ',0)
277      endif
278
279#ifndef MIRROR
280      call ga_distribution(g_a,me,ilo, ihi, jlo, jhi)
281#else
282      lproc = me - ga_cluster_procid(inode,0)
283      call ga_distribution(g_a, lproc, ilo, ihi, jlo, jhi)
284#endif
285      call ga_sync()
286c
287c     Zero the array
288c
289      if (me .eq. 0) then
290         write(6,21)
291 21      format(/'> Checking zero ... ')
292         call ffflush(6)
293      endif
294      call ga_zero(g_a)
295c      call ga_print(g_a)
296c
297c     Check that it is indeed zero
298c
299      lo(1) = 1
300      lo(2) = 1
301      hi(1) = n
302      hi(2) = n
303      call nga_get_field(g_a, lo, hi, 0, 4, breal, n)
304      call nga_get_field(g_a, lo, hi, 4, 4, bimg,  n)
305      call ga_sync()
306      do i = 1, n
307         do j = 1, n
308            if(breal(i,j).ne.0d0) then
309               write(6,*) me,' real zero ', i, j, breal(i,j)
310               call ffflush(6)
311c              call ga_error('... exiting ',0)
312            endif
313            if(bimg(i,j).ne.0d0) then
314               write(6,*) me,' img zero ', i, j, bimg(i,j)
315               call ffflush(6)
316c              call ga_error('... exiting ',0)
317            endif
318         enddo
319      enddo
320      if (me.eq.0) then
321         write(6,*)
322         write(6,*) ' ga_zero is OK'
323         write(6,*)
324      endif
325      call ga_sync()
326c
327c     Each node fills in disjoint sections of the array
328c
329      if (me .eq. 0) then
330         write(6,2)
331 2       format(/'> Checking disjoint put ... ')
332         call ffflush(6)
333      endif
334      call ga_sync()
335c
336      inc = (n-1)/4 + 1
337      ij = 0
338      do j = 1, n, inc
339         do i = 1, n, inc
340#ifndef MIRROR
341            if (mod(ij,nproc) .eq. me) then
342#else
343            if (mod(ij,lprocs) .eq. iproc) then
344#endif
345               ilo = i
346               ihi = min(i+inc, n)
347               jlo = j
348               jhi = min(j+inc, n)
349c               call ga_put(g_a, ilo, ihi, jlo, jhi, a(ilo, jlo), n)
350               lo(1) = ilo
351               lo(2) = jlo
352               hi(1) = ihi
353               hi(2) = jhi
354               call nga_put_field(g_a, lo, hi, 0, 4, areal(ilo, jlo), n)
355               call nga_put_field(g_a, lo, hi, 4, 4, aimg(ilo, jlo), n)
356            endif
357            ij = ij + 1
358         enddo
359      enddo
360      call ga_sync()
361c
362c     All nodes check all of a
363c
364c      call ga_get(g_a, 1, n, 1, n, b, n)
365      lo(1) = 1
366      lo(2) = 1
367      hi(1) = n
368      hi(2) = n
369      call nga_get_field(g_a, lo, hi, 0, 4, breal, n)
370      call nga_get_field(g_a, lo, hi, 4, 4, bimg, n)
371c
372      do i = 1, n
373         do j = 1, n
374            if (breal(i,j) .ne. areal(i,j)) then
375               write(6,*) ' put(real) ', me, i, j, areal(i,j),breal(i,j)
376               call ffflush(6)
377c               call ga_error('... exiting ',0)
378            endif
379            if (bimg(i,j) .ne. aimg(i,j)) then
380               write(6,*) ' put(img) ', me, i, j, aimg(i,j),bimg(i,j)
381               call ffflush(6)
382c               call ga_error('... exiting ',0)
383            endif
384         enddo
385      enddo
386      if (me.eq.0) then
387         write(6,*)
388         write(6,*) ' ga_put is OK'
389         write(6,*)
390      endif
391      call ga_sync()
392c
393c     Now check nloop random gets from each node
394c
395      if (me .eq. 0) then
396         write(6,5) nloop
397 5       format(/'> Checking random get (',i5,' calls)...')
398         call ffflush(6)
399      endif
400      call ga_sync()
401c
402      nwords = 0
403c
404      crap = drand(ga_nodeid()*51 +1) !Different seed for each proc
405      do loop = 1, nloop
406         ilo = iran(loop)
407         ihi = iran(loop)
408         if (ihi.lt. ilo) then
409            itmp = ihi
410            ihi = ilo
411            ilo = itmp
412         endif
413         jlo = iran(loop)
414         jhi = iran(loop)
415         if (jhi.lt. jlo) then
416            itmp = jhi
417            jhi = jlo
418            jlo = itmp
419         endif
420c
421         nwords = nwords + (ihi-ilo+1)*(jhi-jlo+1)
422c
423c         call ga_get(g_a, ilo, ihi, jlo, jhi, b(ilo, jlo), n)
424         lo(1) = ilo
425         lo(2) = jlo
426         hi(1) = ihi
427         hi(2) = jhi
428         call nga_get_field(g_a, lo, hi, 0, 4, breal(ilo, jlo), n)
429         call nga_get_field(g_a, lo, hi, 4, 4, bimg(ilo, jlo), n)
430         if (me .eq. 0 .and. mod(loop-1, max(1,nloop/20)).eq.0) then
431            write(6,1) loop, me, ilo, ihi, jlo, jhi, nwords
432 1          format(' call ',i5, ' node ',i2,' checking get ',4i4,
433     $           ' total ',d9.2)
434            call ffflush(6)
435         endif
436         do j = jlo, jhi
437            do i = ilo, ihi
438               if (breal(i,j) .ne. areal(i,j)) then
439                  write(6,*)'error:', i, j, breal(i,j), areal(i,j)
440                  call ga_error('... exiting ',0)
441               endif
442               if (bimg(i,j) .ne. aimg(i,j)) then
443                  write(6,*)'error:', i, j, bimg(i,j), aimg(i,j)
444                  call ga_error('... exiting ',0)
445               endif
446            enddo
447         enddo
448c
449      enddo
450      if (me.eq.0) then
451         write(6,*)
452         write(6,*) ' ga_get is OK'
453         write(6,*)
454         call ffflush(6)
455      endif
456      call ga_sync()
457c
458c     Check the ga_copy function
459c
460      if (me .eq. 0) then
461         write(6,*)
462         write(6,*)'> Checking copy'
463         write(6,*)
464         call ffflush(6)
465      endif
466      call ga_sync()
467#ifndef MIRROR
468      if(me.eq.0) then
469c         call ga_put(g_a, 1, n, 1, n, a, n)
470         lo(1) = 1
471         lo(2) = 1
472         hi(1) = n
473         hi(2) = n
474         call nga_put_field(g_a, lo, hi, 0, 4, areal, n)
475         call nga_put_field(g_a, lo, hi, 4, 4, aimg, n)
476      endif
477#else
478      if(iproc.eq.0) then
479c         call ga_put(g_a, 1, n, 1, n, a, n)
480         lo(1) = 1
481         lo(2) = 1
482         hi(1) = n
483         hi(2) = n
484         call nga_put_field(g_a, lo, hi, 0, 4, areal, n)
485         call nga_put_field(g_a, lo, hi, 4, 4, aimg, n)
486      endif
487#endif
488      call ga_copy(g_a, g_b)
489c      call ga_get(g_b, 1, n, 1, n, b, n)
490      lo(1) = 1
491      lo(2) = 1
492      hi(1) = n
493      hi(2) = n
494      call nga_get_field(g_b, lo, hi, 0, 4, breal, n)
495      call nga_get_field(g_b, lo, hi, 4, 4, bimg, n)
496      do j = 1, n
497         do i = 1, n
498            if (breal(i,j) .ne. areal(i,j)) then
499               write(6,*) ' copy ', me, i, j, areal(i,j), breal(i,j)
500               call ga_error('... exiting ',0)
501            endif
502            if (bimg(i,j) .ne. aimg(i,j)) then
503               write(6,*) ' copy ', me, i, j, aimg(i,j), bimg(i,j)
504               call ga_error('... exiting ',0)
505            endif
506         enddo
507      enddo
508      if (me.eq.0) then
509         write(6,*)
510         write(6,*) ' copy is OK '
511         write(6,*)
512      endif
513c
514c     Delete the global arrays
515c
516
517      status = ga_destroy(g_b)
518      status = ga_destroy(g_a)
519c
520      end
521c-----------------------------------------------------------------
522
523      subroutine check_complex()
524      implicit none
525#include "mafdecls.fh"
526#include "global.fh"
527#include "testutil.fh"
528c
529      integer n,m
530      parameter (n = 60)
531      parameter (m = 2*n)
532      double precision areal(n,n), breal(n,n)
533      double precision aimg(n,n), bimg(n,n)
534#ifdef MIRROR
535      integer ndim, dims(2), chunk(2), p_mirror
536#else
537#  ifdef NEW_API
538      integer ndim, dims(2), chunk(2), p_mirror
539#  endif
540#endif
541      integer iv(m), jv(m)
542      logical status
543      integer g_a, g_b
544      integer iran, i,j, loop,nloop,maxloop, ilo, ihi, jlo, jhi, itmp
545      integer nproc, me, int, ij, inc, ii, jj, nnodes
546      parameter (maxloop = 100)
547      integer maxproc
548      parameter (maxproc = 4096)
549      double precision crap, real
550      double precision nwords
551      double complex   x, sum1, sum2, factor
552      integer lo(2), hi(2)
553      integer lprocs, inode, iproc, lproc
554#ifdef USE_RESTRICTED
555      integer num_rstrctd
556      integer rstrctd_list(maxproc/2)
557#endif
558#ifdef BLOCK_CYCLIC
559      integer block_size(2), proc_grid(2)
560#endif
561      intrinsic int
562      iran(i) = int(drand(0)*real(i)) + 1
563c
564      nproc = ga_nnodes()
565      me    = ga_nodeid()
566      inode = ga_cluster_nodeid()
567      lprocs = ga_cluster_nprocs(inode)
568      nnodes = ga_cluster_nnodes()
569      iproc = mod(me,lprocs)
570      nloop = Min(maxloop,n)
571#ifdef USE_RESTRICTED
572      num_rstrctd = nproc/2
573      if (num_rstrctd.eq.0) num_rstrctd = 1
574      do i = 1, num_rstrctd
575        rstrctd_list(i) = (num_rstrctd/2) + i-1
576      end do
577#endif
578#ifdef BLOCK_CYCLIC
579      block_size(1) = 32
580      block_size(2) = 32
581#ifdef USE_SCALAPACK_DISTR
582      if (mod(nproc,2).ne.0)
583     +  call ga_error("Available procs must be divisible by 2",0)
584      proc_grid(1) = 2
585      proc_grid(2) = nproc/2
586#endif
587#endif
588c
589c     a() is a local copy of what the global array should start as
590c
591      do j = 1, n
592         do i = 1, n
593#ifndef MIRROR
594            areal(i,j) = dble(i-1)
595            aimg(i,j) = dble((j-1)*n)
596#else
597            areal(i,j) = dble(inode) + dble(i-1)
598            aimg(i,j) = 0.0d00 + dble((j-1)*n)
599#endif
600            breal(i,j) = -1d0
601            bimg(i,j) =  1d0
602         enddo
603      enddo
604c
605c     Create  type
606c
607c
608c     Create a global array
609c
610c     print *,ga_nodeid(), ' creating array'
611      call ffflush(6)
612c     call setdbg(1)
613#ifdef NEW_API
614      ndim = 2
615      dims(1) = n
616      dims(2) = n
617      g_a = ga_create_handle()
618      call ga_set_data(g_a,ndim,dims,MT_DCPL)
619      call ga_set_array_name(g_a,'a')
620#ifdef USE_RESTRICTED
621      call ga_set_restricted(g_a, rstrctd_list, num_rstrctd)
622#endif
623#ifdef BLOCK_CYCLIC
624#ifdef USE_SCALAPACK_DISTR
625      call ga_set_block_cyclic_proc_grid(g_a,block_size,proc_grid)
626#else
627      call ga_set_block_cyclic(g_a,block_size)
628#endif
629#endif
630#  ifdef MIRROR
631      p_mirror = ga_pgroup_get_mirror()
632      call ga_set_pgroup(g_a,p_mirror)
633#  endif
634      status = ga_allocate(g_a)
635#else
636#  ifndef MIRROR
637      status = ga_create(MT_DCPL, n, n, 'a', 0, 0, g_a)
638#  else
639      ndim = 2
640      dims(1) = n
641      dims(2) = n
642      chunk(1) = 0
643      chunk(2) = 0
644      p_mirror = ga_pgroup_get_mirror()
645      status = nga_create_config(MT_DCPL, ndim, dims, 'a', chunk,
646     +                           p_mirror, g_a)
647#  endif
648#endif
649      if (.not. status) then
650         write(6,*) ' ga_create failed'
651         call ga_error('... exiting ',0)
652      endif
653#ifdef NEW_API
654      g_b = ga_create_handle()
655      call ga_set_data(g_b,ndim,dims,MT_DCPL)
656      call ga_set_array_name(g_b,'b')
657#ifdef BLOCK_CYCLIC
658#ifdef USE_SCALAPACK_DISTR
659      call ga_set_block_cyclic_proc_grid(g_b,block_size,proc_grid)
660#else
661      call ga_set_block_cyclic(g_b,block_size)
662#endif
663#endif
664#  ifdef MIRROR
665      call ga_set_pgroup(g_b,p_mirror)
666#  endif
667      if (.not.ga_allocate(g_b)) then
668#else
669#  ifndef MIRROR
670      if (.not. ga_create(MT_DCPL, n, n, 'b', 0, 0, g_b)) then
671#  else
672      if (.not. nga_create_config(MT_DCPL, ndim, dims, 'b', chunk,
673     _                            p_mirror, g_b)) then
674#  endif
675#endif
676         call ga_error('ga_create failed for second array ',0)
677      endif
678
679#ifndef MIRROR
680      call ga_distribution(g_a,me,ilo, ihi, jlo, jhi)
681#else
682      lproc = me - ga_cluster_procid(inode,0)
683      call ga_distribution(g_a, lproc, ilo, ihi, jlo, jhi)
684#endif
685      call ga_sync()
686c
687c     Zero the array
688c
689      if (me .eq. 0) then
690         write(6,21)
691 21      format('> Checking zero ... ')
692         call ffflush(6)
693      endif
694      call ga_zero(g_a)
695c
696c     Check that it is indeed zero
697c
698c      call ga_get(g_a, 1, n, 1, n, b, n)
699      lo(1) = 1
700      lo(2) = 1
701      hi(1) = n
702      hi(2) = n
703      call nga_get_field(g_a, lo, hi, 0, 8, breal, n)
704      call nga_get_field(g_a, lo, hi, 8, 8, bimg, n)
705      call ga_sync()
706      do i = 1, n
707         do j = 1, n
708            if(breal(i,j).ne.(0d0,0d0)) then
709               write(6,*) me,' real zero ', i, j, breal(i,j)
710               call ffflush(6)
711              call ga_error('... exiting ',0)
712            endif
713            if(bimg(i,j).ne.(0d0,0d0)) then
714               write(6,*) me,' img zero ', i, j, bimg(i,j)
715               call ffflush(6)
716              call ga_error('... exiting ',0)
717            endif
718         enddo
719      enddo
720      if (me.eq.0) then
721         write(6,*)
722         write(6,*) ' ga_zero is OK'
723         write(6,*)
724      endif
725      call ga_sync()
726c
727c     Each node fills in disjoint sections of the array
728c
729      if (me .eq. 0) then
730         write(6,2)
731 2       format('> Checking disjoint put ... ')
732         call ffflush(6)
733      endif
734      call ga_sync()
735c
736      inc = (n-1)/20 + 1
737      ij = 0
738      do j = 1, n, inc
739         do i = 1, n, inc
740#ifndef MIRROR
741            if (mod(ij,nproc) .eq. me) then
742#else
743            if (mod(ij,lprocs) .eq. iproc) then
744#endif
745               ilo = i
746               ihi = min(i+inc, n)
747               jlo = j
748               jhi = min(j+inc, n)
749c               call ga_put(g_a, ilo, ihi, jlo, jhi, a(ilo, jlo), n)
750               lo(1) = ilo
751               lo(2) = jlo
752               hi(1) = ihi
753               hi(2) = jhi
754               call nga_put_field(g_a, lo, hi, 0, 8, areal(ilo,jlo),n)
755               call nga_put_field(g_a, lo, hi, 8, 8, aimg(ilo,jlo),n)
756            endif
757            ij = ij + 1
758         enddo
759      enddo
760      call ga_sync()
761c      call ga_print(g_a)
762c
763c     All nodes check all of a
764c
765      call util_qfill(n*n, 0.0d0, breal, 1)
766      call util_qfill(n*n, 0d0, bimg, 1)
767c     call ga_get(g_a, 1, n, 1, n, b, n)
768      lo(1) = 1
769      lo(2) = 1
770      hi(1) = n
771      hi(2) = n
772      call nga_get_field(g_a, lo, hi, 0, 8, breal, n)
773      call nga_get_field(g_a, lo, hi, 8, 8, bimg, n)
774c
775      do i = 1, n
776         do j = 1, n
777            if (breal(i,j) .ne. areal(i,j)) then
778               write(6,*) ' real put ', me, i, j, areal(i,j),breal(i,j)
779               call ffflush(6)
780               call ga_error('... exiting ',0)
781            endif
782            if (bimg(i,j) .ne. aimg(i,j)) then
783               write(6,*) ' img put ', me, i, j, aimg(i,j),bimg(i,j)
784               call ffflush(6)
785               call ga_error('... exiting ',0)
786            endif
787         enddo
788      enddo
789      if (me.eq.0) then
790         write(6,*)
791         write(6,*) ' ga_put is OK'
792         write(6,*)
793      endif
794      call ga_sync()
795c
796c     Now check nloop random gets from each node
797c
798      if (me .eq. 0) then
799         write(6,5) nloop
800 5       format('> Checking random get (',i5,' calls)...')
801         call ffflush(6)
802      endif
803      call ga_sync()
804c
805      nwords = 0
806c
807      crap = drand(ga_nodeid()*51 +1) !Different seed for each proc
808      do loop = 1, nloop
809         ilo = iran(loop)
810         ihi = iran(loop)
811         if (ihi.lt. ilo) then
812            itmp = ihi
813            ihi = ilo
814            ilo = itmp
815         endif
816         jlo = iran(loop)
817         jhi = iran(loop)
818         if (jhi.lt. jlo) then
819            itmp = jhi
820            jhi = jlo
821            jlo = itmp
822         endif
823c
824         nwords = nwords + (ihi-ilo+1)*(jhi-jlo+1)
825c
826         call util_qfill(n*n, 0.0d0, breal, 1)
827         call util_qfill(n*n, 0d0, bimg, 1)
828c         call ga_get(g_a, ilo, ihi, jlo, jhi, b(ilo, jlo), n)
829         lo(1) = ilo
830         lo(2) = jlo
831         hi(1) = ihi
832         hi(2) = jhi
833         call nga_get_field(g_a, lo, hi, 0, 8, breal(ilo,jlo),n)
834         call nga_get_field(g_a, lo, hi, 8, 8, bimg(ilo,jlo),n)
835         if (me .eq. 0 .and. mod(loop-1, max(1,nloop/20)).eq.0) then
836            write(6,1) loop, me, ilo, ihi, jlo, jhi, nwords
837 1          format(' call ',i5, ' node ',i2,' checking get ',4i4,
838     $           ' total ',d9.2)
839            call ffflush(6)
840         endif
841         do j = jlo, jhi
842            do i = ilo, ihi
843               if (breal(i,j) .ne. areal(i,j)) then
844                  write(6,*)'real error:', i, j, breal(i,j), areal(i,j)
845                  call ga_error('... exiting ',0)
846               endif
847               if (bimg(i,j) .ne. aimg(i,j)) then
848                  write(6,*)'img error:', i, j, bimg(i,j), aimg(i,j)
849                  call ga_error('... exiting ',0)
850               endif
851            enddo
852         enddo
853c
854      enddo
855      if (me.eq.0) then
856         write(6,*)
857         write(6,*) ' ga_get is OK'
858         write(6,*)
859         call ffflush(6)
860      endif
861      call ga_sync()
862c
863c     Check the ga_copy function
864c
865      if (me .eq. 0) then
866         write(6,*)
867         write(6,*)'> Checking copy'
868         write(6,*)
869         call ffflush(6)
870      endif
871      call ga_sync()
872#ifndef MIRROR
873      if(me.eq.0) then
874c         call ga_put(g_a, 1, n, 1, n, a, n)
875         lo(1) = 1
876         lo(2) = 1
877         hi(1) = n
878         hi(2) = n
879         call nga_put_field(g_a, lo, hi, 0, 8, areal, n)
880         call nga_put_field(g_a, lo, hi, 8, 8, aimg, n)
881      endif
882#else
883      if(iproc.eq.0) then
884c         call ga_put(g_a, 1, n, 1, n, a, n)
885         lo(1) = 1
886         lo(2) = 1
887         hi(1) = n
888         hi(2) = n
889         call nga_put_field(g_a, lo, hi, 0, 8, areal, n)
890         call nga_put_field(g_a, lo, hi, 8, 8, aimg, n)
891      endif
892#endif
893      call ga_copy(g_a, g_b)
894c      call ga_get(g_b, 1, n, 1, n, b, n)
895      lo(1) = 1
896      lo(2) = 1
897      hi(1) = n
898      hi(2) = n
899      call nga_get_field(g_b, lo, hi, 0, 8, breal, n)
900      call nga_get_field(g_b, lo, hi, 8, 8, bimg, n)
901      do j = 1, n
902         do i = 1, n
903            if (breal(i,j) .ne. areal(i,j)) then
904               write(6,*) ' copy ', me, i, j, areal(i,j), breal(i,j)
905               call ga_error('... exiting ',0)
906            endif
907            if (bimg(i,j) .ne. aimg(i,j)) then
908               write(6,*) ' copy ', me, i, j, aimg(i,j), bimg(i,j)
909               call ga_error('... exiting ',0)
910            endif
911         enddo
912      enddo
913      if (me.eq.0) then
914         write(6,*)
915         write(6,*) ' copy is OK '
916         write(6,*)
917      endif
918c
919c     Delete the global arrays
920c
921      status = ga_destroy(g_b)
922      status = ga_destroy(g_a)
923
924      end
925
926      subroutine util_qfill(n,val,a,ia)
927      implicit none
928      double precision  a(*), val
929      integer n, ia, i
930c
931c     initialise double complex array to scalar value
932c
933      if (ia.eq.1) then
934         do 10 i = 1, n
935            a(i) = val
936 10      continue
937      else
938         do 20 i = 1,(n-1)*ia+1,ia
939            a(i) = val
940 20      continue
941      endif
942c
943      end
944
945