Context.C 25 KB
Newer Older
Francois Gygi committed
1 2
////////////////////////////////////////////////////////////////////////////////
//
Francois Gygi committed
3 4 5 6
// Copyright (c) 2008 The Regents of the University of California
//
// This file is part of Qbox
//
Francois Gygi committed
7 8
// Qbox is distributed under the terms of the GNU General Public License
// as published by the Free Software Foundation, either version 2 of
Francois Gygi committed
9 10 11 12 13 14
// the License, or (at your option) any later version.
// See the file COPYING in the root directory of this distribution
// or <http://www.gnu.org/licenses/>.
//
////////////////////////////////////////////////////////////////////////////////
//
Francois Gygi committed
15 16 17
// Context.C
//
////////////////////////////////////////////////////////////////////////////////
18
// $Id: Context.C,v 1.18 2009-11-30 02:33:49 fgygi Exp $
Francois Gygi committed
19

20
#include "Context.h"
Francois Gygi committed
21 22 23 24
#include <iostream>
#include <iomanip>
#include <cmath>
#include <cassert>
25
#include <vector>
Francois Gygi committed
26 27 28 29 30

#ifdef SCALAPACK
#include "blacs.h"
#endif

31
using namespace std;
Francois Gygi committed
32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73

#ifndef SCALAPACK
void Cblacs_pinfo(int *mypnum, int *nprocs)
{
    *mypnum = 0;
    *nprocs = 1;
    return;
}

void Cblacs_get(int icontxt, int what, int *val)
{
    *val=0;
    return;
}

void Cblacs_gridinit(int *icontxt, char order[], int nprow, int npcol)
{
    *icontxt=0;
    return;
}

void Cblacs_gridmap(int *icontxt, int *pmap, int ldpmap, int nprow, int npcol)
{
    pmap[0]=0;
    return;
}

void Cblacs_barrier(int icontxt, char* scope)
{
    return;
}

void Cblacs_abort(int icontxt, int errornum)
{
    return;
}

void Cblacs_gridexit(int icontxt)
{
    return;
}

74
void Cblacs_gridinfo(int icontxt, int *nprow, int *npcol,
Francois Gygi committed
75 76 77 78 79 80 81 82 83 84 85 86 87 88
                                  int *myprow, int *mypcol)
{
    *nprow=1; *npcol=1; *myprow=0; *mypcol=0; return;
}

int Cblacs_pnum(int icontxt, int prow, int pcol)
{ return 0;}

void Cdgesd2d(int icontxt,int m,int n,double *A,int lda,int rdest,int cdest)
{ return; }

void Cdgerv2d(int icontxt,int m,int n,double *A,int lda,int rdest,int cdest)
{ return; }

89
void Cdgsum2d(int icontxt, char* scope, char* top,
Francois Gygi committed
90 91 92 93 94
              int m, int n, double *a, int lda, int rdest, int cdest)
{ return; }

void Cdgamx2d(int icontxt, char scope[], char top[],int m,int n,double *A,
              int lda, int *ra, int *ca, int rcflag, int rdest, int cdest)
95
{ if ( rcflag > 0 ) { *ca = 0; *ra = 0; } return; }
Francois Gygi committed
96 97 98

void Cdgamn2d(int icontxt, char scope[], char top[],int m,int n,double *A,
              int lda, int *ra, int *ca, int rcflag, int rdest, int cdest)
99
{ if ( rcflag > 0 ) { *ca = 0; *ra = 0; } return; }
Francois Gygi committed
100 101 102 103 104 105 106 107 108 109 110 111 112 113 114

void Cdgebs2d(int icontxt, char scope[], char top[],int m,int n,double *A,
              int lda)
{ return; }

void Cdgebr2d(int icontxt, char scope[], char top[],int m,int n,double *A,
              int lda, int rsrc, int csrc)
{ return; }

void Cigesd2d(int icontxt,int m,int n,int *A,int lda,int rdest,int cdest)
{ return; }

void Cigerv2d(int icontxt,int m,int n,int *A,int lda,int rdest,int cdest)
{ return; }

115
void Cigsum2d(int icontxt, char* scope, char* top,
Francois Gygi committed
116 117 118
              int m, int n, int *a, int lda, int rdest, int cdest)
{ return; }

119
void Cigamx2d(int icontxt, char scope[], char top[],int m,int n,int *A,int lda,
Francois Gygi committed
120
              int *ra, int *ca, int rcflag, int rdest, int cdest)
121
{ if ( rcflag > 0 ) { *ca = 0; *ra = 0; } return; }
Francois Gygi committed
122

