test/EpetraBenchmarkTest/cxx_main.cpp

Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //               Epetra: Linear Algebra Services Package 
00005 //                 Copyright (2001) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ************************************************************************
00027 //@HEADER
00028 #include <iostream>
00029 #include <cstring>
00030 #include "Epetra_Map.h"
00031 #include "Epetra_LocalMap.h"
00032 #include "Epetra_BlockMap.h"
00033 #include "Epetra_Time.h"
00034 #include "Epetra_CrsMatrix.h"
00035 #include "Epetra_VbrMatrix.h"
00036 #include "Epetra_Vector.h"
00037 #include "Epetra_IntVector.h"
00038 #include "Epetra_MultiVector.h"
00039 #include "Epetra_IntSerialDenseVector.h"
00040 #include "Epetra_SerialDenseVector.h"
00041 #ifdef EPETRA_MPI
00042 #include "Epetra_MpiComm.h"
00043 #include "mpi.h"
00044 #else
00045 #include "Epetra_SerialComm.h"
00046 #endif
00047 #include "Epetra_Version.h"
00048 
00049 // prototypes
00050 
00051 void GenerateCrsProblem(int numNodesX, int numNodesY, int numProcsX, int numProcsY, int numPoints, 
00052       int * xoff, int * yoff,
00053       const Epetra_Comm  &comm, 
00054       Epetra_Map *& map, 
00055       Epetra_CrsMatrix *& A, 
00056       Epetra_Vector *& b, 
00057       Epetra_Vector *& bt,
00058       Epetra_Vector *&xexact, bool StaticProfile, bool MakeLocalOnly);
00059 
00060 void GenerateCrsProblem(int numNodesX, int numNodesY, int numProcsX, int numProcsY, int numPoints, 
00061       int * xoff, int * yoff, int nrhs,
00062       const Epetra_Comm  &comm, 
00063       Epetra_Map *& map, 
00064       Epetra_CrsMatrix *& A, 
00065       Epetra_MultiVector *& b, 
00066       Epetra_MultiVector *& bt,
00067       Epetra_MultiVector *&xexact, bool StaticProfile, bool MakeLocalOnly);
00068  
00069 
00070 void GenerateMyGlobalElements(int numNodesX, int numNodesY, int numProcsX, int numProcs,
00071             int myPID, int * & myGlobalElements);
00072 
00073 void runMatrixTests(Epetra_CrsMatrix * A,  Epetra_MultiVector * b, Epetra_MultiVector * bt,
00074         Epetra_MultiVector * xexact, bool StaticProfile);
00075 
00076 int main(int argc, char *argv[])
00077 {
00078 
00079   const string EpetraBenchmarkTest = "Epetra Benchmark Test Version 0.2 08/30/2007"; // Change this number and date when changes are made to functionality
00080   int ierr = 0;
00081   double elapsed_time;
00082   double total_flops;
00083   double MFLOPs;
00084   double global_dimension;
00085   double global_nonzero_count;
00086     
00087 
00088 #ifdef EPETRA_MPI
00089 
00090   // Initialize MPI
00091   MPI_Init(&argc,&argv);
00092   Epetra_MpiComm comm( MPI_COMM_WORLD );
00093 #else
00094   Epetra_SerialComm comm;
00095 #endif
00096 
00097   // Check if we should print only selected timing info to standard out
00098   bool printflops = true, printmflops = true, printtime = true;
00099   if (argc>5) if (argv[5][0]=='-' && argv[5][1]=='f') {printmflops = false; printtime = false;}
00100   bool mflopsonly = false;
00101   if (argc>5) if (argv[5][0]=='-' && argv[5][1]=='m') {printflops = false; printtime = false;}
00102   bool timeonly = false;
00103   if (argc>5) if (argv[5][0]=='-' && argv[5][1]=='t') {printflops = false; printmflops = false;}
00104 
00105   // Check if we should print header to standard out
00106   bool makeheader = false;
00107   if (argc>6) if (argv[6][0]=='-' && argv[6][1]=='h') makeheader = true;
00108 
00109   if(argc < 5) {
00110     cerr << "Usage: " << argv[0]
00111          << " NumNodesX NumNodesY NumProcX NumProcY [-a|-f|-m|-t [-h]]" << endl
00112          << "where:" << endl
00113          << "NumNodesX         - Number of mesh nodes in X direction per processor" << endl
00114          << "NumNodesY         - Number of mesh nodes in Y direction per processor" << endl
00115          << "NumProcX          - Number of processors to use in X direction" << endl
00116          << "NumProcY          - Number of processors to use in Y direction" << endl
00117    << "-a|-f|-m|-t       - Type of information to print: a=all, f=FLOPS, m=MFLOP/s, t=time (sec).  Default is -a."
00118          << "-h                - (Optional) Printer output table header if -h present (typically done on first call)." << endl
00119          << " NOTES: NumProcX*NumProcY must equal the number of processors used to run the problem." << endl << endl
00120    << " Serial example:" << endl << endl
00121          << argv[0] << " 200 300 1 1 -m" << endl
00122    << " Run this program on 1 processor using a 200 by 300 grid, printing only MFLOP/s information without a header."<< endl <<endl
00123    << " MPI example:" << endl << endl
00124          << "mpirun -np 32 " << argv[0] << " 250 200 4 8 -a -h" << endl
00125    << " Run this program on 32 processors putting a 250 by 200 subgrid on each processor using 4 processors "<< endl
00126    << " in the X direction and 8 in the Y direction.  Total grid size is 1000 points in X and 1600 in Y for a total of 1.6M equations."<< endl
00127    << " Print all information. Print header." << endl
00128          << endl;
00129     return(1);
00130 
00131   }
00132     //char tmp;
00133     //if (comm.MyPID()==0) cout << "Press any key to continue..."<< endl;
00134     //if (comm.MyPID()==0) cin >> tmp;
00135     //comm.Barrier();
00136 
00137   if (makeheader && comm.MyPID()==0)
00138     cout << EpetraBenchmarkTest << endl
00139    << "Using " << Epetra_Version() << endl << endl;
00140   if (makeheader) cout << comm <<endl;
00141 
00142 
00143   // Redefine makeheader to only print on PE 0
00144 
00145   if (makeheader && comm.MyPID()!=0) makeheader = false;
00146 
00147   int numNodesX = atoi(argv[1]);
00148   int numNodesY = atoi(argv[2]);
00149   int numProcsX = atoi(argv[3]);
00150   int numProcsY = atoi(argv[4]);
00151   int numPoints = 25;
00152 
00153   if (makeheader) {
00154     cout << " Number of local nodes in X direction  = " << numNodesX << endl
00155    << " Number of local nodes in Y direction  = " << numNodesY << endl
00156    << " Number of global nodes in X direction = " << numNodesX*numProcsX << endl
00157    << " Number of global nodes in Y direction = " << numNodesY*numProcsY << endl
00158    << " Number of local nonzero entries       = " << numNodesX*numNodesY*numPoints << endl
00159    << " Number of global nonzero entries      = " << numNodesX*numNodesY*numPoints*numProcsX*numProcsY << endl
00160    << " Number of Processors in X direction   = " << numProcsX << endl
00161    << " Number of Processors in Y direction   = " << numProcsY << endl
00162    << " Number of Points in stencil           = " << numPoints << endl << endl;
00163     cout << " Timing the following:" <<endl
00164    << " SpMV - Sparse matrix vector product using Epetra_CrsMatrix class" << endl
00165    << " SpMM2- Sparse matrix times 2-column multivector using Epetra_CrsMatrix class" << endl
00166    << " SpMM4- Sparse matrix times 4-column multivector using Epetra_CrsMatrix class" << endl
00167    << " SpMM8- Sparse matrix times 8-column multivector using Epetra_CrsMatrix class" << endl
00168    << " 2-norm of an Epetra_MultiVector" << endl
00169    << " Dot-product of 2 Epetra_MultiVectors" << endl
00170    << " AXPY of 2 Epetra_MultiVectors" << endl << endl;
00171   }
00172 
00173   if (numProcsX*numProcsY!=comm.NumProc()) {
00174     cerr << "Number of processors = " << comm.NumProc() << endl
00175    << " is not the product of " << numProcsX << " and " << numProcsY << endl << endl;
00176     return(1);
00177   }
00178 
00179   if (numNodesX*numNodesY<=0) {
00180     cerr << "Product of number of nodes is <= zero" << endl << endl;
00181     return(1);
00182   }
00183 
00184   Epetra_IntSerialDenseVector Xoff, XLoff, XUoff;
00185   Epetra_IntSerialDenseVector Yoff, YLoff, YUoff;
00186   // Generate a 25-point 2D Finite Difference matrix
00187   Xoff.Size(25);
00188   Yoff.Size(25);
00189   int xi = 0, yi = 0;
00190   int xo = -2, yo = -2;
00191   Xoff[xi++] = xo++;  Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++;
00192   Yoff[yi++] = yo  ;  Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; 
00193   xo = -2, yo++;
00194   Xoff[xi++] = xo++;  Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++;
00195   Yoff[yi++] = yo  ;  Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; 
00196   xo = -2, yo++;
00197   Xoff[xi++] = xo++;  Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++;
00198   Yoff[yi++] = yo  ;  Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; 
00199   xo = -2, yo++;
00200   Xoff[xi++] = xo++;  Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++;
00201   Yoff[yi++] = yo  ;  Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; 
00202   xo = -2, yo++;
00203   Xoff[xi++] = xo++;  Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++; Xoff[xi++] = xo++;
00204   Yoff[yi++] = yo  ;  Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; Yoff[yi++] = yo  ; 
00205   
00206   
00207 
00208   Epetra_Map * map;
00209   Epetra_CrsMatrix * A;
00210   Epetra_MultiVector * b;
00211   Epetra_MultiVector * bt;
00212   Epetra_MultiVector * xexact;
00213   Epetra_SerialDenseVector resvec(0);
00214   global_dimension = (double) (numNodesX*numNodesY); // local grid dimension
00215   global_dimension *= (double) (numProcsX*numProcsY); // times number of processors
00216   global_nonzero_count = global_dimension * (double) numPoints; // number of nonzeros (w/o accounting for boundaries)
00217 
00218   //Timings
00219   Epetra_Time timer(comm);
00220   const int maxtimings = 7;
00221   const int maxtrials = 10;
00222   double results[maxtimings][3]; // timing results stored here
00223 
00224   int nrhs = 1;
00225   int ntimings = 0;
00226   for (int k=0; k<4; k++) {
00227     if (k>0) nrhs = nrhs*2;
00228       
00229     GenerateCrsProblem(numNodesX, numNodesY, numProcsX, numProcsY, numPoints,
00230            Xoff.Values(), Yoff.Values(), nrhs, comm,
00231            map, A, b, bt, xexact, true, false);
00232 
00233       
00234     Epetra_MultiVector z(*b);
00235     Epetra_MultiVector r(*b);
00236     Epetra_SerialDenseVector resvec(b->NumVectors());
00237     
00238     //Timings
00239     Epetra_Time timer(A->Comm());
00240     A->OptimizeStorage();
00241     
00242     timer.ResetStartTime();
00243     
00244     //maxtrials matvecs
00245     for( int i = 0; i < maxtrials; ++i )
00246       A->Multiply(false, *xexact, z); // Compute z = A*xexact or z = A'*xexact
00247     
00248     double elapsed_time = timer.ElapsedTime();
00249     double total_flops = 2.0*global_nonzero_count *((double) maxtrials);
00250     
00251     // Compute residual
00252     r.Update(-1.0, z, 1.0, *b, 0.0); // r = b - z
00253     
00254     r.Norm2(resvec.Values());
00255     double diff = resvec.NormInf();
00256     if (diff>1.0e-8 && comm.MyPID()==0) cerr << "Warning: Residual vector unusually large = " << diff << endl;
00257 
00258     double MFLOPs = total_flops/elapsed_time/1000000.0;
00259 
00260     results[ntimings][0] = total_flops;
00261     results[ntimings][1] = elapsed_time;
00262     results[ntimings++][2] = MFLOPs;
00263 
00264     delete A;
00265     delete b;
00266     delete bt; 
00267     delete xexact;
00268   } // end of k loop
00269 
00270 
00271   // *************** Vector timings *********************
00272 
00273   if (ntimings+3>maxtimings) cerr << "Variable maxtimings = " << maxtimings << " must be at least = " << ntimings+3 << endl;
00274 
00275   Epetra_MultiVector q(*map, nrhs);
00276   Epetra_MultiVector z(q);
00277   Epetra_MultiVector r(q);
00278   
00279   delete map;
00280 
00281   resvec.Resize(nrhs);
00282   
00283   
00284   timer.ResetStartTime();
00285   
00286   //maxtrials norms
00287   for( int i = 0; i < maxtrials; ++i )
00288     q.Norm2( resvec.Values() );
00289   
00290   elapsed_time = timer.ElapsedTime();
00291   total_flops = 2.0*global_dimension *((double) maxtrials);
00292   MFLOPs = total_flops/elapsed_time/1000000.0;
00293   results[ntimings][0] = total_flops;
00294   results[ntimings][1] = elapsed_time;
00295   results[ntimings++][2] = MFLOPs;
00296 
00297 
00298 
00299   timer.ResetStartTime();
00300   
00301   //maxtrials dot's
00302   for( int i = 0; i < maxtrials; ++i )
00303     q.Dot(z, resvec.Values());
00304   
00305   elapsed_time = timer.ElapsedTime();
00306   total_flops = 2.0*global_dimension *((double) maxtrials);
00307   MFLOPs = total_flops/elapsed_time/1000000.0;
00308   results[ntimings][0] = total_flops;
00309   results[ntimings][1] = elapsed_time;
00310   results[ntimings++][2] = MFLOPs;
00311   
00312   timer.ResetStartTime();
00313   
00314   //maxtrials updates
00315   for( int i = 0; i < maxtrials; ++i )
00316     q.Update(1.0, z, 1.0, r, 0.0);
00317   
00318   elapsed_time = timer.ElapsedTime();
00319   total_flops = 2.0*global_dimension *((double) maxtrials);
00320   MFLOPs = total_flops/elapsed_time/1000000.0;
00321   results[ntimings][0] = total_flops;
00322   results[ntimings][1] = elapsed_time;
00323   results[ntimings++][2] = MFLOPs;
00324 
00325   if (makeheader) 
00326     cout << "Metric_\t\t_Procs_\t__SpMV__\t_SpMM2__\t_SpMM4__\t_SpMM8__\t__NORM___\t__DOT___\t__AXPY__" << endl;
00327   if (comm.MyPID()==0) {
00328     cout.setf(std::ios::scientific);
00329     cout.precision(2);
00330     for (int j=0; j<3; j++) {
00331       bool doloop = false;
00332       if (j==0 && printflops) { cout << "FLOPS\t"; doloop = true;}
00333       else if (j==1 && printtime) {cout << "Time(s)\t"; doloop = true;}
00334       else if (j==2 && printmflops) {cout << "MFLOP/s\t"; doloop = true;}
00335       if (doloop) {
00336   cout << "\t" << comm.NumProc();
00337   for (int i=0; i<maxtimings; i++)
00338     cout << "\t" << results[i][j];
00339   cout << endl;
00340       }
00341     }
00342   }
00343 
00344 #ifdef EPETRA_MPI
00345   MPI_Finalize() ;
00346 #endif
00347 
00348 return ierr ;
00349 }
00350 
00351 // Constructs a 2D PDE finite difference matrix using the list of x and y offsets.
00352 // 
00353 // nx      (In) - number of grid points in x direction
00354 // ny      (In) - number of grid points in y direction
00355 //   The total number of equations will be nx*ny ordered such that the x direction changes
00356 //   most rapidly: 
00357 //      First equation is at point (0,0)
00358 //      Second at                  (1,0)
00359 //       ...
00360 //      nx equation at             (nx-1,0)
00361 //      nx+1st equation at         (0,1)
00362 
00363 // numPoints (In) - number of points in finite difference stencil
00364 // xoff    (In) - stencil offsets in x direction (of length numPoints)
00365 // yoff    (In) - stencil offsets in y direction (of length numPoints)
00366 //   A standard 5-point finite difference stencil would be described as:
00367 //     numPoints = 5
00368 //     xoff = [-1, 1, 0,  0, 0]
00369 //     yoff = [ 0, 0, 0, -1, 1]
00370 
00371 // nrhs - Number of rhs to generate. (First interface produces vectors, so nrhs is not needed
00372 
00373 // comm    (In) - an Epetra_Comm object describing the parallel machine (numProcs and my proc ID)
00374 // map    (Out) - Epetra_Map describing distribution of matrix and vectors/multivectors
00375 // A      (Out) - Epetra_CrsMatrix constructed for nx by ny grid using prescribed stencil
00376 //                Off-diagonal values are random between 0 and 1.  If diagonal is part of stencil,
00377 //                diagonal will be slightly diag dominant.
00378 // b      (Out) - Generated RHS.  Values satisfy b = A*xexact
00379 // bt     (Out) - Generated RHS.  Values satisfy b = A'*xexact
00380 // xexact (Out) - Generated exact solution to Ax = b and b' = A'xexact
00381 
00382 // Note: Caller of this function is responsible for deleting all output objects.
00383 
00384 void GenerateCrsProblem(int numNodesX, int numNodesY, int numProcsX, int numProcsY, int numPoints, 
00385       int * xoff, int * yoff,
00386       const Epetra_Comm  &comm, 
00387       Epetra_Map *& map, 
00388       Epetra_CrsMatrix *& A, 
00389       Epetra_Vector *& b, 
00390       Epetra_Vector *& bt,
00391       Epetra_Vector *&xexact, bool StaticProfile, bool MakeLocalOnly) {
00392 
00393   Epetra_MultiVector * b1, * bt1, * xexact1;
00394   
00395   GenerateCrsProblem(numNodesX, numNodesY, numProcsX, numProcsY, numPoints, 
00396          xoff, yoff, 1, comm, 
00397          map, A, b1, bt1, xexact1, StaticProfile, MakeLocalOnly);
00398 
00399   b = dynamic_cast<Epetra_Vector *>(b1);
00400   bt = dynamic_cast<Epetra_Vector *>(bt1);
00401   xexact = dynamic_cast<Epetra_Vector *>(xexact1);
00402 
00403   return;
00404 }
00405 
00406 void GenerateCrsProblem(int numNodesX, int numNodesY, int numProcsX, int numProcsY, int numPoints, 
00407       int * xoff, int * yoff, int nrhs,
00408       const Epetra_Comm  &comm,
00409       Epetra_Map *& map, 
00410       Epetra_CrsMatrix *& A, 
00411       Epetra_MultiVector *& b, 
00412       Epetra_MultiVector *& bt,
00413       Epetra_MultiVector *&xexact, bool StaticProfile, bool MakeLocalOnly) {
00414   
00415   Epetra_Time timer(comm);
00416   // Determine my global IDs
00417   int * myGlobalElements;
00418   GenerateMyGlobalElements(numNodesX, numNodesY, numProcsX, numProcsY, comm.MyPID(), myGlobalElements);
00419 
00420   int numMyEquations = numNodesX*numNodesY;
00421   
00422   map = new Epetra_Map(-1, numMyEquations, myGlobalElements, 0, comm); // Create map with 2D block partitioning.
00423   delete [] myGlobalElements;
00424 
00425   int numGlobalEquations = map->NumGlobalElements();
00426 
00427   int profile = 0; if (StaticProfile) profile = numPoints;
00428 
00429 #ifdef EPETRA_HAVE_STATICPROFILE
00430 
00431   if (MakeLocalOnly) 
00432     A = new Epetra_CrsMatrix(Copy, *map, *map, profile, StaticProfile); // Construct matrix with rowmap=colmap
00433   else
00434     A = new Epetra_CrsMatrix(Copy, *map, profile, StaticProfile); // Construct matrix
00435 
00436 #else
00437 
00438   if (MakeLocalOnly) 
00439     A = new Epetra_CrsMatrix(Copy, *map, *map, profile); // Construct matrix with rowmap=colmap
00440   else
00441     A = new Epetra_CrsMatrix(Copy, *map, profile); // Construct matrix
00442 
00443 #endif
00444 
00445   int * indices = new int[numPoints];
00446   double * values = new double[numPoints];
00447 
00448   double dnumPoints = (double) numPoints;
00449   int nx = numNodesX*numProcsX;
00450 
00451   for (int i=0; i<numMyEquations; i++) {
00452 
00453     int rowID = map->GID(i);
00454     int numIndices = 0;
00455 
00456     for (int j=0; j<numPoints; j++) {
00457       int colID = rowID + xoff[j] + nx*yoff[j]; // Compute column ID based on stencil offsets
00458       if (colID>-1 && colID<numGlobalEquations) {
00459   indices[numIndices] = colID;
00460   double value = - ((double) rand())/ ((double) RAND_MAX);
00461   if (colID==rowID)
00462     values[numIndices++] = dnumPoints - value; // Make diagonal dominant
00463   else
00464     values[numIndices++] = value;
00465       }
00466     }
00467     //cout << "Building row " << rowID << endl;
00468     A->InsertGlobalValues(rowID, numIndices, values, indices);
00469   }
00470 
00471   delete [] indices;
00472   delete [] values;
00473   double insertTime = timer.ElapsedTime();
00474   timer.ResetStartTime();
00475   A->FillComplete();
00476   double fillCompleteTime = timer.ElapsedTime();
00477 
00478   if (nrhs<=1) {  
00479     b = new Epetra_Vector(*map);
00480     bt = new Epetra_Vector(*map);
00481     xexact = new Epetra_Vector(*map);
00482   }
00483   else {
00484     b = new Epetra_MultiVector(*map, nrhs);
00485     bt = new Epetra_MultiVector(*map, nrhs);
00486     xexact = new Epetra_MultiVector(*map, nrhs);
00487   }
00488 
00489   xexact->Random(); // Fill xexact with random values
00490 
00491   A->Multiply(false, *xexact, *b);
00492   A->Multiply(true, *xexact, *bt);
00493 
00494   return;
00495 }
00496 
00497 
00498 
00499 
00500 void GenerateMyGlobalElements(int numNodesX, int numNodesY, int numProcsX, int numProcs,
00501             int myPID, int * & myGlobalElements) {
00502 
00503   myGlobalElements = new int[numNodesX*numNodesY];
00504   int myProcX = myPID%numProcsX;
00505   int myProcY = myPID/numProcsX;
00506   int curGID = myProcY*(numProcsX*numNodesX)*numNodesY+myProcX*numNodesX;
00507   for (int j=0; j<numNodesY; j++) {
00508     for (int i=0; i<numNodesX; i++) {
00509       myGlobalElements[j*numNodesX+i] = curGID+i;
00510     }
00511     curGID+=numNodesX*numProcsX;
00512   }
00513   //for (int i=0; i<numNodesX*numNodesY; i++) cout << "MYPID " << myPID <<" GID "<< myGlobalElements[i] << endl;
00514   
00515   return;
00516 }
00517 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Wed Apr 13 09:58:40 2011 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.6.3