~ubuntu-branches/ubuntu/utopic/nwchem/utopic

« back to all changes in this revision

Viewing changes to src/tools/ga-5-1/global/testing/patch.F

  • Committer: Package Import Robot
  • Author(s): Michael Banck, Daniel Leidert, Andreas Tille, Michael Banck
  • Date: 2013-07-04 12:14:55 UTC
  • mfrom: (1.1.2)
  • Revision ID: package-import@ubuntu.com-20130704121455-5tvsx2qabor3nrui
Tags: 6.3-1
* New upstream release.
* Fixes anisotropic properties (Closes: #696361).
* New features include:
  + Multi-reference coupled cluster (MRCC) approaches
  + Hybrid DFT calculations with short-range HF 
  + New density-functionals including Minnesota (M08, M11) and HSE hybrid
    functionals
  + X-ray absorption spectroscopy (XAS) with TDDFT
  + Analytical gradients for the COSMO solvation model
  + Transition densities from TDDFT 
  + DFT+U and Electron-Transfer (ET) methods for plane wave calculations
  + Exploitation of space group symmetry in plane wave geometry optimizations
  + Local density of states (LDOS) collective variable added to Metadynamics
  + Various new XC functionals added for plane wave calculations, including
    hybrid and range-corrected ones
  + Electric field gradients with relativistic corrections 
  + Nudged Elastic Band optimization method
  + Updated basis sets and ECPs 

[ Daniel Leidert ]
* debian/watch: Fixed.

[ Andreas Tille ]
* debian/upstream: References

[ Michael Banck ]
* debian/upstream (Name): New field.
* debian/patches/02_makefile_flags.patch: Refreshed.
* debian/patches/06_statfs_kfreebsd.patch: Likewise.
* debian/patches/07_ga_target_force_linux.patch: Likewise.
* debian/patches/05_avoid_inline_assembler.patch: Removed, no longer needed.
* debian/patches/09_backported_6.1.1_fixes.patch: Likewise.
* debian/control (Build-Depends): Added gfortran-4.7 and gcc-4.7.
* debian/patches/10_force_gcc-4.7.patch: New patch, explicitly sets
  gfortran-4.7 and gcc-4.7, fixes test suite hang with gcc-4.8 (Closes:
  #701328, #713262).
* debian/testsuite: Added tests for COSMO analytical gradients and MRCC.
* debian/rules (MRCC_METHODS): New variable, required to enable MRCC methods.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
#if HAVE_CONFIG_H
2
 
#   include "config.fh"
3
 
#endif
4
 
C $iD: PATch.F,v 1.21 2002/10/15 23:35:52 vinod Exp $
5
 
#if (defined(CRAY) && !defined(__crayx1)) || defined(KSR)
6
 
#   define xgemm SGEMM
7
 
#   define ygemm CGEMM
8
 
#else
9
 
#   define xgemm TEST_DGEMM
10
 
#   define ygemm TEST_ZGEMM
11
 
#endif
12
 
#if defined(FUJITSU) || defined(CRAY_YMP)
13
 
# define THRESH  1.0d-10
14
 
#else
15
 
# define THRESH  1.0d-20
16
 
#endif
17
 
#define MISMATCH(x,y) abs(x-y)/max(1,abs(x)).gt.THRESH
18
 
c
19
 
c#define BLOCK_CYCLIC
20
 
c#define USE_SCALAPACK_DISTR
21
 
#ifdef USE_SCALAPACK_DISTR
22
 
#define BLOCK_CYCLIC
23
 
#endif
24
 
c
25
 
      program test
26
 
      implicit none
27
 
#include "mafdecls.fh"
28
 
#include "global.fh"
29
 
      logical status
30
 
c
31
 
c***  Initialize a message passing library
32
 
c
33
 
#include "mp3.fh"
34
 
c
35
 
      call ga_initialize()
36
 
      if(ga_nodeid().eq.0)then
37
 
         print *,' GA initialized'
38
 
         call ffflush(6)
39
 
      endif
40
 
c
41
 
      status = ma_init(MT_DBL, 500000, 900000/ga_nnodes())
42
 
      if (.not. status)call ga_error( 'ma_init failed', -1)
43
 
      if(ga_nodeid().eq.0)then
44
 
         print *,' '
45
 
         print *,'CHECKING PATCH OPERATIONS FOR DOUBLES '
46
 
         print *,' '
47
 
      endif
48
 
      call dpatch_test()
49
 
c
50
 
      if(ga_nodeid().eq.0)then
51
 
         print *,' '
52
 
         print *,'CHECKING PATCH OPERATIONS FOR DOUBLE COMPLEX'
53
 
         print *,' '
54
 
      endif
55
 
      call zpatch_test()
56
 
57
 
      if(ga_nodeid().eq.0)then
58
 
         print *,' '
59
 
         print *,'CHECKING PATCH OPERATIONS FOR SINGLE PRECISION'
60
 
         print *,' '
61
 
      endif
62
 
      call spatch_test()
63
 
c     
64
 
      if(ga_nodeid().eq.0) print *,'All tests successful '
65
 
c
66
 
      call ga_terminate()
67
 
c
68
 
      call MP_FINALIZE()
69
 
      end
70
 
 
71
 
 
72
 
 
73
 
      subroutine dpatch_test()
74
 
      implicit none
75
 
#include "mafdecls.fh"
76
 
#include "global.fh"
77
 
#include "testutil.fh"
78
 
c
79
 
      integer n,m
80
 
      parameter (n = 128)
81
 
      parameter (m = n*n)
82
 
      double precision a(n,n), b(n,n), c(n,n), buf(m), val
83
 
      double precision alpha, beta
84
 
      integer nproc, me 
85
 
      integer i, j, ailo, ajlo, bilo, bjlo, base, iran
86
 
      integer aihi, ajhi, bihi, bjhi
87
 
      integer g_a, g_b, g_c 
88
 
      integer rows, cols, loop
89
 
      logical status, dist_same 
90
 
#ifdef BLOCK_CYCLIC
91
 
      integer ndim, dims(2)
92
 
      integer block_size(2), proc_grid(2)
93
 
#endif
94
 
      iran(i) = int(drand(1)*real(i)) + 1
95
 
      dist_same = .false.
96
 
c
97
 
      me = ga_nodeid()
98
 
      nproc = ga_nnodes()
99
 
#ifdef BLOCK_CYCLIC
100
 
      block_size(1) = 32
101
 
      block_size(2) = 32
102
 
#ifdef USE_SCALAPACK_DISTR
103
 
      if (mod(nproc,2).ne.0)
104
 
     +  call ga_error("Available procs must be divisible by 2",0)
105
 
      proc_grid(1) = 2
106
 
      proc_grid(2) = nproc/2
107
 
#endif
108
 
#endif
109
 
 
110
 
c
111
 
      do j = 1, n
112
 
         do i = 1, n
113
 
            a(i,j) = i-1 + (j-1)*n
114
 
            b(i,j) = i+j 
115
 
         enddo
116
 
      enddo
117
 
c
118
 
c***  Create a global array
119
 
c
120
 
#ifndef BLOCK_CYCLIC
121
 
      status = ga_create(MT_DBL, n, n, 'a', 0, 0, g_a)
122
 
#else
123
 
      g_a = ga_create_handle()
124
 
      ndim = 2
125
 
      dims(1) = n
126
 
      dims(2) = n
127
 
      call ga_set_data(g_a,ndim,dims,MT_DBL)
128
 
      call ga_set_array_name(g_a,'a')
129
 
#ifdef USE_SCALAPACK_DISTR
130
 
      call ga_set_block_cyclic_proc_grid(g_a,block_size,proc_grid)
131
 
#else
132
 
      call ga_set_block_cyclic(g_a,block_size)
133
 
#endif
134
 
      status = ga_allocate(g_a)
135
 
#endif
136
 
      if (.not. status) then
137
 
         write(6,*) ' ga_create failed'
138
 
         call ffflush(6)
139
 
         call ga_error('... exiting ',0)
140
 
      endif
141
 
c
142
 
      if(dist_same) then
143
 
         status = ga_duplicate(g_a, g_b, 'a_duplicated')
144
 
         if(.not.ga_compare_distr(g_a, g_b))
145
 
     $           call ga_error("g_b distribution different",0) 
146
 
         status = ga_duplicate(g_a, g_c, 'a_duplicated_again')
147
 
         if(.not.ga_compare_distr(g_a, g_c))
148
 
     $           call ga_error("g_c distribution different",0) 
149
 
      else
150
 
#ifndef BLOCK_CYCLIC
151
 
         status = ga_create(MT_DBL, n, n, 'b', 0, n, g_b)
152
 
#else
153
 
         g_b = ga_create_handle()
154
 
         ndim = 2
155
 
         dims(1) = n
156
 
         dims(2) = n
157
 
         call ga_set_data(g_b,ndim,dims,MT_DBL)
158
 
         call ga_set_array_name(g_b,'b')
159
 
#ifdef USE_SCALAPACK_DISTR
160
 
         call ga_set_block_cyclic_proc_grid(g_b,block_size,proc_grid)
161
 
#else
162
 
         call ga_set_block_cyclic(g_b,block_size)
163
 
#endif
164
 
         status = ga_allocate(g_b)
165
 
#endif
166
 
         if (.not. status) call ga_error('ga_create failed:b',0) 
167
 
#ifndef BLOCK_CYCLIC
168
 
         status = ga_create(MT_DBL, n, n, 'c', n, 0, g_c)
169
 
#else
170
 
         g_c = ga_create_handle()
171
 
         ndim = 2
172
 
         dims(1) = n
173
 
         dims(2) = n
174
 
         call ga_set_data(g_c,ndim,dims,MT_DBL)
175
 
         call ga_set_array_name(g_c,'c')
176
 
#ifdef USE_SCALAPACK_DISTR
177
 
         call ga_set_block_cyclic_proc_grid(g_c,block_size,proc_grid)
178
 
#else
179
 
         call ga_set_block_cyclic(g_c,block_size)
180
 
#endif
181
 
         status = ga_allocate(g_c)
182
 
#endif
183
 
         if (.not. status) call ga_error('ga_create failed:c',0) 
184
 
      endif
185
 
c
186
 
c***
187
 
      if (me .eq. 0) then
188
 
         print *, ' '
189
 
         write(6,*)'> Checking ga_fill_patch ... '
190
 
         call ffflush(6)
191
 
      endif
192
 
c
193
 
      val = 1d0
194
 
      call ga_fill_patch(g_a, 2,n/2, 2,n, val)
195
 
*     call ga_print(g_a, 1)
196
 
      do j = 2+me, n, nproc 
197
 
         call ga_get(g_a, 1,n/2, j,j, buf,n/2) 
198
 
         do i = 2, n/2
199
 
            if(buf(i) .ne.val ) then
200
 
               print *,me, ' error ',i,j, buf(i),val 
201
 
               call ga_error('exiting ...',0)
202
 
            endif
203
 
         enddo
204
 
      enddo
205
 
c
206
 
      call ga_sync()
207
 
      if (me .eq. 0) then
208
 
         write(6,*)'   OK '
209
 
         call ffflush(6)
210
 
      endif
211
 
c  
212
 
      do j = 1+me, n, nproc 
213
 
         call ga_put(g_a,1,n,j,j,a(1,j),n) 
214
 
      enddo
215
 
c***
216
 
      if (me .eq. 0) then
217
 
         print *, ' '
218
 
         write(6,*)'> Checking ga_copy_patch ... '
219
 
         call ffflush(6)
220
 
      endif
221
 
c
222
 
      do loop =1, 10
223
 
              ailo = iran(n/2)
224
 
              ajlo = iran(n/2)
225
 
              aihi = min(n, -1+ailo+n/2)
226
 
              ajhi = min(n, -1+ajlo+n/4)
227
 
              rows = aihi -ailo+1
228
 
              cols = ajhi -ajlo +1
229
 
c
230
 
              bilo = iran(n/3)
231
 
              bjlo = iran(n/3)
232
 
              bihi = bilo + rows -1
233
 
              bjhi = bjlo + cols -1
234
 
              if (me .eq. 0) then
235
 
                write(6,'(2x,1h[,4i4,1h],5h --> ,1h[,4i4,1h])')
236
 
     $                ailo,aihi,ajlo,ajhi,
237
 
     $                bilo,bihi, bjlo, bjhi
238
 
 
239
 
                call ffflush(6)
240
 
              endif
241
 
 
242
 
c
243
 
              call ga_copy_patch('n', g_a, ailo, aihi, ajlo, ajhi,
244
 
     &                                g_b, bilo, bihi, bjlo, bjhi)
245
 
c             call ga_print(g_a,1)
246
 
c             call ga_print(g_b,1)
247
 
              call ga_get(g_b,bilo,bihi, bjlo, bjhi, buf, rows)
248
 
              base = 0
249
 
              do j = ajlo, ajhi
250
 
                 if(Mod(j,nproc).eq.me) then
251
 
                    do i = ailo, aihi
252
 
                       base = base+1
253
 
                       if(buf(base) .ne. a(i,j)) then
254
 
                          print *,me, ' error ',i,j, buf(base), a(i,j)
255
 
                          call ga_error('exiting ...',0)
256
 
                       endif
257
 
                    enddo
258
 
                 else
259
 
                    base = base + rows
260
 
                 endif
261
 
              enddo
262
 
      enddo
263
 
c
264
 
      ailo = iran(n/2)
265
 
      ajlo = iran(n/2)
266
 
      bilo = iran(n/2)
267
 
      bjlo = iran(n/2)
268
 
c
269
 
#ifdef BLOCK_CYCLIC
270
 
      if (me .eq. 0) then
271
 
         write(6,*)'  without transpose:   OK '
272
 
         call ffflush(6)
273
 
      endif
274
 
#else
275
 
      if (me .eq. 0) then
276
 
                write(6,'(2x,1h[,4i4,1h],5h --> ,1h[,4i4,1h])')
277
 
     $                ailo,aihi,ajlo,ajhi,
278
 
     $                bilo,bihi, bjlo, bjhi
279
 
      endif
280
 
 
281
 
      call ga_copy_patch('n', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
282
 
     $                        g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
283
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,1)
284
 
*     call ga_print(g_b, 1)
285
 
      call ga_get(g_b,bilo,bilo+n/3,bjlo, bjlo+n/2,buf,n/3+1) 
286
 
      base = 0
287
 
      do j = ajlo, ajlo+n/3 
288
 
         if(Mod(j,nproc).eq.me) then
289
 
            do i = ailo, ailo+n/2 
290
 
               base = base+1
291
 
               if(buf(base) .ne. a(i,j)) then
292
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)
293
 
                  call ga_error('exiting ...',0)
294
 
               endif
295
 
            enddo
296
 
         else
297
 
            base = base +n/2+1
298
 
         endif
299
 
      enddo
300
 
      call ga_sync()
301
 
c
302
 
      if (me .eq. 0) then
303
 
         write(6,*)'  without transpose:   OK '
304
 
         call ffflush(6)
305
 
      endif
306
 
307
 
      if (me .eq. 0) then
308
 
                write(6,'(2x,1h[,4i4,2h]~,5h --> ,1h[,4i4,1h])')
309
 
     $                ailo,aihi,ajlo,ajhi,
310
 
     $                bilo,bihi, bjlo, bjhi
311
 
      endif
312
 
      call ga_copy_patch('t', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
313
 
     $                        g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
314
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,1)
315
 
*     call ga_print(g_b, 1)
316
 
      call ga_get(g_b,bilo,bilo+n/3,bjlo, bjlo+n/2,buf,n/3+1)
317
 
      base = 0
318
 
      do i = ailo, ailo+n/2
319
 
         if(Mod(j,nproc).eq.me) then
320
 
         do j = ajlo, ajlo+n/3
321
 
               base = base+1
322
 
               if(buf(base) .ne. a(i,j)) then
323
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)
324
 
                  call ga_error('exiting ...',0)