123
void Cigamn2d(int icontxt, char scope[], char top[],int m,int n,int *A,int lda,
Francois Gygi committed
124
              int *ra, int *ca, int rcflag, int rdest, int cdest)
125
{ if ( rcflag > 0 ) { *ca = 0; *ra = 0; } return; }
Francois Gygi committed
126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144

void Cigebs2d(int icontxt, char scope[], char top[],int m,int n,int *A,
              int lda)
{ return; }

void Cigebr2d(int icontxt, char scope[], char top[],int m,int n,int *A,
              int lda, int rsrc, int csrc)
{ return; }

#endif

#if USE_MPI
#include <mpi.h>
#else
typedef int MPI_Comm;
const int MPI_COMM_NULL = 0;
#endif

struct ContextRep
145
{
Francois Gygi committed
146
  private:
147

Francois Gygi committed
148 149 150 151 152 153 154 155 156 157
  int ictxt_;
  int myrow_;
  int mycol_;
  int nprow_;
  int npcol_;
  int size_;
  int myproc_;
  int mype_;
  bool onpe0_;
  bool active_;
158

Francois Gygi committed
159 160
  vector<int> pmap_;
  MPI_Comm comm_;
161

Francois Gygi committed
162 163 164 165 166
  // keep assignment and copy constructors private
  ContextRep& operator=(const Context& c);
  ContextRep(const Context& c);

  public:
167

Francois Gygi committed
168 169 170 171 172
  int ictxt() const { return ictxt_; }
  int myrow() const { return myrow_; }
  int mycol() const { return mycol_; }
  int nprow() const { return nprow_; }
  int npcol() const { return npcol_; }
173

Francois Gygi committed
174 175 176 177 178 179 180
  // number of processes in the context
  // returns -1 if current process is not part of this context
  int size() const { return size_; }
  // position of current process in row-major order
  // returns -1 if current process is not part of this context
  int myproc() const { return myproc_; }
  int mype() const { return mype_; }
181
  int pmap(int irow, int icol) const { return pmap_[irow+nprow_*icol]; }
182

Francois Gygi committed
183 184 185 186 187
  bool onpe0(void) const { return onpe0_; }
  bool active(void) const { return active_; }
  void abort(int ierr) const { Cblacs_abort(ictxt_,ierr); }
  void barrier(void) const { Cblacs_barrier(ictxt_,"A"); }
  void barrier(char scope) const { Cblacs_barrier(ictxt_,&scope); }
188

Francois Gygi committed
189 190 191 192
  void dsend(int m, int n, double* a, int lda, int rdest, int cdest) const
  {
    Cdgesd2d(ictxt_,m,n,a,lda,rdest,cdest);
  }
193

Francois Gygi committed
194 195 196 197
  void drecv(int m, int n, double* a, int lda, int rsrc, int csrc) const
  {
    Cdgerv2d(ictxt_,m,n,a,lda,rsrc,csrc);
  }
198

Francois Gygi committed
199 200 201 202 203
  void dsum(char scope, char topology, int m, int n, double* a, int lda,
            int rdest, int cdest) const
  {
    Cdgsum2d(ictxt_,&scope,&topology,m,n,a,lda,rdest,cdest);
  }
204

Francois Gygi committed
205 206 207 208 209
  void dmax(char scope, char topology, int m, int n, double* a, int lda,
            int rdest, int cdest) const
  {
    Cdgamx2d(ictxt_,&scope,&topology,m,n,a,lda,(int*)0,(int*)0,-1,rdest,cdest);
  }
210

Francois Gygi committed
211 212 213 214 215
  void dmax(char scope, char topology, int m, int n, double* a, int lda,
            int* ra, int* ca, int rcflag, int rdest, int cdest) const
  {
    Cdgamx2d(ictxt_,&scope,&topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest);
  }
216

Francois Gygi committed
217 218 219 220 221
  void dmin(char scope, char topology, int m, int n, double* a, int lda,
            int* ra, int* ca, int rcflag, int rdest, int cdest) const
  {
    Cdgamn2d(ictxt_,&scope,&topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest);
  }
222

Francois Gygi committed
223 224 225 226 227
  void dmin(char scope, char topology, int m, int n, double* a, int lda,
            int rdest, int cdest) const
  {
    Cdgamn2d(ictxt_,&scope,&topology,m,n,a,lda,(int*)0,(int*)0,-1,rdest,cdest);
  }
228 229

  void dbcast_send(char scope, char topology,
Francois Gygi committed
230 231 232 233
                   int m, int n, double* a,int lda) const
  {
    Cdgebs2d(ictxt_,&scope,&topology,m,n,a,lda);
  }
234

Francois Gygi committed
235 236 237 238 239
  void dbcast_recv(char scope, char topology, int m, int n, double* a, int lda,
               int rsrc, int csrc) const
  {
    Cdgebr2d(ictxt_,&scope,&topology,m,n,a,lda,rsrc,csrc);
  }
240

Francois Gygi committed
241 242 243 244
  void isend(int m, int n, int* a, int lda, int rdest, int cdest) const
  {
    Cigesd2d(ictxt_,m,n,a,lda,rdest,cdest);
  }
245

Francois Gygi committed
246 247 248 249
  void irecv(int m, int n, int* a, int lda, int rsrc, int csrc) const
  {
    Cigerv2d(ictxt_,m,n,a,lda,rsrc,csrc);
  }
250

Francois Gygi committed
251 252 253 254 255
  void isum(char scope, char topology, int m, int n, int* a, int lda,
            int rdest, int cdest) const
  {
    Cigsum2d(ictxt_,&scope,&topology,m,n,a,lda,rdest,cdest);
  }
256

Francois Gygi committed
257 258 259 260 261
  void imax(char scope, char topology, int m, int n, int* a, int lda,
            int* ra, int* ca, int rcflag, int rdest, int cdest) const
  {
    Cigamx2d(ictxt_,&scope,&topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest);
  }
262

Francois Gygi committed
263 264 265 266 267
  void imax(char scope, char topology, int m, int n, int* a, int lda,
            int rdest, int cdest) const
  {
    Cigamx2d(ictxt_,&scope,&topology,m,n,a,lda,(int*)0,(int*)0,-1,rdest,cdest);
  }
268

Francois Gygi committed
269 270 271 272 273
  void imin(char scope, char topology, int m, int n, int* a, int lda,
            int* ra, int* ca, int rcflag, int rdest, int cdest) const
  {
    Cigamn2d(ictxt_,&scope,&topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest);
  }
274

Francois Gygi committed
275 276 277 278 279
  void imin(char scope, char topology, int m, int n, int* a, int lda,
            int rdest, int cdest) const
  {
    Cigamn2d(ictxt_,&scope,&topology,m,n,a,lda,(int*)0,(int*)0,-1,rdest,cdest);
  }
280 281

  void ibcast_send(char scope, char topology,
Francois Gygi committed
282 283 284 285
                   int m, int n, int* a,int lda) const
  {
    Cigebs2d(ictxt_,&scope,&topology,m,n,a,lda);
  }
286

Francois Gygi committed
287 288 289 290 291
  void ibcast_recv(char scope, char topology, int m, int n, int* a, int lda,
               int rsrc, int csrc) const
  {
    Cigebr2d(ictxt_,&scope,&topology,m,n,a,lda,rsrc,csrc);
  }
292

Francois Gygi committed
293 294
  void string_send(string& s, int rdest, int cdest) const
  {
295
#if USE_MPI
Francois Gygi committed
296 297 298 299 300 301
    int len = s.size();
    isend(1,1,&len,1,rdest,cdest);
    int ilen = len/(sizeof(int)/sizeof(char));
    if ( len%(sizeof(int)/sizeof(char)) != 0 ) ilen++;
    int* ibuf = new int[ilen];
    s.copy((char*)ibuf,string::npos);
302
    isend(ilen,1,ibuf,ilen,rdest,cdest);
Francois Gygi committed
303
    delete [] ibuf;
304
#endif
Francois Gygi committed
305
  }
306

Francois Gygi committed
307 308
  void string_recv(string& s, int rsrc, int csrc) const
  {
309
#if USE_MPI
Francois Gygi committed
310 311 312 313 314
    int len = -1;
    irecv(1,1,&len,1,rsrc,csrc);
    int ilen = len/(sizeof(int)/sizeof(char));
    if ( len%(sizeof(int)/sizeof(char)) != 0 ) ilen++;
    int* ibuf = new int[ilen];
315
    irecv(ilen,1,ibuf,ilen,rsrc,csrc);
Francois Gygi committed
316 317 318
    s.resize(len);
    s.assign((char*)ibuf,len);
    delete [] ibuf;
319
#endif
Francois Gygi committed
320
  }
321

Francois Gygi committed
322 323
  void string_bcast(string& s, int isrc) const
  {
324
#if USE_MPI
Francois Gygi committed
325 326 327 328 329 330 331
    int len;
    if ( mype() == isrc )
    {
      len = s.length();
    }
    MPI_Bcast(&len,1,MPI_INT,isrc,comm());
    char* buf = new char[len+1];
Francois Gygi committed
332 333 334 335 336 337 338
    // s is initialized only on task isrc
    if ( mype() == isrc )
    {
      s.copy(buf,string::npos);
      buf[len]=0;
      assert(buf[len]=='\0');
    }
Francois Gygi committed
339 340 341
    MPI_Bcast(buf,len+1,MPI_CHAR,isrc,comm());
    s = buf;
    delete [] buf;
342
#endif
Francois Gygi committed
343
  }
344

Francois Gygi committed
345 346
  bool operator==(const ContextRep& c) const
  { return (ictxt_==c.ictxt());}
347

Francois Gygi committed
348 349 350
  // MPI communicator for this context. Returns MPI_COMM_NULL if
  // this process is not part of the context
  MPI_Comm comm(void) const { return comm_; }
351

Francois Gygi committed
352 353 354 355
  // Constructors

  // default global context: construct a single-row global ContextRep
  explicit ContextRep();
356

357 358
  // global ContextRep of size nprow * npcol with column major order
  explicit ContextRep(int nprow, int npcol);
359

360
  // construct a ContextRep of size nprow*npcol using the processes
361 362
  // in context c starting at process (irstart,icstart)
  explicit ContextRep(const ContextRep &c, int nprow, int npcol,
363
    int irstart, int icstart);
364

Francois Gygi committed
365
  ~ContextRep();
366

Francois Gygi committed
367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385
  void print(ostream& os) const;
};

