~ubuntu-branches/ubuntu/saucy/nwchem/saucy

« back to all changes in this revision

Viewing changes to src/tools/ga-5-1/global/src/global.nalg.c

  • Committer: Package Import Robot
  • Author(s): Michael Banck, Michael Banck, Daniel Leidert
  • Date: 2012-02-09 20:02:41 UTC
  • mfrom: (1.1.1)
  • Revision ID: package-import@ubuntu.com-20120209200241-jgk03qfsphal4ug2
Tags: 6.1-1
* New upstream release.

[ Michael Banck ]
* debian/patches/02_makefile_flags.patch: Updated.
* debian/patches/02_makefile_flags.patch: Use internal blas and lapack code.
* debian/patches/02_makefile_flags.patch: Define GCC4 for LINUX and LINUX64
  (Closes: #632611 and LP: #791308).
* debian/control (Build-Depends): Added openssh-client.
* debian/rules (USE_SCALAPACK, SCALAPACK): Removed variables (Closes:
  #654658).
* debian/rules (LIBDIR, USE_MPIF4, ARMCI_NETWORK): New variables.
* debian/TODO: New file.
* debian/control (Build-Depends): Removed libblas-dev, liblapack-dev and
  libscalapack-mpi-dev.
* debian/patches/04_show_testsuite_diff_output.patch: New patch, shows the
  diff output for failed tests.
* debian/patches/series: Adjusted.
* debian/testsuite: Optionally run all tests if "all" is passed as option.
* debian/rules: Run debian/testsuite with "all" if DEB_BUILD_OPTIONS
  contains "checkall".

[ Daniel Leidert ]
* debian/control: Used wrap-and-sort. Added Vcs-Svn and Vcs-Browser fields.
  (Priority): Moved to extra according to policy section 2.5.
  (Standards-Version): Bumped to 3.9.2.
  (Description): Fixed a typo.
* debian/watch: Added.
* debian/patches/03_hurd-i386_define_path_max.patch: Added.
  - Define MAX_PATH if not defines to fix FTBFS on hurd.
* debian/patches/series: Adjusted.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#if HAVE_CONFIG_H
 
2
#   include "config.h"
 
3
#endif
 
4
 
 
5
/*************************************************************************\
 
6
 Purpose:   File global.nalg.c contains a set of linear algebra routines 
 
7
            that operate on n-dim global arrays in the SPMD mode. 
 
8
 
 
9
 Date: 10.22.98
 
10
 Author: Jarek Nieplocha
 
11
\************************************************************************/
 
12
 
 
13
 
 
14
#if HAVE_STDIO_H
 
15
#   include <stdio.h>
 
16
#endif
 
17
#if HAVE_STRING_H
 
18
#   include <string.h>
 
19
#endif
 
20
#include "message.h"
 
21
#include "globalp.h"
 
22
#include "armci.h"
 
23
#include "ga-papi.h"
 
24
#include "ga-wapi.h"
 
25
 
 
26
#ifdef MPI
 
27
extern ARMCI_Group* ga_get_armci_group_(int);
 
28
#endif
 
29
 
 
30
/* work arrays used in all routines */
 
31
static Integer dims[MAXDIM], ld[MAXDIM-1];
 
32
static Integer lo[MAXDIM],hi[MAXDIM];
 
33
static Integer one_arr[MAXDIM]={1,1,1,1,1,1,1};
 
34
 
 
35
#define GET_ELEMS(ndim,lo,hi,ld,pelems){\
 
36
int _i;\
 
37
      for(_i=0, *pelems = hi[ndim-1]-lo[ndim-1]+1; _i< ndim-1;_i++) {\
 
38
         if(ld[_i] != (hi[_i]-lo[_i]+1)) pnga_error("layout problem",_i);\
 
39
         *pelems *= hi[_i]-lo[_i]+1;\
 
40
      }\
 
41
}
 
42
 
 
43
#define GET_ELEMS_W_GHOSTS(ndim,lo,hi,ld,pelems){\
 
44
int _i;\
 
45
      for(_i=0, *pelems = hi[ndim-1]-lo[ndim-1]+1; _i< ndim-1;_i++) {\
 
46
         if(ld[_i] < (hi[_i]-lo[_i]+1))\
 
47
           pnga_error("layout problem with ghosts",_i);\
 
48
         *pelems *= hi[_i]-lo[_i]+1;\
 
49
      }\
 
50
}
 
51
 
 
52
 
 
53
#if HAVE_SYS_WEAK_ALIAS_PRAGMA
 
54
#   pragma weak wnga_zero = pnga_zero
 
55
#endif
 
56
void pnga_zero(Integer g_a)
 
57
{
 
58
  Integer ndim, type, me, elems, p_handle;
 
59
  Integer num_blocks;
 
60
  void *ptr;
 
61
  /*register Integer i;*/
 
62
  int local_sync_begin,local_sync_end;
 
63
 
 
64
  local_sync_begin = _ga_sync_begin; local_sync_end = _ga_sync_end;
 
65
  _ga_sync_begin = 1; _ga_sync_end=1; /*remove any previous masking*/
 
66
  p_handle = pnga_get_pgroup(g_a);
 
67
 
 
68
  if(local_sync_begin) pnga_pgroup_sync(p_handle);
 
69
 
 
70
  me = pnga_pgroup_nodeid(p_handle);
 
71
 
 
72
  pnga_check_handle(g_a, "ga_zero");
 
73
  GA_PUSH_NAME("ga_zero");
 
74
 
 
75
  num_blocks = pnga_total_blocks(g_a);
 
76
 
 
77
  pnga_inquire(g_a, &type, &ndim, dims);
 
78
  if (num_blocks < 0) {
 
79
    pnga_distribution(g_a, me, lo, hi);
 
80
 
 
81
    if ( lo[0]> 0 ){ /* base index is 1: we get 0 if no elements stored on p */
 
82
 
 
83
      if (pnga_has_ghosts(g_a)) {
 
84
        pnga_zero_patch(g_a,lo,hi);
 
85
        return;
 
86
      }
 
87
      pnga_access_ptr(g_a, lo, hi, &ptr, ld);
 
88
      GET_ELEMS(ndim,lo,hi,ld,&elems);
 
89
 
 
90
      /* switch (type){ */
 
91
/*         int *ia; */
 
92
/*         double *da; */
 
93
/*         float *fa; */
 
94
/*         long *la; */
 
95
/*         long long *lla; */
 
96
/*         case C_INT: */
 
97
/*         ia = (int*)ptr; */
 
98
/*         for(i=0;i<elems;i++) ia[i]  = 0; */
 
99
/*         break; */
 
100
/*         case C_DCPL: */
 
101
/*         elems *=2; */
 
102
/*         case C_DBL: */
 
103
/*         da = (double*)ptr; */
 
104
/*         for(i=0;i<elems;i++) da[i] = 0; */
 
105
/*         break; */
 
106
/*         case C_SCPL: */
 
107
/*         elems *=2; */
 
108
/*         case C_FLOAT: */
 
109
/*         fa = (float*)ptr; */
 
110
/*         for(i=0;i<elems;i++) fa[i]  = 0; */
 
111
/*         break; */
 
112
/*         case C_LONG: */
 
113
/*         la = (long*)ptr; */
 
114
/*         for(i=0;i<elems;i++) la[i]  = 0; */
 
115
/*         break; */
 
116
/*         case C_LONGLONG: */
 
117
/*         lla = (long long*)ptr; */
 
118
/*         for(i=0;i<elems;i++) lla[i]  = 0; */
 
119
/*         break; */
 
120
/*         default: pnga_error(" wrong data type ",type); */
 
121
/*       } */
 
122
      memset(ptr, 0, GAsizeofM(type)*elems);
 
123
 
 
124
      /* release access to the data */
 
125
      pnga_release_update(g_a, lo, hi);
 
126
    } 
 
127
  } else {
 
128
    pnga_access_block_segment_ptr(g_a, me, &ptr, &elems);
 
129
/*     switch (type){ */
 
130
/*       int *ia; */
 
131
/*       double *da; */
 
132
/*       float *fa; */
 
133
/*       long *la; */
 
134
/*       long long *lla; */
 
135
/*       case C_INT: */
 
136
/*         ia = (int*)ptr; */
 
137
/*         for(i=0;i<elems;i++) ia[i]  = 0; */
 
138
/*         break; */
 
139
/*       case C_DCPL: */
 
140
/*         elems *=2; */
 
141
/*       case C_DBL: */
 
142
/*         da = (double*)ptr; */
 
143
/*         for(i=0;i<elems;i++) da[i] = 0; */
 
144
/*         break; */
 
145
/*       case C_SCPL: */
 
146
/*         elems *=2; */
 
147
/*       case C_FLOAT: */
 
148
/*         fa = (float*)ptr; */
 
149
/*         for(i=0;i<elems;i++) fa[i]  = 0; */
 
150
/*         break; */
 
151
/*       case C_LONG: */
 
152
/*       la = (long*)ptr; */
 
153
/*       for(i=0;i<elems;i++) la[i]  = 0; */
 
154
/*       break; */
 
155
/*       case C_LONGLONG: */
 
156
/*       lla = (long long*)ptr; */
 
157
/*       for(i=0;i<elems;i++) lla[i]  = 0; */
 
158
/*       break; */
 
159
/*       default: pnga_error(" wrong data type ",type); */
 
160
/*     } */
 
161
      memset(ptr, 0, GAsizeofM(type)*elems);
 
162
 
 
163
    /* release access to the data */
 
164
    pnga_release_update_block_segment(g_a, me);
 
165
  }
 
166
  if(local_sync_end)pnga_pgroup_sync(p_handle);
 
167
  GA_POP_NAME;
 
168
}
 