325
 
               endif
326
 
            enddo
327
 
         else
328
 
            base = base +n/3+1
329
 
         endif
330
 
      enddo
331
 
      call ga_sync()
332
 
c
333
 
      if (me .eq. 0) then
334
 
         write(6,*)'  transposed:   OK '
335
 
         call ffflush(6)
336
 
      endif
337
 
#endif
338
 
c
339
 
c***
340
 
      if (me .eq. 0) then
341
 
         print *, ' '
342
 
         write(6,*)'> Checking ga_scale_patch ... '
343
 
         call ffflush(6)
344
 
      endif
345
 
      call ga_copy_patch('n', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
346
 
     $                        g_b, bilo,bilo+n/2, bjlo, bjlo+n/3)
347
 
      val = 1.d0
348
 
      call ga_scale_patch(g_b, bilo,bilo+n/2,bjlo, bjlo+n/3,val)
349
 
      call ga_get(g_b,bilo,bilo+n/2,bjlo, bjlo+n/3,buf,n/2+1)
350
 
      base = 0
351
 
      do j = ajlo, ajlo+n/3
352
 
         if(Mod(j,nproc).eq.me) then
353
 
            do i = ailo, ailo+n/2
354
 
               base = base+1
355
 
               if(buf(base) .ne. a(i,j)*val) then
356
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)*val
357
 
                  call ga_error('exiting ...',0)
358
 
               endif
359
 
            enddo
360
 
         else
361
 
            base = base +n/2+1
362
 
         endif