////////////////////////////////////////////////////////////////////////////////
ContextRep::ContextRep() : ictxt_(-1), myrow_(-1), mycol_(-1)
{
  // construct a single row global context
  int nprocs;
  char order='R';
  Cblacs_pinfo( &mype_, &nprocs );
  nprow_ = 1;
  npcol_ = nprocs;

  Cblacs_get(0, 0, &ictxt_ );
  Cblacs_gridinit( &ictxt_, &order, nprow_, npcol_ );

  // get values of nprow_, npcol_, myrow_ and mycol_ in the new context
  if ( ictxt_ >= 0 )
    Cblacs_gridinfo(ictxt_, &nprow_, &npcol_, &myrow_, &mycol_);
386

Francois Gygi committed
387 388 389 390 391 392 393 394 395 396 397 398 399 400
  size_ = nprow_ * npcol_;
  myproc_ = myrow_ < 0 ? -1 : mycol_ + npcol_ * myrow_;
  onpe0_ = ( mype_ == 0 );
  active_ = ( ictxt_ >= 0 );

  pmap_.resize(size_);
  for ( int i = 0; i < size_; i++ )
    pmap_[i] = i;

#if USE_MPI
  MPI_Group group_world, subgroup;
  MPI_Comm_group(MPI_COMM_WORLD,&group_world);
  MPI_Group_incl(group_world,size_,&pmap_[0],&subgroup);
  MPI_Comm_create(MPI_COMM_WORLD,subgroup,&comm_);
401 402
  MPI_Group_free(&group_world);
  MPI_Group_free(&subgroup);
Francois Gygi committed
403 404 405 406 407 408
#else
  comm_ = 0;
#endif
}