169
 
 
170
 
 
171
 
 
172
#if 0
 
173
/*\ COPY ONE GLOBAL ARRAY INTO ANOTHER
 
174
\*/
 
175
static void snga_copy_old(Integer g_a, Integer g_b)
 
176
{
 
177
Integer  ndim, ndimb, type, typeb, me, elems=0, elemsb=0;
 
178
Integer dimsb[MAXDIM];
 
179
void *ptr_a, *ptr_b;
 
180
 
 
181
   me = pnga_nodeid();
 
182
 
 
183
   GA_PUSH_NAME("ga_copy");
 
184
 
 
185
   if(g_a == g_b) pnga_error("arrays have to be different ", 0L);
 
186
 
 
187
   pnga_inquire(g_a,  &type, &ndim, dims);
 
188
   pnga_inquire(g_b,  &typeb, &ndimb, dimsb);
 
189
 
 
190
   if(type != typeb) pnga_error("types not the same", g_b);
 
191
 
 
192
   if(!pnga_compare_distr(g_a,g_b))
 
193
 
 
194
      pnga_copy_patch("n",g_a, one_arr, dims, g_b, one_arr, dimsb);
 
195
 
 
196
   else {
 
197
 
 
198
     pnga_sync();
 
199
 
 
200
     pnga_distribution(g_a, me, lo, hi);
 
201
     if(lo[0]>0){
 
202
        pnga_access_ptr(g_a, lo, hi, &ptr_a, ld);
 
203
        if (pnga_has_ghosts(g_a)) {
 
204
          GET_ELEMS_W_GHOSTS(ndim,lo,hi,ld,&elems);
 
205
        } else {
 
206
          GET_ELEMS(ndim,lo,hi,ld,&elems);
 
207
        }
 
208
     }
 
209
 
 
210
     pnga_distribution(g_b, me, lo, hi);
 
211
     if(lo[0]>0){
 
212
        pnga_access_ptr(g_b, lo, hi, &ptr_b, ld);
 
213
        if (pnga_has_ghosts(g_b)) {
 
214
          GET_ELEMS_W_GHOSTS(ndim,lo,hi,ld,&elems);
 
215
        } else {
 
216
          GET_ELEMS(ndim,lo,hi,ld,&elems);
 
217
        }
 
218
     }
 
219
  
 
220
     if(elems!= elemsb)pnga_error("inconsistent number of elements",elems-elemsb);
 
221
 
 
222
     if(elems>0){
 
223
        ARMCI_Copy(ptr_a, ptr_b, (int)elems*GAsizeofM(type));
 
224
        pnga_release(g_a,lo,hi);
 
225
        pnga_release(g_b,lo,hi);
 
226
     }
 
227
 
 
228
     pnga_sync();
 
229
   }
 
230
 
 
231
   GA_POP_NAME;
 
232
}
 
233
#endif
 
234
 
 
235
 
 
236
 
 
237
/*\ COPY ONE GLOBAL ARRAY INTO ANOTHER
 
238
\*/
 
239
#if HAVE_SYS_WEAK_ALIAS_PRAGMA
 
240
#   pragma weak wnga_copy = pnga_copy
 
241
#endif
 
242
void pnga_copy(Integer g_a, Integer g_b)
 