363
 
      enddo
364
 
      call ga_sync()
365
 
c
366
 
      if (me .eq. 0) then
367
 
         write(6,*)'  OK '
368
 
         call ffflush(6)
369
 
      endif
370
 
c
371
 
c***
372
 
      if (me .eq. 0) then
373
 
         print *, ' '
374
 
         write(6,*)'> Checking ga_add_patch ... '
375
 
         call ffflush(6)
376
 
      endif
377
 
      alpha = .1d0
378
 
      beta = .2d0
379
 
      call ga_zero(g_c)
380
 
      call ga_add_patch(alpha, g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
381
 
     $                   beta, g_b, bilo,bilo+n/2, bjlo, bjlo+n/3,
382
 
     $                         g_c, bilo,bilo+n/2, bjlo, bjlo+n/3)
383
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/3,buf,n/2+1)
384
 
      base = 0
385
 
      val = val*beta + alpha
386
 
      do j = ajlo, ajlo+n/3
387
 
         if(Mod(j,nproc).eq.me) then
388
 
            do i = ailo, ailo+n/2
389
 
               base = base+1
390
 
               if(ABS(buf(base)- a(i,j)*val).gt.1d-5) then
391
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)*val
392
 
                  call ga_error('exiting ...',0)
393
 
               endif
394
 
            enddo
395
 
         else
396
 
            base = base +n/2+1
397
 
         endif
398
 
      enddo
399
 
      call ga_sync()
400
 
c
401
 
      if (me .eq. 0) then
402
 
         write(6,*)'  OK '
403
 
         call ffflush(6)
404
 
      endif
405
 
c
406
 
c***
407
 
      if (me .eq. 0) then
408
 
         print *, ' '
409
 
         write(6,*)'> Checking ga_ddot_patch ... '
410
 
         call ffflush(6)
411
 
      endif
412
 
      alpha= ga_ddot_patch(g_a,'n', ailo,ailo+n/2, ajlo, ajlo+n/3,
413
 
     $                     g_c,'n', bilo,bilo+n/2, bjlo, bjlo+n/3)
414
 
      beta = 0d0
415
 
      do j = ajlo, ajlo+n/3
416
 
            do i = ailo, ailo+n/2
417
 
               beta = beta + a(i,j)**2
418
 
            enddo
419
 
      enddo
420
 
 
421
 
      if(ABS(beta*val- alpha).gt.1d-6*alpha) then
422
 
             print *,me, ' error ', beta*val, alpha
423
 
             call ga_error('exiting ...',0)
424
 
      endif
425
 
      call ga_sync()
426
 
c
427
 
      if (me .eq. 0) then
428
 
         write(6,*)'  OK '
429
 
         call ffflush(6)
430
 
      endif
431
 
c
432
 
c......................................................
433
 
      if (me .eq. 0) then
434
 
         print *, ' '
435
 
         write(6,*)'> Checking ga_matmul_patch ... '
436
 
         call ffflush(6)
437
 
      endif
438
 
      do j = 1+me, n, nproc
439
 
         call ga_put(g_b,1,n,j,j,b(1,j),n)
440
 
      enddo
441
 
      call ga_sync()
442
 
      call ga_matmul_patch('n','n', 1d0, 0d0, 
443
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
444
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
445
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
446
 
      call xgemm('n','n',n/2+1,n/2+1,n/3+1,1d0,a(ailo,ajlo), n,
447
 
     $            b(bilo,bjlo),n, 0d0, c, n)
448
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3)
449
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
450
 
*     call ga_print(g_c, 1)
451
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
452
 
      base = 0
453
 
      do j = 1, 1+n/2
454
 
         if(Mod(j,nproc).eq.me) then
455
 
            do i = 1, 1+n/2
456
 
               base = base+1
457
 
               if(ABS(buf(base)- c(i,j)).gt.1d-8) then
458
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
459
 
                  call ga_error('exiting ...',0)
460
 
               endif
461
 
            enddo
462
 
         else
463
 
            base = base +n/2+1
464
 
         endif
465
 
      enddo
466
 
c
467
 
      call ga_sync()
468
 
      if (me .eq. 0) then
469
 
         write(6,*)'  a*b: OK '
470
 
         call ffflush(6)
471
 
      endif
472
 
c
473
 
      call ga_sync()