////////////////////////////////////////////////////////////////////////////////
409
ContextRep::ContextRep(int nprow, int npcol) :
Francois Gygi committed
410 411
  ictxt_(-1), myrow_(-1), mycol_(-1), nprow_(nprow), npcol_(npcol)
{
412 413
  int nprocs;
  char order = 'C';
Francois Gygi committed
414 415 416
  Cblacs_pinfo( &mype_, &nprocs );
  if ( nprocs < nprow * npcol )
  {
417 418 419
    cout << " nprocs=" << nprocs << endl;
    cout << " Context nprow*npcol > nprocs" << endl;
    Cblacs_abort(ictxt_, 1);
Francois Gygi committed
420 421 422 423 424 425 426
  }
  Cblacs_get(0, 0, &ictxt_ );
  Cblacs_gridinit( &ictxt_, &order, nprow, npcol );

  // get values of nprow_, npcol_, myrow_ and mycol_ in the new context
  if ( ictxt_ >= 0 )
    Cblacs_gridinfo(ictxt_, &nprow_, &npcol_, &myrow_, &mycol_);
427

Francois Gygi committed
428
  size_ = nprow_ * npcol_;
429
  myproc_ = Cblacs_pnum(ictxt_,myrow_,mycol_);
Francois Gygi committed
430 431
  onpe0_ = ( mype_ == 0 );
  active_ = ( ictxt_ >= 0 );
432

Francois Gygi committed
433
  pmap_.resize(size_);
434 435 436 437 438 439 440 441
  // column-major order
  int i = 0;
  for ( int ic = 0; ic < npcol; ic++ )
    for ( int ir = 0; ir < nprow; ir++ )
    {
      pmap_[ir+nprow*ic] = i;
      i++;
    }
Francois Gygi committed
442 443 444 445 446 447

#if USE_MPI
  MPI_Group group_world, subgroup;
  MPI_Comm_group(MPI_COMM_WORLD,&group_world);
  MPI_Group_incl(group_world,size_,&pmap_[0],&subgroup);
  MPI_Comm_create(MPI_COMM_WORLD,subgroup,&comm_);
448 449
  MPI_Group_free(&group_world);
  MPI_Group_free(&subgroup);
Francois Gygi committed
450 451 452 453
#else
  comm_ = 0;
#endif
}
454