243
{
 
244
Integer  ndim, ndimb, type, typeb, me_a, me_b;
 
245
Integer dimsb[MAXDIM],i;
 
246
Integer a_grp, b_grp, anproc, bnproc;
 
247
Integer num_blocks_a, num_blocks_b;
 
248
Integer blocks[MAXDIM], block_dims[MAXDIM];
 
249
void *ptr_a, *ptr_b;
 
250
int local_sync_begin,local_sync_end,use_put;
 
251
 
 
252
   GA_PUSH_NAME("ga_copy");
 
253
 
 
254
   local_sync_begin = _ga_sync_begin; local_sync_end = _ga_sync_end;
 
255
   _ga_sync_begin = 1; _ga_sync_end=1; /*remove any previous masking*/
 
256
   a_grp = pnga_get_pgroup(g_a);
 
257
   b_grp = pnga_get_pgroup(g_b);
 
258
   me_a = pnga_pgroup_nodeid(a_grp);
 
259
   me_b = pnga_pgroup_nodeid(b_grp);
 
260
   anproc = pnga_get_pgroup_size(a_grp);
 
261
   bnproc = pnga_get_pgroup_size(b_grp);
 
262
   num_blocks_a = pnga_total_blocks(g_a);
 
263
   num_blocks_b = pnga_total_blocks(g_b);
 
264
   if (anproc <= bnproc) {
 
265
     use_put = 1;
 
266
   } else {
 
267
     use_put = 0;
 
268
   }
 
269
   /*if (a_grp != b_grp)
 
270
     pnga_error("Both arrays must be defined on same group",0L); */
 
271
   if(local_sync_begin) {
 
272
     if (anproc <= bnproc) {
 
273
       pnga_pgroup_sync(a_grp);
 
274
     } else if (a_grp == pnga_pgroup_get_world() &&
 
275
                b_grp == pnga_pgroup_get_world()) {
 
276
       pnga_sync();
 
277
     } else {
 
278
       pnga_pgroup_sync(b_grp);
 
279
     }
 
280
   }
 
281
 
 
282
   if(g_a == g_b) pnga_error("arrays have to be different ", 0L);
 
283
 
 
284
   pnga_inquire(g_a,  &type, &ndim, dims);
 
285
   pnga_inquire(g_b,  &typeb, &ndimb, dimsb);
 
286
 
 
287
   if(type != typeb) pnga_error("types not the same", g_b);
 
288
   if(ndim != ndimb) pnga_error("dimensions not the same", ndimb);
 
289
 
 
290
   for(i=0; i< ndim; i++)if(dims[i]!=dimsb[i]) 
 
291
                          pnga_error("dimensions not the same",i);
 
292
 
 
293
   if ((pnga_is_mirrored(g_a) && pnga_is_mirrored(g_b)) ||
 
294
       (!pnga_is_mirrored(g_a) && !pnga_is_mirrored(g_b))) {
 
295
     /* Both global arrays are mirrored or both global arrays are not mirrored.
 
296
        Copy operation is straightforward */
 
297
 
 
298
     if (use_put) {
 
299
       if (num_blocks_a < 0) {
 
300
         pnga_distribution(g_a, me_a, lo, hi);
 
301
         if(lo[0]>0){
 
302
           pnga_access_ptr(g_a, lo, hi, &ptr_a, ld);
 
303
           pnga_put(g_b, lo, hi, ptr_a, ld);
 
304
         }
 
305
       } else {
 
306
         if (!pnga_uses_proc_grid(g_a)) {
 
307
           for (i=me_a; i<num_blocks_a; i += anproc) {
 
308
             pnga_distribution(g_a, i, lo, hi);
 
309
             if (lo[0]>0) {
 
310
               pnga_access_block_ptr(g_a, i, &ptr_a, ld);
 
311
               pnga_put(g_b, lo, hi, ptr_a, ld);
 
312
             }
 
313
           }
 
314
         } else {
 
315
           Integer proc_index[MAXDIM], index[MAXDIM];
 
316
           Integer topology[MAXDIM], chk;
 
317
           pnga_get_proc_index(g_a, me_a, proc_index);
 
318
           pnga_get_proc_index(g_a, me_a, index);
 
319
           pnga_get_block_info(g_a, blocks, block_dims);
 
320
           pnga_get_proc_grid(g_a, topology);
 
321
           while (index[ndim-1] < blocks[ndim-1]) {
 
322
             /* find bounding coordinates of block */
 
323
             chk = 1;
 
324
             for (i = 0; i < ndim; i++) {
 
325
               lo[i] = index[i]*block_dims[i]+1;
 
326
               hi[i] = (index[i] + 1)*block_dims[i];
 
327
               if (hi[i] > dims[i]) hi[i] = dims[i];
 
328
               if (hi[i] < lo[i]) chk = 0;
 
329
             }
 
330
             if (chk) {
 
331
               pnga_access_block_grid_ptr(g_a, index, &ptr_a, ld);
 
332
               pnga_put(g_b, lo, hi, ptr_a, ld);
 
333
             }
 
334
             /* increment index to get next block on processor */
 
335
             index[0] += topology[0];
 
336
             for (i = 0; i < ndim; i++) {
 
337
               if (index[i] >= blocks[i] && i<ndim-1) {
 
338
                 index[i] = proc_index[i];
 
339
                 index[i+1] += topology[i+1];
 
340
               }
 
341
             }
 
342
           }
 
343
         }
 
344
       }
 
345
     } else {
 
346
       if (num_blocks_b < 0) {
 
347
         pnga_distribution(g_b, me_b, lo, hi);
 
348
         if(lo[0]>0){
 
349
           pnga_access_ptr(g_b, lo, hi, &ptr_b, ld);
 
350
           pnga_get(g_a, lo, hi, ptr_b, ld);
 
351
         }
 
352
       } else {
 
353
         if (!pnga_uses_proc_grid(g_a)) {
 
354
           for (i=me_b; i<num_blocks_b; i += bnproc) {
 
355
             pnga_distribution(g_b, i, lo, hi);
 
356
             if (lo[0]>0) {
 
357
               pnga_access_block_ptr(g_b, i, &ptr_b, ld);
 
358
               pnga_get(g_a, lo, hi, ptr_b, ld);
 
359
             }
 
360
           }
 
361
         } else {
 
362
           Integer proc_index[MAXDIM], index[MAXDIM];
 
363
           Integer topology[MAXDIM], chk;
 
364
           pnga_get_proc_index(g_b, me_b, proc_index);
 
365
           pnga_get_proc_index(g_b, me_b, index);
 
366
           pnga_get_block_info(g_b, blocks, block_dims);
 
367
           pnga_get_proc_grid(g_b, topology);
 
368
           while (index[ndim-1] < blocks[ndim-1]) {
 
369
             /* find bounding coordinates of block */
 
370
             chk = 1;
 
371
             for (i = 0; i < ndim; i++) {
 
372
               lo[i] = index[i]*block_dims[i]+1;
 
373
               hi[i] = (index[i] + 1)*block_dims[i];
 
374
               if (hi[i] > dims[i]) hi[i] = dims[i];
 
375
               if (hi[i] < lo[i]) chk = 0;
 
376
             }
 
377
             if (chk) {
 
378
               pnga_access_block_grid_ptr(g_b, index, &ptr_b, ld);
 
379
               pnga_get(g_a, lo, hi, ptr_b, ld);
 
380
             }
 
381
             /* increment index to get next block on processor */
 
382
             index[0] += topology[0];
 
383
             for (i = 0; i < ndim; i++) {
 
384
               if (index[i] >= blocks[i] && i<ndim-1) {
 
385
                 index[i] = proc_index[i];
 
386
                 index[i+1] += topology[i+1];
 
387
               }
 
388
             }
 
389
           }
 
390
         }
 
391
       }
 
392
     }
 
393
   } else {
 
394
     /* One global array is mirrored and the other is not */
 
395
     if (pnga_is_mirrored(g_a)) {
 
396
       /* Source array is mirrored and destination
 
397
          array is distributed. Assume source array is consistent */
 
398
       pnga_distribution(g_b, me_b, lo, hi);
 
399
       if (lo[0]>0) {
 
400
         pnga_access_ptr(g_b, lo, hi, &ptr_b, ld);
 
401
         pnga_get(g_a, lo, hi, ptr_b, ld);
 
402
       } 
 
403
     } else {
 
404
       /* source array is distributed and destination
 
405
          array is mirrored */
 
406
       pnga_zero(g_b);
 
407
       pnga_distribution(g_a, me_a, lo, hi);
 
408
       if (lo[0] > 0) {
 
409
         pnga_access_ptr(g_a, lo, hi, &ptr_a, ld);
 
410
         pnga_put(g_b, lo, hi, ptr_a, ld);
 
411
       }
 
412
       pnga_merge_mirrored(g_b);
 
413
     }
 
414
   }
 
415
 
 
416
   if(local_sync_end) {
 
417
     if (anproc <= bnproc) {
 
418
       pnga_pgroup_sync(a_grp);
 
419
     } else if (a_grp == pnga_pgroup_get_world() &&
 
420
                b_grp == pnga_pgroup_get_world()) {
 
421
       pnga_sync();
 
422
     } else {
 
423
       pnga_pgroup_sync(b_grp);
 
424
     }
 
425
   }
 
426
   GA_POP_NAME;
 
427
}
 