474
 
      call ga_matmul_patch('t','n', 1d0, 0d0,
475
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
476
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
477
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
478
 
      call xgemm('t','n',n/2+1,n/2+1,n/3+1,1d0,a(ajlo,ailo), n,
479
 
     $            b(bilo,bjlo),n, 0d0, c, n)
480
 
*     call ga_print(g_a,1) 
481
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
482
 
*     call ga_print(g_c, 1)
483
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
484
 
      base = 0
485
 
      do j = 1, 1+n/2
486
 
         if(Mod(j,nproc).eq.me) then
487
 
            do i = 1, 1+n/2
488
 
               base = base+1
489
 
               if(ABS(buf(base)- c(i,j)).gt.1d-8) then
490
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
491
 
                  call ga_error('exiting ...',0)
492
 
               endif
493
 
            enddo
494
 
         else
495
 
            base = base +n/2+1
496
 
         endif
497
 
      enddo
498
 
c
499
 
      call ga_sync()
500
 
      if (me .eq. 0) then
501
 
         write(6,*)'  trans(a)*b: OK '
502
 
         call ffflush(6)
503
 
      endif
504
 
c
505
 
      call ga_sync()
506
 
      call ga_matmul_patch('n','t', 1d0, 0d0,
507
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
508
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
509
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
510
 
      call xgemm('n','t',n/2+1,n/2+1,n/3+1,1d0,a(ailo,ajlo), n,
511
 
     $            b(bjlo,bilo),n, 0d0, c, n)
512
 
*     call ga_print(g_a,1)
513
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
514
 
*     call ga_print(g_c, 1)
515
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
516
 
      base = 0
517
 
      do j = 1, 1+n/2
518
 
         if(Mod(j,nproc).eq.me) then
519
 
            do i = 1, 1+n/2
520
 
               base = base+1
521
 
               if(ABS(buf(base)- c(i,j)).gt.1d-8) then
522
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
523
 
                  call ga_error('exiting ...',0)
524
 
               endif
525
 
            enddo
526
 
         else
527
 
            base = base +n/2+1
528
 
         endif
529
 
      enddo
530
 
c
531
 
      call ga_sync()
532
 
      if (me .eq. 0) then
533
 
         write(6,*)'  a*trans(b): OK '
534
 
         call ffflush(6)
535
 
      endif
536
 
c
537
 
      call ga_sync()
538
 
      call ga_matmul_patch('t','t', 1d0, 0d0,
539
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
540
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
541
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
542
 
      call xgemm('t','t',n/2+1,n/2+1,n/3+1,1d0,a(ajlo,ailo), n,
543
 
     $            b(bjlo,bilo),n, 0d0, c, n)
544
 
*     call ga_print(g_a,1)
545
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
546
 
*     call ga_print(g_c, 1)
547
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
548
 
      base = 0
549
 
      do j = 1, 1+n/2
550
 
         if(Mod(j,nproc).eq.me) then
551
 
            do i = 1, 1+n/2
552
 
               base = base+1
553
 
               if(ABS(buf(base)- c(i,j)).gt.1d-8) then
554
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
555
 
                  call ga_error('exiting ...',0)
556
 
               endif
557
 
            enddo
558
 
         else
559
 
            base = base +n/2+1
560
 
         endif
561
 
      enddo
562
 
c
563
 
      call ga_sync()
564
 
      if (me .eq. 0) then
565
 
         write(6,*)'  trans(a)*trans(b): OK '
566
 
         call ffflush(6)
567
 
      endif
568
 
      status = ga_destroy(g_a)
569
 
      status = status .and. ga_destroy(g_b)
570
 
      status = status .and. ga_destroy(g_c)
571
 
      if(.not. status) print *, 'ga_destroy failed'
572
 
c
573
 
      end
574
 
 
575
 
 
576
 
 
577
 
 
578
 
 
579
 
      subroutine zpatch_test()
580
 
      implicit none
581
 
#include "mafdecls.fh"
582
 
#include "global.fh"
583
 
#include "testutil.fh"
584
 
c
585
 
      integer n,m
586
 
      parameter (n = 128)
587
 
      parameter (m = n*n)
588
 
      double complex a(n,n), b(n,n), c(n,n), buf(m),  val
589
 
      double complex alpha, beta
590
 
      integer nproc, me 
591
 
      integer i, j, ailo, ajlo, bilo, bjlo, base, iran
592
 
      integer aihi, ajhi, bihi, bjhi
593
 
      integer g_a, g_b, g_c 
594
 
      integer rows, cols, loop
595
 
      logical status, dist_same 
596
 
#ifdef BLOCK_CYCLIC
597
 
      integer ndim, dims(2)
598
 
      integer block_size(2), proc_grid(2)
599
 
#endif
600
 
      iran(i) = int(drand(1)*real(i)) + 1
601
 
      dist_same = .false.
602
 
c
603
 
      me = ga_nodeid()
604
 
      nproc = ga_nnodes()
605
 
#ifdef BLOCK_CYCLIC
606
 
      block_size(1) = 32
607
 
      block_size(2) = 32
608
 
#ifdef USE_SCALAPACK_DISTR
609
 
      if (mod(nproc,2).ne.0)
610
 
     +  call ga_error("Available procs must be divisible by 2",0)
611
 
      proc_grid(1) = 2
612
 
      proc_grid(2) = nproc/2
613
 
#endif
614
 
#endif
615
 
c
616
 
      do j = 1, n
617
 
         do i = 1, n
618
 
            a(i,j) = cmplx(dble(i-1), dble((j-1)*n))
619
 
            b(i,j) = cmplx(dble(i+j),1d0)
620
 
         enddo
621
 
      enddo
622
 
c
623
 
c***  Create a global array
624
 
c
625
 
#ifndef BLOCK_CYCLIC
626
 
      status = ga_create(MT_DCPL, n, n, 'a', 0, 0, g_a)
627
 
#else
628
 
      g_a = ga_create_handle()
629
 
      ndim = 2
630
 
      dims(1) = n
631
 
      dims(2) = n
632
 
      call ga_set_data(g_a,ndim,dims,MT_DCPL)
633
 
      call ga_set_array_name(g_a,'a')
634
 
#ifdef USE_SCALAPACK_DISTR
635
 
      call ga_set_block_cyclic_proc_grid(g_a,block_size,proc_grid)
636
 
#else
637
 
      call ga_set_block_cyclic(g_a,block_size)
638
 
#endif
639
 
      status = ga_allocate(g_a)
640
 
#endif
641
 
      if (.not. status) then
642
 
         write(6,*) ' ga_create failed'
643
 
         call ffflush(6)
644
 
         call ga_error('... exiting ',0)
645
 
      endif
646
 
c
647
 
      if(dist_same) then
648
 
         status = ga_duplicate(g_a, g_b, 'a_duplicated')
649
 
         if(.not.ga_compare_distr(g_a, g_b))
650
 
     $           call ga_error("g_b distribution different",0) 
651
 
         status = ga_duplicate(g_a, g_c, 'a_duplicated_again')
652
 
         if(.not.ga_compare_distr(g_a, g_c))
653
 
     $           call ga_error("g_c distribution different",0) 
654
 
      else
655
 
#ifndef BLOCK_CYCLIC
656
 
         status = ga_create(MT_DCPL, n, n, 'b', 0, n, g_b)
657
 
#else
658
 
         g_b = ga_create_handle()
659
 
         ndim = 2
660
 
         dims(1) = n
661
 
         dims(2) = n
662
 
         call ga_set_data(g_b,ndim,dims,MT_DCPL)
663
 
         call ga_set_array_name(g_b,'b')
664
 
#ifdef USE_SCALAPACK_DISTR
665
 
         call ga_set_block_cyclic_proc_grid(g_b,block_size,proc_grid)
666
 
#else
667
 
         call ga_set_block_cyclic(g_b,block_size)
668
 
#endif
669
 
         status = ga_allocate(g_b)
670
 
#endif
671
 
         if (.not. status) call ga_error('ga_create failed:b',0) 
672
 
#ifndef BLOCK_CYCLIC
673
 
         status = ga_create(MT_DCPL, n, n, 'c', n, 0, g_c)
674
 
#else
675
 
         g_c = ga_create_handle()
676
 
         ndim = 2
677
 
         dims(1) = n
678
 
         dims(2) = n
679
 
         call ga_set_data(g_c,ndim,dims,MT_DCPL)
680
 
         call ga_set_array_name(g_c,'c')
681
 
#ifdef USE_SCALAPACK_DISTR
682
 
         call ga_set_block_cyclic_proc_grid(g_c,block_size,proc_grid)
683
 
#else
684
 
         call ga_set_block_cyclic(g_c,block_size)
685
 
#endif
686
 
         status = ga_allocate(g_c)
687
 
#endif
688
 
         if (.not. status) call ga_error('ga_create failed:c',0) 
689
 
      endif
690
 
c
691
 
c***
692
 
      if (me .eq. 0) then
693
 
         print *, ' '
694
 
         write(6,*)'> Checking ga_fill_patch ... '
695
 
         call ffflush(6)
696
 
      endif
697
 
c
698
 
      val = (1d0,-1d0)
699
 
      call ga_fill_patch(g_a, 2,n/2, 2,n, val)
700
 
*     call ga_print(g_a, 1)
701
 
      do j = 2+me, n, nproc 
702
 
         call ga_get(g_a, 1,n/2, j,j, buf,n/2) 
703
 
         do i = 2, n/2
704
 
            if(buf(i) .ne.val ) then
705
 
               print *,me, ' error ',i,j, buf(i),val 
706
 
               call ga_error('exiting ...',0)
707
 
            endif
708
 
         enddo
709
 
      enddo
710
 
c
711
 
      call ga_sync()
712
 
      if (me .eq. 0) then
713
 
         write(6,*)'   OK '
714
 
         call ffflush(6)
715
 
      endif
716
 
c  
717
 
      do j = 1+me, n, nproc 
718
 
         call ga_put(g_a,1,n,j,j,a(1,j),n) 
719
 
      enddo
720
 
c***
721
 
      if (me .eq. 0) then
722
 
         print *, ' '
723
 
         write(6,*)'> Checking ga_copy_patch ... '
724
 
         call ffflush(6)
725
 
      endif
726
 
c
727
 
      do loop =1, 10
728
 
              ailo = iran(n/2)
729
 
              ajlo = iran(n/2)
730
 
              aihi = min(n, -1+ailo+n/2)
731
 
              ajhi = min(n, -1+ajlo+n/4)
732
 
              rows = aihi -ailo+1
733
 
              cols = ajhi -ajlo +1
734
 
c
735
 
              bilo = iran(n/3)
736
 
              bjlo = iran(n/3)
737
 
              bihi = bilo + rows -1
738
 
              bjhi = bjlo + cols -1
739
 
              if (me .eq. 0) then
740
 
                write(6,'(2x,1h[,4i4,1h],5h --> ,1h[,4i4,1h])')
741
 
     $                ailo,aihi,ajlo,ajhi,
742
 
     $                bilo,bihi, bjlo, bjhi
743
 
 
744
 
                call ffflush(6)
745
 
              endif
746
 
 
747
 
c
748
 
              call ga_copy_patch('n', g_a, ailo, aihi, ajlo, ajhi,
749
 
     &                                g_b, bilo, bihi, bjlo, bjhi)
750
 
c             call ga_print(g_a,1)
751
 
c             call ga_print(g_b,1)
752
 
              call ga_get(g_b,bilo,bihi, bjlo, bjhi, buf, rows)
753
 
              base = 0
754
 
              do j = ajlo, ajhi
755
 
                 if(Mod(j,nproc).eq.me) then
756
 
                    do i = ailo, aihi
757
 
                       base = base+1
758
 
                       if(buf(base) .ne. a(i,j)) then
759
 
                          print *,me, ' error ',i,j, buf(base), a(i,j)
760
 
                          call ga_error('exiting ...',0)
761
 
                       endif
762
 
                    enddo
763
 
                 else
764
 
                    base = base + rows
765
 
                 endif
766
 
              enddo
767
 
      enddo
768
 
c
769
 
      ailo = iran(n/2)
770
 
      ajlo = iran(n/2)
771
 
      bilo = iran(n/2)
772
 
      bjlo = iran(n/2)
773
 
c
774
 
#ifdef BLOCK_CYCLIC
775
 
      if (me .eq. 0) then
776
 
         write(6,*)'  without transpose:   OK '
777
 
         call ffflush(6)
778
 
      endif
779
 
#else
780
 
      if (me .eq. 0) then
781
 
                write(6,'(2x,1h[,4i4,1h],5h --> ,1h[,4i4,1h])')
782
 
     $                ailo,aihi,ajlo,ajhi,
783
 
     $                bilo,bihi, bjlo, bjhi
784
 
      endif
785
 
 
786
 
      call ga_copy_patch('n', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
787
 
     $                        g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
788
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,1)
789
 
*     call ga_print(g_b, 1)
790
 
      call ga_get(g_b,bilo,bilo+n/3,bjlo, bjlo+n/2,buf,n/3+1) 
791
 
      base = 0
792
 
      do j = ajlo, ajlo+n/3 
793
 
         if(Mod(j,nproc).eq.me) then
794
 
            do i = ailo, ailo+n/2 
795
 
               base = base+1
796
 
               if(buf(base) .ne. a(i,j)) then
797
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)
798
 
                  call ga_error('exiting ...',0)
799
 
               endif
800
 
            enddo
801
 
         else
802
 
            base = base +n/2+1
803
 
         endif
804
 
      enddo
805
 
      call ga_sync()
806
 
c
807
 
      if (me .eq. 0) then
808
 
         write(6,*)'  without transpose:   OK '
809
 
         call ffflush(6)
810
 
      endif
811
 
812
 
      if (me .eq. 0) then
813
 
                write(6,'(2x,1h[,4i4,2h]~,5h --> ,1h[,4i4,1h])')
814
 
     $                ailo,aihi,ajlo,ajhi,
815
 
     $                bilo,bihi, bjlo, bjhi
816
 
      endif
817
 
      call ga_copy_patch('t', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
818
 
     $                        g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
819
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,1)
820
 
*     call ga_print(g_b, 1)
821
 
      call ga_get(g_b,bilo,bilo+n/3,bjlo, bjlo+n/2,buf,n/3+1)
822
 
      base = 0
823
 
      do i = ailo, ailo+n/2
824
 
         if(Mod(j,nproc).eq.me) then
825
 
         do j = ajlo, ajlo+n/3
826
 
               base = base+1
827
 
               if(buf(base) .ne. a(i,j)) then
828
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)
829
 
                  call ga_error('exiting ...',0)
830
 
               endif
831
 
            enddo
832
 
         else
833
 
            base = base +n/3+1
834
 
         endif
835
 
      enddo
836
 
      call ga_sync()
837
 
c
838
 
      if (me .eq. 0) then
839
 
         write(6,*)'  transposed:   OK '
840
 
         call ffflush(6)
841
 
      endif
842
 
#endif
843
 
c
844
 
c***
845
 
      if (me .eq. 0) then
846
 
         print *, ' '
847
 
         write(6,*)'> Checking ga_scale_patch ... '
848
 
         call ffflush(6)
849
 
      endif
850
 
      call ga_copy_patch('n', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
851
 
     $                        g_b, bilo,bilo+n/2, bjlo, bjlo+n/3)
852
 
      val = 1.d0
853
 
      call ga_scale_patch(g_b, bilo,bilo+n/2,bjlo, bjlo+n/3,val)
854
 
      call ga_get(g_b,bilo,bilo+n/2,bjlo, bjlo+n/3,buf,n/2+1)
855
 
      base = 0
856
 
      do j = ajlo, ajlo+n/3
857
 
         if(Mod(j,nproc).eq.me) then
858
 
            do i = ailo, ailo+n/2
859
 
               base = base+1
860
 
               if(buf(base) .ne. a(i,j)*val) then
861
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)*val
862
 
                  call ga_error('exiting ...',0)