Francois Gygi committed
455
////////////////////////////////////////////////////////////////////////////////
456
ContextRep::ContextRep(const ContextRep& c, int nprow, int npcol,
457
  int irstart, int icstart) :
Francois Gygi committed
458 459 460
ictxt_(-1), myrow_(-1), mycol_(-1), nprow_(nprow), npcol_(npcol)
{
  assert(c.active());
461
  vector<int> gmap_;
462

Francois Gygi committed
463 464 465
  int nprocs;
  Cblacs_pinfo( &mype_, &nprocs );
  // construct a (nprow*npcol) context using the processes in c
466
  // located at (irstart,icstart)
467
  if ( irstart < 0 || icstart < 0 ||
468
       irstart+nprow > c.nprow() || icstart+npcol > c.npcol() )
Francois Gygi committed
469
  {
470
    cout << " Context::Context: cut rectangle: invalid parameters" << endl;
471
    cout << irstart << " " << icstart << " " << nprow << " " << npcol << endl;
Francois Gygi committed
472 473 474
    Cblacs_abort(ictxt_, 1);
  }
  pmap_.resize(nprow*npcol);
475
  gmap_.resize(nprow*npcol);
Francois Gygi committed
476 477
  // build pmap
  int i = 0;
478 479 480 481 482 483 484
  for ( int ic = icstart; ic < icstart+npcol; ic++ )
    for ( int ir = irstart; ir < irstart+nprow; ir++ )
    {
      pmap_[i] = c.pmap(ir,ic);
      gmap_[i] = ic + c.npcol()*ir;
      i++;
    }
485

Francois Gygi committed
486
  Cblacs_get(c.ictxt(), 10, &ictxt_ );
487
  Cblacs_gridmap(&ictxt_,&gmap_[0],nprow,nprow,npcol);
488

Francois Gygi committed
489 490 491
  // get values of nprow_, npcol_, myrow_ and mycol_ in the new context
  if ( ictxt_ >= 0 )
    Cblacs_gridinfo(ictxt_, &nprow_, &npcol_, &myrow_, &mycol_);
492

Francois Gygi committed
493 494 495 496
  size_ = nprow_ * npcol_;
  myproc_ = myrow_ < 0 ? -1 : mycol_ + npcol_ * myrow_;
  onpe0_ = ( mype_ == 0 );
  active_ = ( ictxt_ >= 0 );
497

Francois Gygi committed
498 499 500 501 502
#if USE_MPI
  MPI_Group c_group, subgroup;
  MPI_Comm_group(c.comm(),&c_group);
  MPI_Group_incl(c_group,size_,&pmap_[0],&subgroup);
  MPI_Comm_create(c.comm(),subgroup,&comm_);
503 504
  MPI_Group_free(&c_group);
  MPI_Group_free(&subgroup);
Francois Gygi committed
505 506 507 508 509 510 511 512
#else
  comm_ = 0;
#endif
}

////////////////////////////////////////////////////////////////////////////////
ContextRep::~ContextRep()
{
513
  if ( myrow_ != -1 )
Francois Gygi committed
514 515 516
  {
    // cout << " ContextRep destructor called on ictxt = " << ictxt_ << endl;
    Cblacs_gridexit( ictxt_ );
517
#if USE_MPI
Francois Gygi committed
518
    MPI_Comm_free(&comm_);
519
#endif
Francois Gygi committed
520 521 522 523 524 525 526 527 528 529 530 531 532 533 534
  }
}

////////////////////////////////////////////////////////////////////////////////
void ContextRep::print(ostream& os) const
{
  if ( active_ )
  {
    os << " " << nprow_ << "x" << npcol_ << " ";
    os << "{ ";
    for ( int ir = 0; ir < nprow_; ir++ )
    {
      os << "{ ";
      for ( int ic = 0; ic < npcol_; ic++ )
      {
535
        os << pmap_[ir+nprow_*ic] << " ";
Francois Gygi committed
536 537 538
      }
      os << "} ";
    }
539 540 541 542 543 544
    //os << "} (ictxt=" << ictxt_ << ")" << endl;
    os << "} (ictxt=" << ictxt_ << ")";
    int handle;
    Cblacs_get(ictxt_,10,&handle);
    os << "(handle=" << handle << ")";
    os << "(comm=" << comm_<< ")" << endl;
Francois Gygi committed
545 546 547 548 549 550 551 552 553 554 555
  }
  else
  {
    os << " (inactive)" << endl;
  }
}