428
 
 
429
 
 
430
 
 
431
/*\ internal version of dot product
 
432
\*/
 
433
#if HAVE_SYS_WEAK_ALIAS_PRAGMA
 
434
#   pragma weak wnga_dot = pnga_dot
 
435
#endif
 
436
void pnga_dot(int Type, Integer g_a, Integer g_b, void *value)
 
437
{
 
438
Integer  ndim=0, type=0, atype=0, me=0, elems=0, elemsb=0;
 
439
register Integer i=0;
 
440
int isum=0;
 
441
long lsum=0;
 
442
long long llsum=0;
 
443
DoubleComplex zsum ={0.,0.};
 
444
SingleComplex csum ={0.,0.};
 
445
float fsum=0.0;
 
446
void *ptr_a=NULL, *ptr_b=NULL;
 
447
int alen=0;
 
448
Integer a_grp=0, b_grp=0;
 
449
Integer num_blocks_a=0, num_blocks_b=0;
 
450
 
 
451
Integer andim, adims[MAXDIM];
 
452
Integer bndim, bdims[MAXDIM];
 
453
 
 
454
   _ga_sync_begin = 1; _ga_sync_end=1; /*remove any previous masking*/
 
455
 
 
456
   GA_PUSH_NAME("ga_dot");
 
457
   a_grp = pnga_get_pgroup(g_a);
 
458
   b_grp = pnga_get_pgroup(g_b);
 
459
   if (a_grp != b_grp)
 
460
     pnga_error("Both arrays must be defined on same group",0L);
 
461
   me = pnga_pgroup_nodeid(a_grp);
 
462
 
 
463
   /* Check to see if either GA is block cyclic distributed */
 
464
   num_blocks_a = pnga_total_blocks(g_a);
 
465
   num_blocks_b = pnga_total_blocks(g_b);
 
466
   if (num_blocks_a >= 0 || num_blocks_b >= 0) {
 
467
     pnga_inquire(g_a, &type, &andim, adims);
 
468
     pnga_inquire(g_b, &type, &bndim, bdims);
 
469
     pnga_dot_patch(g_a, "n", one_arr, adims, g_b, "n", one_arr, bdims,
 
470
         value);
 
471
     GA_POP_NAME;
 
472
     return;
 
473
   }
 
474
 
 
475
   if(pnga_compare_distr(g_a,g_b) == FALSE ||
 
476
      pnga_has_ghosts(g_a) || pnga_has_ghosts(g_b)) {
 
477
       /* distributions not identical */
 
478
       pnga_inquire(g_a, &type, &andim, adims);
 
479
       pnga_inquire(g_b, &type, &bndim, bdims);
 
480
 
 
481
       pnga_dot_patch(g_a, "n", one_arr, adims, g_b, "n", one_arr, bdims,
 
482
                      value);
 
483
       
 
484
       GA_POP_NAME;
 
485
       return;
 
486
   }
 
487
   
 
488
   pnga_pgroup_sync(a_grp);
 
489
   pnga_inquire(g_a,  &type, &ndim, dims);
 
490
   if(type != Type) pnga_error("type not correct", g_a);
 
491
   pnga_distribution(g_a, me, lo, hi);
 
492
   if(lo[0]>0){
 
493
      pnga_access_ptr(g_a, lo, hi, &ptr_a, ld);
 
494
      if (pnga_has_ghosts(g_a)) {
 
495
        GET_ELEMS_W_GHOSTS(ndim,lo,hi,ld,&elems);
 
496
      } else {
 
497
        GET_ELEMS(ndim,lo,hi,ld,&elems);
 
498
      }
 
499
   }
 
500
 
 
501
   if(g_a == g_b){
 
502
     elemsb = elems;
 
503
     ptr_b = ptr_a;
 
504
   }else {  
 
505
     pnga_inquire(g_b,  &type, &ndim, dims);
 
506
     if(type != Type) pnga_error("type not correct", g_b);
 
507
     pnga_distribution(g_b, me, lo, hi);
 
508
     if(lo[0]>0){
 
509
        pnga_access_ptr(g_b, lo, hi, &ptr_b, ld);
 
510
        if (pnga_has_ghosts(g_b)) {
 
511
          GET_ELEMS_W_GHOSTS(ndim,lo,hi,ld,&elemsb);
 
512
        } else {
 
513
          GET_ELEMS(ndim,lo,hi,ld,&elemsb);
 
514
        }
 
515
     }
 
516
   }
 
517
 
 
518
   if(elems!= elemsb)pnga_error("inconsistent number of elements",elems-elemsb); 
 
519
 
 
520
 
 
521
      /* compute "local" contribution to the dot product */
 
522
      switch (type){
 
523
        int *ia, *ib;
 
524
        double *da,*db;
 
525
        float *fa, *fb;
 
526
        long *la,*lb;
 
527
        long long *lla,*llb;
 
528
        case C_INT:
 
529
           ia = (int*)ptr_a;
 
530
           ib = (int*)ptr_b;
 
531
           for(i=0;i<elems;i++) 
 
532
                 isum += ia[i]  * ib[i];
 
533
           *(int*)value = isum; 
 
534
           type = C_INT;
 
535
           alen = 1;
 
536
           break;
 
537
 
 
538
        case C_DCPL:
 
539
           for(i=0;i<elems;i++){
 
540
               DoubleComplex a = ((DoubleComplex*)ptr_a)[i];
 
541
               DoubleComplex b = ((DoubleComplex*)ptr_b)[i];
 
542
               zsum.real += a.real*b.real  - b.imag * a.imag;
 
543
               zsum.imag += a.imag*b.real  + b.imag * a.real;
 
544
           }
 
545
           *(DoubleComplex*)value = zsum; 
 
546
           type = C_DCPL;
 
547
           alen = 2;
 
548
           break;
 
549
 
 
550
        case C_SCPL:
 
551
           for(i=0;i<elems;i++){
 
552
               SingleComplex a = ((SingleComplex*)ptr_a)[i];
 
553
               SingleComplex b = ((SingleComplex*)ptr_b)[i];
 
554
               csum.real += a.real*b.real  - b.imag * a.imag;
 
555
               csum.imag += a.imag*b.real  + b.imag * a.real;
 
556
           }
 
557
           *(SingleComplex*)value = csum; 
 
558
           type = C_SCPL;
 
559
           alen = 2;
 
560
           break;
 
561
 
 
562
        case C_DBL:
 
563
           da = (double*)ptr_a;
 
564
           db = (double*)ptr_b;
 
565
           for(i=0;i<elems;i++) 
 
566
                 zsum.real += da[i]  * db[i];
 
567
           *(double*)value = zsum.real; 
 
568
           type = C_DBL;
 
569
           alen = 1;
 
570
           break;
 
571
        case C_FLOAT:
 
572
           fa = (float*)ptr_a;
 
573
           fb = (float*)ptr_b;
 
574
           for(i=0;i<elems;i++)
 
575
                 fsum += fa[i]  * fb[i];
 
576
           *(float*)value = fsum;
 
577
           type = C_FLOAT;
 
578
           alen = 1;
 
579
           break;         
 
580
        case C_LONG:
 
581
           la = (long*)ptr_a;
 
582
           lb = (long*)ptr_b;
 
583
           for(i=0;i<elems;i++)
 
584
                lsum += la[i]  * lb[i];
 
585
           *(long*)value = lsum;
 
586
           type = C_LONG;
 
587
           alen = 1;
 
588
           break;               
 
589
        case C_LONGLONG:
 
590
           lla = (long long*)ptr_a;
 
591
           llb = (long long*)ptr_b;
 
592
           for(i=0;i<elems;i++)
 
593
                llsum += lla[i]  * llb[i];
 
594
           *(long long*)value = llsum;
 
595
           type = C_LONGLONG;
 
596
           alen = 1;
 
597
           break;               
 
598
        default: pnga_error(" wrong data type ",type);
 
599
      }
 
600
   
 
601
      /* release access to the data */
 
602
      if(elems>0){
 
603
         pnga_release(g_a, lo, hi);
 
604
         if(g_a != g_b)pnga_release(g_b, lo, hi);
 
605
      }
 
606
 
 
607
    /*convert from C data type to ARMCI type */
 
608
    switch(type) {
 
609
      case C_FLOAT: atype=ARMCI_FLOAT; break;
 
610
      case C_DBL: atype=ARMCI_DOUBLE; break;
 
611
      case C_INT: atype=ARMCI_INT; break;
 
612
      case C_LONG: atype=ARMCI_LONG; break;
 
613
      case C_LONGLONG: atype=ARMCI_LONG_LONG; break;
 
614
      case C_DCPL: atype=ARMCI_DOUBLE; break;
 
615
      case C_SCPL: atype=ARMCI_FLOAT; break;
 
616
      default: pnga_error("pnga_dot: type not supported",type);
 
617
    }
 
618
 
 
619
   if (pnga_is_mirrored(g_a) && pnga_is_mirrored(g_b)) {
 
620
     armci_msg_gop_scope(SCOPE_NODE,value,alen,"+",atype);
 
621
   } else {
 
622
     if (a_grp == -1) {
 
623
       armci_msg_gop_scope(SCOPE_ALL,value,alen,"+",atype);
 
624
#ifdef MPI
 
625
     } else {
 
626
       armci_msg_group_gop_scope(SCOPE_ALL,value,alen,"+",atype,
 
627
           ga_get_armci_group_((int)a_grp));
 
628
#endif
 
629
     }
 
630
   }
 
631
    
 
632
   GA_POP_NAME;
 
633
 
 
634
}
 