863
 
               endif
864
 
            enddo
865
 
         else
866
 
            base = base +n/2+1
867
 
         endif
868
 
      enddo
869
 
      call ga_sync()
870
 
c
871
 
      if (me .eq. 0) then
872
 
         write(6,*)'  OK '
873
 
         call ffflush(6)
874
 
      endif
875
 
c
876
 
c***
877
 
      if (me .eq. 0) then
878
 
         print *, ' '
879
 
         write(6,*)'> Checking ga_add_patch ... '
880
 
         call ffflush(6)
881
 
      endif
882
 
      alpha = .1d0
883
 
      beta = .2d0
884
 
      call ga_zero(g_c)
885
 
      call ga_add_patch(alpha, g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
886
 
     $                   beta, g_b, bilo,bilo+n/2, bjlo, bjlo+n/3,
887
 
     $                         g_c, bilo,bilo+n/2, bjlo, bjlo+n/3)
888
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/3,buf,n/2+1)
889
 
      base = 0
890
 
      val = val*beta + alpha
891
 
      do j = ajlo, ajlo+n/3
892
 
         if(Mod(j,nproc).eq.me) then
893
 
            do i = ailo, ailo+n/2
894
 
               base = base+1
895
 
               if(ABS(buf(base)- a(i,j)*val).gt.1d-5) then
896
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)*val
897
 
                  call ga_error('exiting ...',0)
898
 
               endif
899
 
            enddo
900
 
         else
901
 
            base = base +n/2+1
902
 
         endif
903
 
      enddo
904
 
      call ga_sync()
905
 
c
906
 
      if (me .eq. 0) then
907
 
         write(6,*)'  OK '
908
 
         call ffflush(6)
909
 
      endif
910
 
c
911
 
c***
912
 
      if (me .eq. 0) then
913
 
         print *, ' '
914
 
         write(6,*)'> Checking ga_zdot_patch ... '
915
 
         call ffflush(6)
916
 
      endif
917
 
      alpha= ga_zdot_patch(g_a,'n', ailo,ailo+n/2, ajlo, ajlo+n/3,
918
 
     $                     g_c,'n', bilo,bilo+n/2, bjlo, bjlo+n/3)
919
 
      beta = (0d0,0d0)
920
 
      do j = ajlo, ajlo+n/3
921
 
            do i = ailo, ailo+n/2
922
 
               beta = beta + a(i,j)*a(i,j)
923
 
            enddo
924
 
      enddo
925
 
      if(ABS(beta*val- alpha)/abs(alpha).gt.1d-6) then
926
 
             print *,me, ' error ', beta*val, alpha
927
 
             call ga_error('exiting ...',0)
928
 
      endif
929
 
      call ga_sync()
930
 
c
931
 
      if (me .eq. 0) then
932
 
         write(6,*)'  OK '
933
 
         call ffflush(6)
934
 
      endif
935
 
c
936
 
c......................................................
937
 
      if (me .eq. 0) then
938
 
         print *, ' '
939
 
         write(6,*)'> Checking ga_matmul_patch ... '
940
 
         call ffflush(6)
941
 
      endif
942
 
      do j = 1+me, n, nproc
943
 
         call ga_put(g_b,1,n,j,j,b(1,j),n)
944
 
      enddo
945
 
      call ga_sync()
946
 
      call ga_matmul_patch('n','n', (1d0,0d0), (0d0,0d0),
947
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
948
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
949
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
950
 
      call ygemm('n','n',n/2+1,n/2+1,n/3+1,(1d0,0d0),
951
 
     $            a(ailo,ajlo), n,
952
 
     $            b(bilo,bjlo),n, (0d0,0d0), c, n)
953
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3)
954
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
955
 
*     call ga_print(g_c, 1)
956
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
957
 
      base = 0
958
 
      do j = 1, 1+n/2
959
 
         if(Mod(j,nproc).eq.me) then
960
 
            do i = 1, 1+n/2
961
 
               base = base+1
962
 
               if(ABS(buf(base)- c(i,j))/abs(c(i,j)).gt.1d-8) then
963
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
964
 
                  call ga_error('exiting ...',0)
965
 
               endif
966
 
            enddo
967
 
         else
968
 
            base = base +n/2+1
969
 
         endif
970
 
      enddo
971
 
c
972
 
      call ga_sync()
973
 
      if (me .eq. 0) then
974
 
         write(6,*)'  a*b: OK '
975
 
         call ffflush(6)
976
 
      endif
977
 
c
978
 
      call ga_sync()