////////////////////////////////////////////////////////////////////////////////
class ContextImpl
{
  private:
556

Francois Gygi committed
557 558
  ContextRep* rep;
  int* pcount;
559

Francois Gygi committed
560
  public:
561

Francois Gygi committed
562 563 564 565 566
  ContextRep* operator->() { return rep; }
  ContextRep* get_rep(void) { return rep; }
  ContextImpl(ContextRep* pp) : rep(pp), pcount(new int(1)) {}
  ContextImpl(const ContextImpl& r) : rep(r.rep), pcount(r.pcount)
  { (*pcount)++; }
567

Francois Gygi committed
568 569 570 571 572 573 574 575 576 577 578 579 580
  ContextImpl& operator=(const ContextImpl& r)
  {
    if ( rep == r.rep ) return *this;
    if ( --(*pcount) == 0 )
    {
      delete rep;
      delete pcount;
    }
    rep = r.rep;
    pcount = r.pcount;
    (*pcount)++;
    return *this;
  }
581

Francois Gygi committed
582 583 584 585 586 587 588 589 590 591 592 593 594 595
  ~ContextImpl(void)
  {
    if ( --(*pcount) == 0 )
    {
      delete rep;
      delete pcount;
    }
  }
};

////////////////////////////////////////////////////////////////////////////////
Context::Context(void) : pimpl_(new ContextImpl(new ContextRep())) {}

////////////////////////////////////////////////////////////////////////////////
596
Context::Context(int nprow, int npcol):
597
pimpl_(new ContextImpl(new ContextRep(nprow,npcol))) {}
598

Francois Gygi committed
599
////////////////////////////////////////////////////////////////////////////////
600 601
Context::Context(const Context &c, int nprow, int npcol,
        int irstart, int icstart) :
602 603
pimpl_(new ContextImpl(new ContextRep(*(c.pimpl_)->get_rep(),nprow,npcol,
irstart,icstart))) {}
604

Francois Gygi committed
605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662
////////////////////////////////////////////////////////////////////////////////
Context::~Context() { delete pimpl_; }

////////////////////////////////////////////////////////////////////////////////
Context::Context(const Context& c) : pimpl_(new ContextImpl(*(c.pimpl_))) {}

////////////////////////////////////////////////////////////////////////////////
Context& Context::operator=(const Context& rhs)
{
  if ( &rhs == this ) return *this;
  *pimpl_ = *(rhs.pimpl_);
  return *this;
}

////////////////////////////////////////////////////////////////////////////////
int Context::ictxt() const { return (*pimpl_)->ictxt(); }

////////////////////////////////////////////////////////////////////////////////
int Context::myrow() const { return (*pimpl_)->myrow(); }

////////////////////////////////////////////////////////////////////////////////
int Context::mycol() const { return (*pimpl_)->mycol(); }

////////////////////////////////////////////////////////////////////////////////
int Context::nprow() const { return (*pimpl_)->nprow(); }

////////////////////////////////////////////////////////////////////////////////
int Context::npcol() const { return (*pimpl_)->npcol(); }

////////////////////////////////////////////////////////////////////////////////
int Context::size() const { return (*pimpl_)->size(); }

////////////////////////////////////////////////////////////////////////////////
int Context::myproc() const { return (*pimpl_)->myproc(); }

////////////////////////////////////////////////////////////////////////////////
int Context::mype() const { return (*pimpl_)->mype(); }

////////////////////////////////////////////////////////////////////////////////
int Context::pmap(int irow, int icol) const
{ return (*pimpl_)->pmap(irow,icol); }

////////////////////////////////////////////////////////////////////////////////
bool Context::onpe0(void) const { return (*pimpl_)->onpe0(); }

////////////////////////////////////////////////////////////////////////////////
bool Context::active(void) const { return (*pimpl_)->active(); }

////////////////////////////////////////////////////////////////////////////////
void Context::abort(int ierr) const { (*pimpl_)->abort(ierr); }

////////////////////////////////////////////////////////////////////////////////
void Context::barrier(void) const { (*pimpl_)->barrier(); }

////////////////////////////////////////////////////////////////////////////////
void Context::barrier(char scope) const { (*pimpl_)->barrier(scope); }

////////////////////////////////////////////////////////////////////////////////
663
void Context::dsend(int m, int n, double* a,
Francois Gygi committed
664 665 666
  int lda, int rdest, int cdest) const
{ (*pimpl_)->dsend(m,n,a,lda,rdest,cdest); }