635
 
 
636
 
 
637
#if HAVE_SYS_WEAK_ALIAS_PRAGMA
 
638
#   pragma weak wnga_scale = pnga_scale
 
639
#endif
 
640
void pnga_scale(Integer g_a, void* alpha)
 
641
{
 
642
  Integer ndim, type, me, elems, grp_id;
 
643
  register Integer i;
 
644
  Integer num_blocks;
 
645
  void *ptr;
 
646
  int local_sync_begin,local_sync_end;
 
647
 
 
648
  local_sync_begin = _ga_sync_begin; local_sync_end = _ga_sync_end;
 
649
  _ga_sync_begin = 1; _ga_sync_end=1; /*remove any previous masking*/
 
650
  grp_id = pnga_get_pgroup(g_a);
 
651
  if(local_sync_begin)pnga_pgroup_sync(grp_id);
 
652
 
 
653
  me = pnga_pgroup_nodeid(grp_id);
 
654
 
 
655
  pnga_check_handle(g_a, "ga_scale");
 
656
  GA_PUSH_NAME("ga_scale");
 
657
  num_blocks = pnga_total_blocks(g_a);
 
658
 
 
659
  pnga_inquire(g_a, &type, &ndim, dims);
 
660
  if (num_blocks < 0) {
 
661
    pnga_distribution(g_a, me, lo, hi);
 
662
    if (pnga_has_ghosts(g_a)) {
 
663
      pnga_scale_patch(g_a, lo, hi, alpha);
 
664
      return;
 
665
    }
 
666
 
 
667
    if ( lo[0]> 0 ){ /* base index is 1: we get 0 if no elements stored on p */
 
668
 
 
669
      pnga_access_ptr(g_a, lo, hi, &ptr, ld);
 
670
      GET_ELEMS(ndim,lo,hi,ld,&elems);
 
671
 
 
672
      switch (type){
 
673
        int *ia;
 
674
        double *da;
 
675
        DoubleComplex *ca, scale;
 
676
        SingleComplex *cfa, cfscale;
 
677
        long *la;
 
678
        long long *lla;
 
679
        float *fa;
 
680
        case C_INT:
 
681
        ia = (int*)ptr;
 
682
        for(i=0;i<elems;i++) ia[i]  *= *(int*)alpha;
 
683
        break;
 
684
        case C_LONG:
 
685
        la = (long*)ptr;
 
686
        for(i=0;i<elems;i++) la[i]  *= *(long*)alpha;
 
687
        break;
 
688
        case C_LONGLONG:
 
689
        lla = (long long*)ptr;
 
690
        for(i=0;i<elems;i++) lla[i]  *= *(long long*)alpha;
 
691
        break;
 
692
        case C_DCPL:
 
693
        ca = (DoubleComplex*)ptr;
 
694
        scale= *(DoubleComplex*)alpha;
 
695
        for(i=0;i<elems;i++){
 
696
          DoubleComplex val = ca[i]; 
 
697
          ca[i].real = scale.real*val.real  - val.imag * scale.imag;
 
698
          ca[i].imag = scale.imag*val.real  + val.imag * scale.real;
 
699
        }
 
700
        break;
 
701
        case C_SCPL:
 
702
        cfa = (SingleComplex*)ptr;
 
703
        cfscale= *(SingleComplex*)alpha;
 
704
        for(i=0;i<elems;i++){
 
705
          SingleComplex val = cfa[i]; 
 
706
          cfa[i].real = cfscale.real*val.real  - val.imag * cfscale.imag;
 
707
          cfa[i].imag = cfscale.imag*val.real  + val.imag * cfscale.real;
 
708
        }
 
709
        break;
 
710
        case C_DBL:
 
711
        da = (double*)ptr;
 
712
        for(i=0;i<elems;i++) da[i] *= *(double*)alpha;
 
713
        break;
 
714
        case C_FLOAT:
 
715
        fa = (float*)ptr;
 
716
        for(i=0;i<elems;i++) fa[i]  *= *(float*)alpha;
 
717
        break;       
 
718
        default: pnga_error(" wrong data type ",type);
 
719
      }
 
720
 
 
721
      /* release access to the data */
 
722
      pnga_release_update(g_a, lo, hi);
 
723
    }
 
724
  } else {
 
725
    pnga_access_block_segment_ptr(g_a, me, &ptr, &elems);
 
726
    switch (type){
 
727
      int *ia;
 
728
      double *da;
 
729
      DoubleComplex *ca, scale;
 
730
      SingleComplex *cfa, cfscale;
 
731
      long *la;
 
732
      long long *lla;
 
733
      float *fa;
 
734
      case C_INT:
 
735
      ia = (int*)ptr;
 
736
      for(i=0;i<elems;i++) ia[i]  *= *(int*)alpha;
 
737
      break;
 
738
      case C_LONG:
 
739
      la = (long*)ptr;
 
740
      for(i=0;i<elems;i++) la[i]  *= *(long*)alpha;
 
741
      break;
 
742
      case C_LONGLONG:
 
743
      lla = (long long*)ptr;
 
744
      for(i=0;i<elems;i++) lla[i]  *= *(long long*)alpha;
 
745
      break;
 
746
      case C_DCPL:
 
747
      ca = (DoubleComplex*)ptr;
 
748
      scale= *(DoubleComplex*)alpha;
 
749
      for(i=0;i<elems;i++){
 
750
        DoubleComplex val = ca[i]; 
 
751
        ca[i].real = scale.real*val.real  - val.imag * scale.imag;
 
752
        ca[i].imag = scale.imag*val.real  + val.imag * scale.real;
 
753
      }
 
754
      break;
 
755
      case C_SCPL:
 
756
      cfa = (SingleComplex*)ptr;
 
757
      cfscale= *(SingleComplex*)alpha;
 
758
      for(i=0;i<elems;i++){
 
759
        SingleComplex val = cfa[i]; 
 
760
        cfa[i].real = cfscale.real*val.real  - val.imag * cfscale.imag;
 
761
        cfa[i].imag = cfscale.imag*val.real  + val.imag * cfscale.real;
 
762
      }
 
763
      break;
 
764
      case C_DBL:
 
765
      da = (double*)ptr;
 
766
      for(i=0;i<elems;i++) da[i] *= *(double*)alpha;
 
767
      break;
 
768
      case C_FLOAT:
 
769
      fa = (float*)ptr;
 
770
      for(i=0;i<elems;i++) fa[i]  *= *(float*)alpha;
 
771
      break;       
 
772
      default: pnga_error(" wrong data type ",type);
 
773
    }
 
774
    /* release access to the data */
 
775
    pnga_release_update_block_segment(g_a, me);
 
776
  }
 
777
  GA_POP_NAME;
 
778
  if(local_sync_end)pnga_pgroup_sync(grp_id); 
 
779
}
 