979
 
      call ga_matmul_patch('t','n', (1d0,0d0), (0d0,0d0),
980
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
981
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
982
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
983
 
      call ygemm('t','n',n/2+1,n/2+1,n/3+1,(1d0,0d0),
984
 
     $            a(ajlo,ailo), n,
985
 
     $            b(bilo,bjlo),n, (0d0,0d0), c, n)
986
 
*     call ga_print(g_a,1) 
987
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
988
 
*     call ga_print(g_c, 1)
989
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
990
 
      base = 0
991
 
      do j = 1, 1+n/2
992
 
         if(Mod(j,nproc).eq.me) then
993
 
            do i = 1, 1+n/2
994
 
               base = base+1
995
 
               if(ABS(buf(base)- c(i,j))/abs(c(i,j)).gt.1d-8) then
996
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
997
 
                  call ga_error('exiting ...',0)
998
 
               endif
999
 
            enddo
1000
 
         else
1001
 
            base = base +n/2+1
1002
 
         endif
1003
 
      enddo
1004
 
c
1005
 
      call ga_sync()
1006
 
      if (me .eq. 0) then
1007
 
         write(6,*)'  trans(a)*b: OK '
1008
 
         call ffflush(6)
1009
 
      endif
1010
 
c
1011
 
      call ga_sync()
1012
 
      call ga_matmul_patch('n','t', (1d0,0d0), (0d0,0d0),
1013
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1014
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
1015
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
1016
 
      call ygemm('n','t',n/2+1,n/2+1,n/3+1,(1d0,0d0),
1017
 
     $            a(ailo,ajlo), n,
1018
 
     $            b(bjlo,bilo),n, (0d0,0d0), c, n)
1019
 
*     call ga_print(g_a,1)
1020
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1021
 
*     call ga_print(g_c, 1)
1022
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
1023
 
      base = 0
1024
 
      do j = 1, 1+n/2
1025
 
         if(Mod(j,nproc).eq.me) then
1026
 
            do i = 1, 1+n/2
1027
 
               base = base+1
1028
 
               if(ABS(buf(base)- c(i,j))/abs(c(i,j)).gt.1d-8) then
1029
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
1030
 
                  call ga_error('exiting ...',0)
1031
 
               endif
1032
 
            enddo
1033
 
         else
1034
 
            base = base +n/2+1
1035
 
         endif
1036
 
      enddo
1037
 
c
1038
 
      call ga_sync()
1039
 
      if (me .eq. 0) then
1040
 
         write(6,*)'  a*trans(b): OK '
1041
 
         call ffflush(6)
1042
 
      endif
1043
 
c
1044
 
      call ga_sync()
1045
 
      call ga_matmul_patch('t','t', (1d0,0d0), (0d0,0d0),
1046
 
     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1047
 
     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
1048
 
     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
1049
 
      call ygemm('t','t',n/2+1,n/2+1,n/3+1,(1d0,0d0),
1050
 
     $            a(ajlo,ailo), n,
1051
 
     $            b(bilo,bjlo),n, (0d0,0d0), c, n)
1052
 
*     call ga_print(g_a,1)
1053
 
*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1054
 
*     call ga_print(g_c, 1)
1055
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
1056
 
      base = 0
1057
 
      do j = 1, 1+n/2
1058
 
         if(Mod(j,nproc).eq.me) then
1059
 
            do i = 1, 1+n/2
1060
 
               base = base+1
1061
 
               if(ABS(buf(base)- c(i,j))/abs(c(i,j)).gt.1d-8) then
1062
 
                  print *,me, ' error ',i,j, buf(base), c(i,j)
1063
 
                  call ga_error('exiting ...',0)
1064
 
               endif
1065
 
            enddo
1066
 
         else
1067
 
            base = base +n/2+1
1068
 
         endif
1069
 
      enddo
1070
 
c
1071
 
      call ga_sync()
1072
 
      if (me .eq. 0) then
1073
 
         write(6,*)'  trans(a)*trans(b): OK '
1074
 
         call ffflush(6)
1075
 
      endif
1076
 
      status = ga_destroy(g_a)
1077
 
      status = status .and. ga_destroy(g_b)
1078
 
      status = status .and. ga_destroy(g_c)
1079
 
      if(.not. status) print *, 'ga_destroy failed'
1080
 
c
1081
 
      end
1082
 
 
1083
 
 
1084
 
 
1085
 
      subroutine spatch_test()
1086
 
      implicit none
1087
 
#include "mafdecls.fh"
1088
 
#include "global.fh"
1089
 
#include "testutil.fh"
1090
 
c
1091
 
      integer n,m
1092
 
      parameter (n = 128)
1093
 
      parameter (m = n*n)
1094
 
      real a(n,n), b(n,n), c(n,n), buf(m), val
1095
 
      real alpha, beta
1096
 
      integer nproc, me
1097
 
      integer i, j, ailo, ajlo, bilo, bjlo, base, iran
1098
 
      integer aihi, ajhi, bihi, bjhi
1099
 
      integer g_a, g_b, g_c
1100
 
      integer rows, cols, loop
1101
 
      logical status, dist_same
1102
 
#ifdef BLOCK_CYCLIC
1103
 
      integer ndim, dims(2)
1104
 
      integer block_size(2), proc_grid(2)
1105
 
#endif
1106
 
      iran(i) = int(drand(1)*real(i)) + 1
1107
 
      dist_same = .false.
1108
 
c
1109
 
      me = ga_nodeid()
1110
 
      nproc = ga_nnodes()
1111
 
#ifdef BLOCK_CYCLIC
1112
 
      block_size(1) = 32
1113
 
      block_size(2) = 32
1114
 
#ifdef USE_SCALAPACK_DISTR
1115
 
      if (mod(nproc,2).ne.0)
1116
 
     +  call ga_error("Available procs must be divisible by 2",0)
1117
 
      proc_grid(1) = 2
1118
 
      proc_grid(2) = nproc/2
1119
 
#endif
1120
 
#endif
1121
 
c
1122
 
      do j = 1, n
1123
 
         do i = 1, n
1124
 
            a(i,j) = i-1 + (j-1)*n
1125
 
            b(i,j) = i+j
1126
 
         enddo
1127
 
      enddo
1128
 
c
1129
 
c***  Create a global array
1130
 
c
1131
 
#ifndef BLOCK_CYCLIC
1132
 
      status = ga_create(MT_REAL, n, n, 'a', 0, 0, g_a)
1133
 
#else
1134
 
      g_a = ga_create_handle()
1135
 
      ndim = 2
1136
 
      dims(1) = n
1137
 
      dims(2) = n
1138
 
      call ga_set_data(g_a,ndim,dims,MT_REAL)
1139
 
      call ga_set_array_name(g_a,'a')
1140
 
#ifdef USE_SCALAPACK_DISTR
1141
 
      call ga_set_block_cyclic_proc_grid(g_a,block_size,proc_grid)
1142
 
#else
1143
 
      call ga_set_block_cyclic(g_a,block_size)
1144
 
#endif
1145
 
      status = ga_allocate(g_a)
1146
 
#endif
1147
 
      if (.not. status) then
1148
 
         write(6,*) ' ga_create failed'
1149
 
         call ffflush(6)
1150
 
         call ga_error('... exiting ',0)
1151
 
      endif
1152
 
c
1153
 
      if(dist_same) then
1154
 
         status = ga_duplicate(g_a, g_b, 'a_duplicated')
1155
 
         if(.not.ga_compare_distr(g_a, g_b))
1156
 
     $           call ga_error("g_b distribution different",0)
1157
 
         status = ga_duplicate(g_a, g_c, 'a_duplicated_again')
1158
 
         if(.not.ga_compare_distr(g_a, g_c))
1159
 
     $           call ga_error("g_c distribution different",0)
1160
 
      else
1161
 
#ifndef BLOCK_CYCLIC
1162
 
         status = ga_create(MT_REAL, n, n, 'b', 0, n, g_b)
1163
 
#else
1164
 
         g_b = ga_create_handle()
1165
 
         ndim = 2
1166
 
         dims(1) = n
1167
 
         dims(2) = n
1168
 
         call ga_set_data(g_b,ndim,dims,MT_REAL)
1169
 
         call ga_set_array_name(g_b,'b')
1170
 
#ifdef USE_SCALAPACK_DISTR
1171
 
         call ga_set_block_cyclic_proc_grid(g_b,block_size,proc_grid)
1172
 
#else
1173
 
         call ga_set_block_cyclic(g_b,block_size)
1174
 
#endif
1175
 
         status = ga_allocate(g_b)
1176
 
#endif
1177
 
         if (.not. status) call ga_error('ga_create failed:b',0)
1178
 
#ifndef BLOCK_CYCLIC
1179
 
         status = ga_create(MT_REAL, n, n, 'c', n, 0, g_c)
1180
 
#else
1181
 
         g_c = ga_create_handle()
1182
 
         ndim = 2
1183
 
         dims(1) = n
1184
 
         dims(2) = n
1185
 
         call ga_set_data(g_c,ndim,dims,MT_REAL)
1186
 
         call ga_set_array_name(g_c,'c')
1187
 
#ifdef USE_SCALAPACK_DISTR
1188
 
         call ga_set_block_cyclic_proc_grid(g_c,block_size,proc_grid)
1189
 
#else
1190
 
         call ga_set_block_cyclic(g_c,block_size)
1191
 
#endif
1192
 
         status = ga_allocate(g_c)
1193
 
#endif
1194
 
         if (.not. status) call ga_error('ga_create failed:c',0)
1195
 
      endif
1196
 
c
1197
 
c***
1198
 
      if (me .eq. 0) then
1199
 
         print *, ' '
1200
 
         write(6,*)'> Checking ga_fill_patch ... '
1201
 
         call ffflush(6)
1202
 
      endif
1203
 
c
1204
 
      val = 1.0
1205
 
      call ga_fill_patch(g_a, 2,n/2, 2,n, val)
1206
 
*     call ga_print(g_a, 1)
1207
 
      do j = 2+me, n, nproc
1208
 
         call ga_get(g_a, 1,n/2, j,j, buf,n/2)
1209
 
         do i = 2, n/2
1210
 
            if(buf(i) .ne.val ) then
1211
 
               print *,me, ' error ',i,j, buf(i),val
1212
 
               call ga_error('exiting ...',0)
1213
 
            endif
1214
 
         enddo
1215
 
      enddo
1216
 
c
1217
 
      call ga_sync()
1218
 
      if (me .eq. 0) then
1219
 
         write(6,*)'   OK '
1220
 
         call ffflush(6)
1221
 
      endif
1222
 
1223
 
      do j = 1+me, n, nproc
1224
 
         call ga_put(g_a,1,n,j,j,a(1,j),n)
1225
 
      enddo
1226
 
c***
1227
 
      if (me .eq. 0) then
1228
 
         print *, ' '
1229
 
         write(6,*)'> Checking ga_copy_patch ... '
1230
 
         call ffflush(6)
1231
 
      endif
1232
 
c
1233
 
      do loop =1, 10
1234
 
              ailo = iran(n/2)
1235
 
              ajlo = iran(n/2)
1236
 
              aihi = min(n, -1+ailo+n/2)
1237
 
              ajhi = min(n, -1+ajlo+n/4)
1238
 
              rows = aihi -ailo+1
1239
 
              cols = ajhi -ajlo +1
1240
 
c
1241
 
              bilo = iran(n/3)
1242
 
              bjlo = iran(n/3)
1243
 
              bihi = bilo + rows -1
1244
 
              bjhi = bjlo + cols -1
1245
 
              if (me .eq. 0) then
1246
 
                write(6,'(2x,1h[,4i4,1h],5h --> ,1h[,4i4,1h])')
1247
 
     $                ailo,aihi,ajlo,ajhi,
1248
 
     $                bilo,bihi, bjlo, bjhi
1249
 
 
1250
 
                call ffflush(6)
1251
 
              endif
1252
 
 
1253
 
c
1254
 
              call ga_copy_patch('n', g_a, ailo, aihi, ajlo, ajhi,
1255
 
     &                                g_b, bilo, bihi, bjlo, bjhi)
1256
 
c             call ga_print(g_a,1)
1257
 
c             call ga_print(g_b,1)
1258
 
              call ga_get(g_b,bilo,bihi, bjlo, bjhi, buf, rows)
1259
 
              base = 0
1260
 
              do j = ajlo, ajhi
1261
 
                 if(Mod(j,nproc).eq.me) then
1262
 
                    do i = ailo, aihi
1263
 
                       base = base+1
1264
 
                       if(buf(base) .ne. a(i,j)) then
1265
 
                          print *,me, ' error ',i,j, buf(base), a(i,j)
1266
 
                          call ga_error('exiting ...',0)
1267
 
                       endif
1268
 
                    enddo
1269
 
                 else
1270
 
                    base = base + rows
1271
 
                 endif
1272
 
              enddo
1273
 
      enddo
1274
 
c
1275
 
      ailo = iran(n/2)
1276
 
      ajlo = iran(n/2)
1277
 
      bilo = iran(n/2)
1278
 
      bjlo = iran(n/2)
1279
 
c
1280
 
#ifdef BLOCK_CYCLIC
1281
 
      if (me .eq. 0) then
1282
 
         write(6,*)'  without transpose:   OK '
1283
 
         call ffflush(6)
1284
 
      endif
1285
 
#else
1286
 
      if (me .eq. 0) then
1287
 
                write(6,'(2x,1h[,4i4,1h],5h --> ,1h[,4i4,1h])')
1288
 
     $                ailo,aihi,ajlo,ajhi,
1289
 
     $                bilo,bihi, bjlo, bjhi
1290
 
      endif
1291
 
 
1292
 
      call ga_copy_patch('n', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1293
 
     $                        g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1294
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,1)
1295
 
*     call ga_print(g_b, 1)
1296
 
      call ga_get(g_b,bilo,bilo+n/3,bjlo, bjlo+n/2,buf,n/3+1)
1297
 
      base = 0
1298
 
      do j = ajlo, ajlo+n/3
1299
 
         if(Mod(j,nproc).eq.me) then
1300
 
            do i = ailo, ailo+n/2
1301
 
               base = base+1
1302
 
               if(buf(base) .ne. a(i,j)) then
1303
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)
1304
 
                  call ga_error('exiting ...',0)
1305
 
               endif
1306
 
            enddo
1307
 
         else
1308
 
            base = base +n/2+1
1309
 
         endif
1310
 
      enddo
1311
 
      call ga_sync()
1312
 
c
1313
 
      if (me .eq. 0) then
1314
 
         write(6,*)'  without transpose:   OK '
1315
 
         call ffflush(6)
1316
 
      endif
1317
 
c
1318
 
      if (me .eq. 0) then
1319
 
                write(6,'(2x,1h[,4i4,2h]~,5h --> ,1h[,4i4,1h])')
1320
 
     $                ailo,aihi,ajlo,ajhi,
1321
 
     $                bilo,bihi, bjlo, bjhi
1322
 
      endif
1323
 
      call ga_copy_patch('t', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1324
 
     $                        g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1325
 
*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,1)
1326
 
*     call ga_print(g_b, 1)
1327
 
      call ga_get(g_b,bilo,bilo+n/3,bjlo, bjlo+n/2,buf,n/3+1)
1328
 
      base = 0
1329
 
      do i = ailo, ailo+n/2
1330
 
         if(Mod(j,nproc).eq.me) then
1331
 
         do j = ajlo, ajlo+n/3
1332
 
               base = base+1
1333
 
               if(buf(base) .ne. a(i,j)) then
1334
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)
1335
 
                  call ga_error('exiting ...',0)
1336
 
               endif
1337
 
            enddo
1338
 
         else
1339
 
            base = base +n/3+1
1340
 
         endif
1341
 
      enddo
1342
 
      call ga_sync()
1343
 
c
1344
 
      if (me .eq. 0) then
1345
 
         write(6,*)'  transposed:   OK '
1346
 
         call ffflush(6)
1347
 
      endif
1348
 
#endif
1349
 
c
1350
 
c***
1351
 
      if (me .eq. 0) then
1352
 
         print *, ' '
1353
 
         write(6,*)'> Checking ga_scale_patch ... '
1354
 
         call ffflush(6)
1355
 
      endif
1356
 
 
1357
 
      call ga_copy_patch('n', g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1358
 
     $                        g_b, bilo,bilo+n/2, bjlo, bjlo+n/3)
1359
 
      val = 1.0
1360
 
      call ga_scale_patch(g_b, bilo,bilo+n/2,bjlo, bjlo+n/3,val)
1361
 
      call ga_get(g_b,bilo,bilo+n/2,bjlo, bjlo+n/3,buf,n/2+1)
1362
 
      base = 0
1363
 
      do j = ajlo, ajlo+n/3
1364
 
         if(Mod(j,nproc).eq.me) then
1365
 
            do i = ailo, ailo+n/2
1366
 
               base = base+1
1367
 
               if(buf(base) .ne. a(i,j)*val) then
1368
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)*val
1369
 
                  call ga_error('exiting ...',0)