667
void Context::drecv(int m, int n, double* a,
Francois Gygi committed
668 669 670
  int lda, int rsrc, int csrc) const
{ (*pimpl_)->drecv(m,n,a,lda,rsrc,csrc); }

671
void Context::dsum(char scope, char topology,
Francois Gygi committed
672 673 674 675 676 677 678 679 680
  int m, int n, double* a, int lda, int rdest, int cdest) const
{ (*pimpl_)->dsum(scope,topology,m,n,a,lda,rdest,cdest); }

void Context::dsum(char scope, int m, int n, double* a, int lda) const
{ (*pimpl_)->dsum(scope,' ',m,n,a,lda,-1,-1); }

void Context::dsum(int m, int n, double* a, int lda) const
{ (*pimpl_)->dsum('A',' ',m,n,a,lda,-1,-1); }

681 682
void Context::dmax(char scope, char topology,
  int m, int n, double* a, int lda, int* ra, int* ca, int rcflag,
Francois Gygi committed
683 684 685
  int rdest, int cdest) const
{ (*pimpl_)->dmax(scope,topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest); }

686
void Context::dmax(char scope, char topology,
Francois Gygi committed
687 688 689 690 691 692 693 694 695
  int m, int n, double* a, int lda, int rdest, int cdest) const
{ (*pimpl_)->dmax(scope,topology,m,n,a,lda,rdest,cdest); }

void Context::dmax(char scope, int m, int n, double* a, int lda) const
{ (*pimpl_)->dmax(scope,' ',m,n,a,lda,-1,-1); }

void Context::dmax(int m, int n, double* a, int lda) const
{ (*pimpl_)->dmax('A',' ',m,n,a,lda,-1,-1); }

696 697
void Context::dmin(char scope, char topology,
  int m, int n, double* a, int lda, int* ra, int* ca, int rcflag,
Francois Gygi committed
698 699 700
  int rdest, int cdest) const
{ (*pimpl_)->dmin(scope,topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest); }

701
void Context::dmin(char scope, char topology,
Francois Gygi committed
702 703 704 705 706 707 708 709 710
  int m, int n, double* a, int lda, int rdest, int cdest) const
{ (*pimpl_)->dmin(scope,topology,m,n,a,lda,rdest,cdest); }

void Context::dmin(char scope, int m, int n, double* a, int lda) const
{ (*pimpl_)->dmin(scope,' ',m,n,a,lda,-1,-1); }

void Context::dmin(int m, int n, double* a, int lda) const
{ (*pimpl_)->dmin('A',' ',m,n,a,lda,-1,-1); }

711
void Context::dbcast_send(char scope, char topology,
Francois Gygi committed
712 713 714 715 716 717 718 719 720
  int m, int n, double* a, int lda) const
{ (*pimpl_)->dbcast_send(scope,topology,m,n,a,lda); }

void Context::dbcast_send(char scope, int m, int n, double* a, int lda) const
{ (*pimpl_)->dbcast_send(scope,' ',m,n,a,lda); }

void Context::dbcast_send(int m, int n, double* a, int lda) const
{ (*pimpl_)->dbcast_send('A',' ',m,n,a,lda); }

721
void Context::dbcast_recv(char scope, char topology,
Francois Gygi committed
722 723 724
  int m, int n, double* a, int lda, int rsrc,int csrc) const
{ (*pimpl_)->dbcast_recv(scope,topology,m,n,a,lda,rsrc,csrc); }

725
void Context::dbcast_recv(char scope,
Francois Gygi committed
726 727 728
  int m, int n, double* a, int lda, int rsrc, int csrc) const
{ (*pimpl_)->dbcast_recv(scope,' ',m,n,a,lda,rsrc,csrc); }

729
void Context::dbcast_recv(int m, int n, double* a, int lda,
Francois Gygi committed
730 731 732 733
  int rsrc,int csrc) const
{ (*pimpl_)->dbcast_recv('A',' ',m,n,a,lda,rsrc,csrc); }

////////////////////////////////////////////////////////////////////////////////
734
void Context::isend(int m, int n, int* a,
Francois Gygi committed
735 736 737
  int lda, int rdest, int cdest) const
{ (*pimpl_)->isend(m,n,a,lda,rdest,cdest); }

738
void Context::irecv(int m, int n, int* a,
Francois Gygi committed
739 740 741
  int lda, int rsrc, int csrc) const
{ (*pimpl_)->irecv(m,n,a,lda,rsrc,csrc); }