780
 
 
781
 
 
782
#if HAVE_SYS_WEAK_ALIAS_PRAGMA
 
783
#   pragma weak wnga_add = pnga_add
 
784
#endif
 
785
void pnga_add(void *alpha, Integer g_a, void* beta, Integer g_b, Integer g_c)
 
786
{
 
787
Integer  ndim, type, typeC, me, elems=0, elemsb=0, elemsa=0;
 
788
register Integer i;
 
789
void *ptr_a, *ptr_b, *ptr_c;
 
790
Integer a_grp, b_grp, c_grp;
 
791
int local_sync_begin,local_sync_end;
 
792
 
 
793
 Integer andim, adims[MAXDIM];
 
794
 Integer bndim, bdims[MAXDIM];
 
795
 Integer cndim, cdims[MAXDIM];
 
796
 
 
797
   local_sync_begin = _ga_sync_begin; local_sync_end = _ga_sync_end;
 
798
   _ga_sync_begin = 1; _ga_sync_end=1; /*remove any previous masking*/
 
799
 
 
800
   GA_PUSH_NAME("ga_add");
 
801
   a_grp = pnga_get_pgroup(g_a);
 
802
   b_grp = pnga_get_pgroup(g_b);
 
803
   c_grp = pnga_get_pgroup(g_c);
 
804
   if (a_grp != b_grp || b_grp != c_grp)
 
805
     pnga_error("All three arrays must be on same group for ga_add",0L);
 
806
 
 
807
   if(local_sync_begin)pnga_pgroup_sync(a_grp);
 
808
 
 
809
   me = pnga_pgroup_nodeid(a_grp);
 
810
   if((pnga_compare_distr(g_a,g_b) == FALSE) ||
 
811
      (pnga_compare_distr(g_a,g_c) == FALSE) ||
 
812
       pnga_has_ghosts(g_a) || pnga_has_ghosts(g_b) || pnga_has_ghosts(g_c) ||
 
813
       pnga_total_blocks(g_a) > 0 || pnga_total_blocks(g_b) > 0 ||
 
814
       pnga_total_blocks(g_c) > 0) {
 
815
       /* distributions not identical */
 
816
       pnga_inquire(g_a, &type, &andim, adims);
 
817
       pnga_inquire(g_b, &type, &bndim, bdims);
 
818
       pnga_inquire(g_b, &type, &cndim, cdims);
 
819
 
 
820
       pnga_add_patch(alpha, g_a, one_arr, adims, beta, g_b, one_arr, bdims,
 
821
                      g_c, one_arr, cdims);
 
822
       
 
823
       GA_POP_NAME;
 
824
       return;
 
825
   }
 
826
 
 
827
   pnga_pgroup_sync(a_grp);
 
828
   pnga_inquire(g_c,  &typeC, &ndim, dims);
 
829
   pnga_distribution(g_c, me, lo, hi);
 
830
   if (  lo[0]>0 ){
 
831
     pnga_access_ptr(g_c, lo, hi, &ptr_c, ld);
 
832
     GET_ELEMS(ndim,lo,hi,ld,&elems);
 
833
   }
 
834
 
 
835
   if(g_a == g_c){
 
836
     ptr_a  = ptr_c;
 
837
     elemsa = elems;
 
838
   }else { 
 
839
     pnga_inquire(g_a,  &type, &ndim, dims);
 
840
     if(type != typeC) pnga_error("types not consistent", g_a);
 
841
     pnga_distribution(g_a, me, lo, hi);
 
842
     if (  lo[0]>0 ){
 
843
       pnga_access_ptr(g_a, lo, hi, &ptr_a, ld);
 
844
       GET_ELEMS(ndim,lo,hi,ld,&elemsa);
 
845
     }
 
846
   }
 
847
 
 
848
   if(g_b == g_c){
 
849
     ptr_b  = ptr_c;
 
850
     elemsb = elems;
 
851
   }else {
 
852
     pnga_inquire(g_b,  &type, &ndim, dims);
 
853
     if(type != typeC) pnga_error("types not consistent", g_b);
 
854
     pnga_distribution(g_b, me, lo, hi);
 
855
     if (  lo[0]>0 ){
 
856
       pnga_access_ptr(g_b, lo, hi, &ptr_b, ld);
 
857
       GET_ELEMS(ndim,lo,hi,ld,&elemsb);
 
858
     }
 
859
   }
 
860
 
 
861
   if(elems!= elemsb)pnga_error("inconsistent number of elements a",elems-elemsb);
 
862
   if(elems!= elemsa)pnga_error("inconsistent number of elements b",elems-elemsa);
 
863
 
 
864
   if (  lo[0]>0 ){
 
865
 
 
866
       /* operation on the "local" piece of data */
 
867
       switch(type){
 
868
         int *ia, *ib, *ic;
 
869
         double *da,*db,*dc;
 
870
         float *fa, *fb, *fc;
 
871
         long *la,*lb,*lc;
 
872
         long long *lla,*llb,*llc;
 
873
         case C_DBL:
 
874
                  da = (double*)ptr_a;
 
875
                  db = (double*)ptr_b;
 
876
                  dc = (double*)ptr_c;
 
877
                  for(i=0; i<elems; i++)
 
878
                      dc[i] = *(double*)alpha *da[i] +
 
879
                              *(double*)beta * db[i];
 
880
              break;
 
881
         case C_DCPL:
 
882
                  for(i=0; i<elems; i++){
 
883
                     DoubleComplex a = ((DoubleComplex*)ptr_a)[i];
 
884
                     DoubleComplex b = ((DoubleComplex*)ptr_b)[i];
 
885
                     DoubleComplex *ac = (DoubleComplex*)ptr_c;
 
886
                     DoubleComplex x= *(DoubleComplex*)alpha;
 
887
                     DoubleComplex y= *(DoubleComplex*)beta;
 
888
                     /* c = x*a + y*b */
 
889
                     ac[i].real = x.real*a.real - 
 
890
                              x.imag*a.imag + y.real*b.real - y.imag*b.imag;
 
891
                     ac[i].imag = x.real*a.imag + 
 
892
                              x.imag*a.real + y.real*b.imag + y.imag*b.real;
 
893
                  }
 
894
              break;
 
895
         case C_SCPL:
 
896
                  for(i=0; i<elems; i++){
 
897
                     SingleComplex a = ((SingleComplex*)ptr_a)[i];
 
898
                     SingleComplex b = ((SingleComplex*)ptr_b)[i];
 
899
                     SingleComplex *ac = (SingleComplex*)ptr_c;
 
900
                     SingleComplex x= *(SingleComplex*)alpha;
 
901
                     SingleComplex y= *(SingleComplex*)beta;
 
902
                     /* c = x*a + y*b */
 
903
                     ac[i].real = x.real*a.real - 
 
904
                              x.imag*a.imag + y.real*b.real - y.imag*b.imag;
 
905
                     ac[i].imag = x.real*a.imag + 
 
906
                              x.imag*a.real + y.real*b.imag + y.imag*b.real;
 
907
                  }
 
908
              break;
 
909
         case C_FLOAT:
 
910
                  fa = (float*)ptr_a;
 
911
                  fb = (float*)ptr_b;
 
912
                  fc = (float*)ptr_c;
 
913
                  for(i=0; i<elems; i++)
 
914
                      fc[i] = *(float*)alpha *fa[i] + *(float*)beta *fb[i]; 
 
915
              break;
 
916
         case C_INT:
 
917
                  ia = (int*)ptr_a;
 
918
                  ib = (int*)ptr_b;
 
919
                  ic = (int*)ptr_c;
 
920
                  for(i=0; i<elems; i++) 
 
921
                      ic[i] = *(int*)alpha *ia[i] + *(int*)beta *ib[i];
 
922
              break;    
 
923
         case C_LONG:
 
924
                  la = (long*)ptr_a;
 
925
                  lb = (long*)ptr_b;
 
926
                  lc = (long*)ptr_c;
 
927
                  for(i=0; i<elems; i++)
 
928
                      lc[i] = *(long*)alpha *la[i] + *(long*)beta *lb[i];
 
929
              break;
 
930
         case C_LONGLONG:
 
931
                  lla = (long long*)ptr_a;
 
932
                  llb = (long long*)ptr_b;
 
933
                  llc = (long long*)ptr_c;
 
934
                  for(i=0; i<elems; i++)
 
935
                      llc[i] = ( *(long long*)alpha *lla[i] +
 
936
                                 *(long long*)beta  * llb[i] );
 
937
                
 
938
       }
 
939
 
 
940
       /* release access to the data */
 
941
       pnga_release_update(g_c, lo, hi);
 
942
       if(g_c != g_a)pnga_release(g_a, lo, hi);
 
943
       if(g_c != g_b)pnga_release(g_b, lo, hi);
 
944
   }
 
945
 
 
946
 
 
947
   GA_POP_NAME;
 
948
   if(local_sync_end)pnga_pgroup_sync(a_grp);
 
949
}
 