1370
 
               endif
1371
 
            enddo
1372
 
         else
1373
 
            base = base +n/2+1
1374
 
         endif
1375
 
      enddo
1376
 
      call ga_sync()
1377
 
c
1378
 
      if (me .eq. 0) then
1379
 
         write(6,*)'  OK '
1380
 
         call ffflush(6)
1381
 
      endif
1382
 
c
1383
 
c***
1384
 
      if (me .eq. 0) then
1385
 
         print *, ' '
1386
 
         write(6,*)'> Checking ga_add_patch ... '
1387
 
         call ffflush(6)
1388
 
      endif
1389
 
      alpha = 0.1
1390
 
      beta = 0.2
1391
 
      call ga_zero(g_c)
1392
 
 
1393
 
      call ga_add_patch(alpha, g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1394
 
     $                   beta, g_b, bilo,bilo+n/2, bjlo, bjlo+n/3,
1395
 
     $                         g_c, bilo,bilo+n/2, bjlo, bjlo+n/3)
1396
 
      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/3,buf,n/2+1)
1397
 
      base = 0
1398
 
      val = val*beta + alpha
1399
 
      do j = ajlo, ajlo+n/3
1400
 
         if(Mod(j,nproc).eq.me) then
1401
 
            do i = ailo, ailo+n/2
1402
 
               base = base+1
1403
 
               if(ABS(buf(base)- a(i,j)*val).gt.1e-3) then
1404
 
                  print *,me, ' error ',i,j, buf(base), a(i,j)*val
1405
 
                  call ga_error('exiting ...',0)
1406
 
               endif
1407
 
            enddo
1408
 
         else
1409
 
            base = base +n/2+1
1410
 
         endif
1411
 
      enddo
1412
 
      call ga_sync()
1413
 
c
1414
 
      if (me .eq. 0) then
1415
 
         write(6,*)'  OK '
1416
 
         call ffflush(6)
1417
 
      endif
1418
 