742
void Context::isum(char scope, char topology,
Francois Gygi committed
743 744 745 746 747 748 749 750 751
  int m, int n, int* a, int lda, int rdest, int cdest) const
{ (*pimpl_)->isum(scope,topology,m,n,a,lda,rdest,cdest); }

void Context::isum(char scope, int m, int n, int* a, int lda) const
{ (*pimpl_)->isum(scope,' ',m,n,a,lda,-1,-1); }

void Context::isum(int m, int n, int* a, int lda) const
{ (*pimpl_)->isum('A',' ',m,n,a,lda,-1,-1); }

752 753
void Context::imax(char scope, char topology,
  int m, int n, int* a, int lda, int* ra, int* ca, int rcflag,
Francois Gygi committed
754 755 756
  int rdest, int cdest) const
{ (*pimpl_)->imax(scope,topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest); }

757
void Context::imax(char scope, char topology,
Francois Gygi committed
758 759 760 761 762 763 764 765 766
  int m, int n, int* a, int lda, int rdest, int cdest) const
{ (*pimpl_)->imax(scope,topology,m,n,a,lda,rdest,cdest); }

void Context::imax(char scope, int m, int n, int* a, int lda) const
{ (*pimpl_)->imax(scope,' ',m,n,a,lda,-1,-1); }

void Context::imax(int m, int n, int* a, int lda) const
{ (*pimpl_)->imax('A',' ',m,n,a,lda,-1,-1); }

767 768
void Context::imin(char scope, char topology,
  int m, int n, int* a, int lda, int* ra, int* ca, int rcflag,
Francois Gygi committed
769 770 771
  int rdest, int cdest) const
{ (*pimpl_)->imin(scope,topology,m,n,a,lda,ra,ca,rcflag,rdest,cdest); }

772
void Context::imin(char scope, char topology,
Francois Gygi committed
773 774 775 776 777 778 779 780 781
  int m, int n, int* a, int lda, int rdest, int cdest) const
{ (*pimpl_)->imin(scope,topology,m,n,a,lda,rdest,cdest); }

void Context::imin(char scope, int m, int n, int* a, int lda) const
{ (*pimpl_)->imin(scope,' ',m,n,a,lda,-1,-1); }

void Context::imin(int m, int n, int* a, int lda) const
{ (*pimpl_)->imin('A',' ',m,n,a,lda,-1,-1); }

782
void Context::ibcast_send(char scope, char topology,
Francois Gygi committed
783 784 785 786 787 788 789 790 791
  int m, int n, int* a, int lda) const
{ (*pimpl_)->ibcast_send(scope,topology,m,n,a,lda); }

void Context::ibcast_send(char scope, int m, int n, int* a, int lda) const
{ (*pimpl_)->ibcast_send(scope,' ',m,n,a,lda); }

void Context::ibcast_send(int m, int n, int* a, int lda) const
{ (*pimpl_)->ibcast_send('A',' ',m,n,a,lda); }

792
void Context::ibcast_recv(char scope, char topology,
Francois Gygi committed
793 794 795
  int m, int n, int* a, int lda, int rsrc,int csrc) const
{ (*pimpl_)->ibcast_recv(scope,topology,m,n,a,lda,rsrc,csrc); }

796
void Context::ibcast_recv(char scope,
Francois Gygi committed
797 798 799
  int m, int n, int* a, int lda, int rsrc, int csrc) const
{ (*pimpl_)->ibcast_recv(scope,' ',m,n,a,lda,rsrc,csrc); }

800
void Context::ibcast_recv(int m, int n, int* a, int lda,
Francois Gygi committed
801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818
  int rsrc,int csrc) const
{ (*pimpl_)->ibcast_recv('A',' ',m,n,a,lda,rsrc,csrc); }

void Context::string_send(string& s, int rdest, int cdest) const
{ (*pimpl_)->string_send(s,rdest,cdest); }

void Context::string_recv(string& s, int rsrc, int csrc) const
{ (*pimpl_)->string_recv(s,rsrc,csrc); }

void Context::string_bcast(string& s, int isrc) const
{ (*pimpl_)->string_bcast(s,isrc); }

////////////////////////////////////////////////////////////////////////////////
bool Context::operator==(const Context& ctxt) const
{ return ( (*pimpl_)->ictxt() == ctxt.ictxt() ); }

////////////////////////////////////////////////////////////////////////////////
MPI_Comm Context::comm(void) const { return (*pimpl_)->comm(); }
819

Francois Gygi committed
820 821 822 823 824 825 826 827 828 829
////////////////////////////////////////////////////////////////////////////////
void Context::print(ostream& os) const { (*pimpl_)->print(os);}

////////////////////////////////////////////////////////////////////////////////
ostream& operator<<(ostream& os, const Context& c)
{
  c.print(os);
  return os;
}