950
 
 
951
 
 
952
static 
 
953
void snga_local_transpose(Integer type, char *ptra, Integer n, Integer stride, char *ptrb)
 
954
{
 
955
int i;
 
956
    switch(type){
 
957
 
 
958
       case C_INT:
 
959
            for(i = 0; i< n; i++, ptrb+= stride) 
 
960
               *(int*)ptrb= ((int*)ptra)[i];
 
961
            break;
 
962
       case C_DCPL:
 
963
            for(i = 0; i< n; i++, ptrb+= stride) 
 
964
               *(DoubleComplex*)ptrb= ((DoubleComplex*)ptra)[i];
 
965
            break;
 
966
       case C_SCPL:
 
967
            for(i = 0; i< n; i++, ptrb+= stride) 
 
968
               *(SingleComplex*)ptrb= ((SingleComplex*)ptra)[i];
 
969
            break;
 
970
       case C_DBL:
 
971
            for(i = 0; i< n; i++, ptrb+= stride) 
 
972
               *(double*)ptrb= ((double*)ptra)[i];
 
973
            break;
 
974
       case C_FLOAT:
 
975
            for(i = 0; i< n; i++, ptrb+= stride)
 
976
               *(float*)ptrb= ((float*)ptra)[i];
 
977
            break;      
 
978
       case C_LONG:
 
979
            for(i = 0; i< n; i++, ptrb+= stride)
 
980
               *(long*)ptrb= ((long*)ptra)[i];
 
981
            break;                                 
 
982
       case C_LONGLONG:
 
983
            for(i = 0; i< n; i++, ptrb+= stride)
 
984
               *(long long*)ptrb= ((long long*)ptra)[i];
 
985
            break;                                 
 
986
       default: pnga_error("bad type:",type);
 
987
    }
 
988
}
 
989
 
 
990
 
 
991
#if HAVE_SYS_WEAK_ALIAS_PRAGMA
 
992
#   pragma weak wnga_transpose = pnga_transpose
 
993
#endif
 
994
void pnga_transpose(Integer g_a, Integer g_b)
 