c
1419
 
c***
1420
 
      if (me .eq. 0) then
1421
 
         print *, ' '
1422
 
         write(6,*)'> Checking ga_sdot_patch ... '
1423
 
         call ffflush(6)
1424
 
      endif
1425
 
      alpha= ga_sdot_patch(g_a,'n', ailo,ailo+n/2, ajlo, ajlo+n/3,
1426
 
     $                     g_c,'n', bilo,bilo+n/2, bjlo, bjlo+n/3)
1427
 
      beta = 0.0
1428
 
 
1429
 
      do j = ajlo, ajlo+n/3
1430
 
            do i = ailo, ailo+n/2
1431
 
               beta = beta + a(i,j)**2
1432
 
            enddo
1433
 
      enddo
1434
 
      if(ABS(beta*val- alpha).gt.(1.0e-2*alpha)) then
1435
 
             print *,me, ' error ', beta*val, alpha
1436
 
             call ga_error('exiting ...',0)
1437
 
      endif
1438
 
      call ga_sync()
1439
 
c
1440
 
      if (me .eq. 0) then
1441
 
         write(6,*)'  OK '
1442
 
         call ffflush(6)
1443
 
      endif
1444
 
c......................................................
1445
 
ccc      if (me .eq. 0) then
1446
 
ccc         print *, ' '
1447
 
ccc         write(6,*)'> Checking ga_matmul_patch ... '
1448
 
ccc         call ffflush(6)
1449
 
ccc      endif
1450
 
ccc      do j = 1+me, n, nproc
1451
 
ccc         call ga_put(g_b,1,n,j,j,b(1,j),n)
1452
 
ccc      enddo
1453
 
ccc      call ga_sync()
1454
 
ccc      call ga_matmul_patch('n','n', 1.0, 0.0,
1455
 
ccc     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1456
 
ccc     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
1457
 
ccc     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
1458
 
ccc      call xgemm('n','n',n/2+1,n/2+1,n/3+1,1d0,a(ailo,ajlo), n,
1459
 
ccc     $            b(bilo,bjlo),n, 0d0, c, n)
1460
 
ccc*     call ga_print_patch(g_a, ailo,ailo+n/2, ajlo, ajlo+n/3)
1461
 
ccc*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1462
 
ccc*     call ga_print(g_c, 1)
1463
 
ccc      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
1464
 
ccc      base = 0
1465
 
ccc      do j = 1, 1+n/2
1466
 
ccc         if(Mod(j,nproc).eq.me) then
1467
 
ccc            do i = 1, 1+n/2
1468
 
ccc               base = base+1
1469
 
ccc               if(ABS(buf(base)- c(i,j)).gt.1e-8) then
1470
 
ccc                  print *,me, ' error ',i,j, buf(base), c(i,j)
1471
 
ccc                  call ga_error('exiting ...',0)
1472
 
ccc               endif
1473
 
ccc            enddo
1474
 
ccc         else
1475
 
ccc            base = base +n/2+1
1476
 
ccc         endif
1477
 
ccc      enddo
1478
 
cccc
1479
 
ccc      call ga_sync()
1480
 
ccc      if (me .eq. 0) then
1481
 
ccc         write(6,*)'  a*b: OK '
1482
 
ccc         call ffflush(6)
1483
 
ccc      endif
1484
 
cccc
1485
 
ccc      call ga_sync()
1486
 
ccc      call ga_matmul_patch('t','n', 1d0, 0d0,
1487
 
ccc     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1488
 
ccc     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
1489
 
ccc     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
1490
 
ccc      call xgemm('t','n',n/2+1,n/2+1,n/3+1,1.0,a(ajlo,ailo), n,
1491
 
ccc     $            b(bilo,bjlo),n, 0.0, c, n)
1492
 
ccc*     call ga_print(g_a,1)
1493
 
ccc*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1494
 
ccc*     call ga_print(g_c, 1)
1495
 
ccc      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
1496
 
ccc      base = 0
1497
 
ccc      do j = 1, 1+n/2
1498
 
ccc         if(Mod(j,nproc).eq.me) then
1499
 
ccc            do i = 1, 1+n/2
1500
 
ccc               base = base+1
1501
 
ccc               if(ABS(buf(base)- c(i,j)).gt.1e-8) then
1502
 
ccc                  print *,me, ' error ',i,j, buf(base), c(i,j)
1503
 
ccc                  call ga_error('exiting ...',0)
1504
 
ccc               endif
1505
 
ccc            enddo
1506
 
ccc         else
1507
 
ccc            base = base +n/2+1
1508
 
ccc         endif
1509
 
ccc      enddo
1510
 
cccc
1511
 
ccc      call ga_sync()
1512
 
ccc      if (me .eq. 0) then
1513
 
ccc         write(6,*)'  trans(a)*b: OK '
1514
 
ccc         call ffflush(6)
1515
 
ccc      endif
1516
 
cccc
1517
 
ccc      call ga_sync()
1518
 
ccc      call ga_matmul_patch('n','t', 1.0, 0.0,
1519
 
ccc     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1520
 
ccc     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
1521
 
ccc     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
1522
 
ccc      call xgemm('n','t',n/2+1,n/2+1,n/3+1,1.0,a(ailo,ajlo), n,
1523
 
ccc     $            b(bjlo,bilo),n, 0.0, c, n)
1524
 
ccc*     call ga_print(g_a,1)
1525
 
ccc*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1526
 
ccc*     call ga_print(g_c, 1)
1527
 
ccc      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
1528
 
ccc      base = 0
1529
 
ccc      do j = 1, 1+n/2
1530
 
ccc         if(Mod(j,nproc).eq.me) then
1531
 
ccc            do i = 1, 1+n/2
1532
 
ccc               base = base+1
1533
 
ccc               if(ABS(buf(base)- c(i,j)).gt.1e-8) then
1534
 
ccc                  print *,me, ' error ',i,j, buf(base), c(i,j)
1535
 
ccc                  call ga_error('exiting ...',0)
1536
 
ccc               endif
1537
 
ccc            enddo
1538
 
ccc         else
1539
 
ccc            base = base +n/2+1
1540
 
ccc         endif
1541
 
ccc      enddo
1542
 
cccc
1543
 
cccc
1544
 
ccc      call ga_sync()
1545
 
ccc      if (me .eq. 0) then
1546
 
ccc         write(6,*)'  a*trans(b): OK '
1547
 
ccc         call ffflush(6)
1548
 
ccc      endif
1549
 
cccc
1550
 
ccc      call ga_sync()
1551
 
ccc      call ga_matmul_patch('t','t', 1d0, 0d0,
1552
 
ccc     $                      g_a, ailo,ailo+n/2, ajlo, ajlo+n/3,
1553
 
ccc     $                      g_b, bilo,bilo+n/3, bjlo, bjlo+n/2,
1554
 
ccc     $                      g_c, bilo,bilo+n/2, bjlo, bjlo+n/2)
1555
 
ccc      call xgemm('t','t',n/2+1,n/2+1,n/3+1,1.0,a(ajlo,ailo), n,
1556
 
ccc     $            b(bjlo,bilo),n, 0.0, c, n)
1557
 
ccc*     call ga_print(g_a,1)
1558
 
ccc*     call ga_print_patch(g_b, bilo,bilo+n/3, bjlo, bjlo+n/2)
1559
 
ccc*     call ga_print(g_c, 1)
1560
 
ccc      call ga_get(g_c,bilo,bilo+n/2,bjlo, bjlo+n/2,buf,n/2+1)
1561
 
ccc      base = 0
1562
 
ccc      do j = 1, 1+n/2
1563
 
ccc         if(Mod(j,nproc).eq.me) then
1564
 
ccc            do i = 1, 1+n/2
1565
 
ccc               base = base+1
1566
 
ccc               if(ABS(buf(base)- c(i,j)).gt.1e-8) then
1567
 
ccc                  print *,me, ' error ',i,j, buf(base), c(i,j)
1568
 
ccc                  call ga_error('exiting ...',0)
1569
 
ccc               endif
1570
 
ccc            enddo
1571
 
ccc         else
1572
 
ccc            base = base +n/2+1
1573
 
ccc         endif
1574
 
ccc      enddo
1575
 
cccc
1576
 
ccc      call ga_sync()
1577
 
ccc      if (me .eq. 0) then
1578
 
ccc         write(6,*)'  trans(a)*trans(b): OK '
1579
 
ccc         call ffflush(6)
1580
 
ccc      endif
1581
 
ccc      status = ga_destroy(g_a)
1582
 
ccc      status = status .and. ga_destroy(g_b)
1583
 
ccc      status = status .and. ga_destroy(g_c)
1584
 
ccc      if(.not. status) print *, 'ga_destroy failed'
1585
 
cccc
1586
 
      end
1587
 
 
1588