995
{
 
996
Integer me = pnga_nodeid();
 
997
Integer nproc = pnga_nnodes(); 
 
998
Integer atype, btype, andim, adims[MAXDIM], bndim, bdims[MAXDIM];
 
999
Integer lo[2],hi[2];
 
1000
int local_sync_begin,local_sync_end;
 
1001
Integer num_blocks_a;
 
1002
char *ptr_tmp, *ptr_a;
 
1003
 
 
1004
    GA_PUSH_NAME("ga_transpose");
 
1005
    
 
1006
    local_sync_begin = _ga_sync_begin; local_sync_end = _ga_sync_end;
 
1007
    _ga_sync_begin = 1; _ga_sync_end=1; /*remove any previous masking*/
 
1008
    if(local_sync_begin)pnga_sync();
 
1009
 
 
1010
    if(g_a == g_b) pnga_error("arrays have to be different ", 0L);
 
1011
 
 
1012
    pnga_inquire(g_a, &atype, &andim, adims);
 
1013
    pnga_inquire(g_b, &btype, &bndim, bdims);
 
1014
 
 
1015
    if(bndim != 2 || andim != 2) pnga_error("dimension must be 2",0);
 
1016
    if(atype != btype ) pnga_error("array type mismatch ", 0L);
 
1017
 
 
1018
    num_blocks_a = pnga_total_blocks(g_a);
 
1019
 
 
1020
    if (num_blocks_a < 0) {
 
1021
      pnga_distribution(g_a, me, lo, hi);
 
1022
 
 
1023
      if(lo[0]>0){
 
1024
        Integer nelem, lob[2], hib[2], nrow, ncol;
 
1025
        int i, size=GAsizeofM(atype);
 
1026
 
 
1027
        nrow   = hi[0] -lo[0]+1;
 
1028
        ncol   = hi[1] -lo[1]+1; 
 
1029
        nelem  = nrow*ncol;
 
1030
        lob[0] = lo[1]; lob[1] = lo[0];
 
1031
        hib[0] = hi[1]; hib[1] = hi[0];
 
1032
 
 
1033
        /* allocate memory for transposing elements locally */
 
1034
        ptr_tmp = (char *) ga_malloc(nelem, atype, "transpose_tmp");
 
1035
 
 
1036
        /* get access to local data */
 
1037
        pnga_access_ptr(g_a, lo, hi, &ptr_a, ld);
 
1038
 
 
1039
        for(i = 0; i < ncol; i++){
 
1040
          char *ptr = ptr_tmp + i*size;
 
1041
 
 
1042
          snga_local_transpose(atype, ptr_a, nrow, ncol*size, ptr);
 
1043
          ptr_a += ld[0]*size;
 
1044
        }
 
1045
 
 
1046
        pnga_release(g_a, lo, hi); 
 
1047
 
 
1048
        pnga_put(g_b, lob, hib, ptr_tmp ,&ncol);
 
1049
 
 
1050
        ga_free(ptr_tmp);
 
1051
      }
 
1052
    } else {
 
1053
      Integer idx;
 
1054
      Integer blocks[MAXDIM], block_dims[MAXDIM];
 
1055
      Integer nelem, lob[2], hib[2], nrow, ncol;
 
1056
      int i, size=GAsizeofM(atype);
 
1057
 
 
1058
      /* allocate memory for transposing elements locally */
 
1059
      pnga_get_block_info(g_a, blocks, block_dims);
 
1060
 
 
1061
      /* Simple block-cyclic data distribution */
 
1062
      nelem = block_dims[0]*block_dims[1];
 
1063
      ptr_tmp = (char *) ga_malloc(nelem, atype, "transpose_tmp");
 
1064
      if (!pnga_uses_proc_grid(g_a)) {
 
1065
        for (idx = me; idx < num_blocks_a; idx += nproc) {
 
1066
          pnga_distribution(g_a, idx, lo, hi);
 
1067
          pnga_access_block_ptr(g_a, idx, &ptr_a, ld);
 
1068
 
 
1069
          nrow   = hi[0] -lo[0]+1;
 
1070
          ncol   = hi[1] -lo[1]+1; 
 
1071
          nelem  = nrow*ncol;
 
1072
          lob[0] = lo[1]; lob[1] = lo[0];
 
1073
          hib[0] = hi[1]; hib[1] = hi[0];
 
1074
          for(i = 0; i < ncol; i++){
 
1075
            char *ptr = ptr_tmp + i*size;
 
1076
 
 
1077
            snga_local_transpose(atype, ptr_a, nrow, ncol*size, ptr);
 
1078
            ptr_a += ld[0]*size;
 
1079
          }
 
1080
          pnga_put(g_b, lob, hib, ptr_tmp ,&ncol);
 
1081
 
 
1082
          pnga_release_update_block(g_a, idx);
 
1083
        }
 
1084
      } else {
 
1085
        /* Uses scalapack block-cyclic data distribution */
 
1086
        Integer chk;
 
1087
        Integer proc_index[MAXDIM], index[MAXDIM];
 
1088
        Integer topology[MAXDIM], ichk;
 
1089
 
 
1090
        pnga_get_proc_index(g_a, me, proc_index);
 
1091
        pnga_get_proc_index(g_a, me, index);
 
1092
        pnga_get_block_info(g_a, blocks, block_dims);
 
1093
        pnga_get_proc_grid(g_a, topology);
 
1094
        /* Verify that processor has data */
 
1095
        ichk = 1;
 
1096
        for (i=0; i<andim; i++) {
 
1097
          if (index[i]<0 || index[i] >= blocks[i]) {
 
1098
            ichk = 0;
 
1099
          }
 
1100
        }
 
1101
 
 
1102
        if (ichk) {
 
1103
          pnga_access_block_grid_ptr(g_a, index, &ptr_a, ld);
 
1104
          while (index[andim-1] < blocks[andim-1]) {
 
1105
            /* find bounding coordinates of block */
 
1106
            chk = 1;
 
1107
            for (i = 0; i < andim; i++) {
 
1108
              lo[i] = index[i]*block_dims[i]+1;
 
1109
              hi[i] = (index[i] + 1)*block_dims[i];
 
1110
              if (hi[i] > adims[i]) hi[i] = adims[i];
 
1111
              if (hi[i] < lo[i]) chk = 0;
 
1112
            }
 
1113
            if (chk) {
 
1114
              pnga_access_block_grid_ptr(g_a, index, &ptr_a, ld);
 
1115
              nrow   = hi[0] -lo[0]+1;
 
1116
              ncol   = hi[1] -lo[1]+1; 
 
1117
              nelem  = nrow*ncol;
 
1118
              lob[0] = lo[1]; lob[1] = lo[0];
 
1119
              hib[0] = hi[1]; hib[1] = hi[0];
 
1120
              for(i = 0; i < ncol; i++){
 
1121
                char *ptr = ptr_tmp + i*size;
 
1122
                snga_local_transpose(atype, ptr_a, nrow, block_dims[0]*size, ptr);
 
1123
                ptr_a += ld[0]*size;
 
1124
              }
 
1125
              pnga_put(g_b, lob, hib, ptr_tmp ,&block_dims[0]);
 
1126
              pnga_release_update_block_grid(g_a, index);
 
1127
            }
 
1128
            /* increment index to get next block on processor */
 
1129
            index[0] += topology[0];
 
1130
            for (i = 0; i < andim; i++) {
 
1131
              if (index[i] >= blocks[i] && i<andim-1) {
 
1132
                index[i] = proc_index[i];
 
1133
                index[i+1] += topology[i+1];
 
1134
              }
 
1135
            }
 
1136
          }
 
1137
        }
 
1138
      }
 
1139
      ga_free(ptr_tmp);
 
1140
    }
 
1141
 
 
1142
    if(local_sync_end)pnga_sync();
 
1143
    GA_POP_NAME;
 
1144
}