EpetraExt Package Browser (Single Doxygen Collection) Development
example/zoltan/cxx_main.cpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ***********************************************************************
00003 // 
00004 //     EpetraExt: Epetra Extended - 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 
00029 // EpetraExt::CrsGraph_Zoltan Example routine
00030 #include <Epetra_ConfigDefs.h>
00031 
00032 #ifdef EPETRA_MPI
00033 #include "Epetra_MpiComm.h"
00034 #include <mpi.h>
00035 #endif
00036 #include "Epetra_SerialComm.h"
00037 #include "Epetra_Time.h"
00038 #include "Epetra_BlockMap.h"
00039 #include "Epetra_CrsGraph.h"
00040 #include "Epetra_CrsMatrix.h"
00041 #include "Epetra_Vector.h"
00042 #include "Epetra_LinearProblem.h"
00043 
00044 #include "Trilinos_Util.h"
00045 
00046 #include "EpetraExt_Zoltan_CrsGraph.h"
00047 #include "EpetraExt_SymmRCM_CrsGraph.h"
00048 //#include "EpetraExt_ZoltanOrder_CrsGraph.h"
00049 #include "EpetraExt_LPTrans_From_GraphTrans.h"
00050 #include "EpetraExt_Version.h"
00051 
00052 #define perror(str) { fprintf(stderr,"%s\n",str);  exit(-1); }
00053 #define perror1(str,ierr) { fprintf(stderr,"%s %d\n",str,ierr);  exit(ierr); }
00054 
00055 int main(int argc, char *argv[]) {
00056 
00057   int i, ierr=0, returnierr=0;
00058 
00059 #ifdef EPETRA_MPI
00060 
00061   // Initialize MPI
00062 
00063   MPI_Init(&argc,&argv);
00064   int size, rank; // Number of MPI processes, My process ID
00065 
00066   MPI_Comm_size(MPI_COMM_WORLD, &size);
00067   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00068 
00069 #else
00070 
00071   int size = 1; // Serial case (not using MPI)
00072   int rank = 0;
00073 
00074 #endif
00075 
00076   bool verbose = false;
00077 
00078   // Check if we should print results to standard out
00079   if (argc>2) if (argv[2][0]=='-' && argv[2][1]=='v') verbose = true;
00080   if (argc<2) perror("error: enter name of data file on cmd line");
00081 
00082 #ifdef EPETRA_MPI
00083   Epetra_MpiComm Comm(MPI_COMM_WORLD);
00084 #else
00085   Epetra_SerialComm Comm;
00086 #endif
00087   if (!verbose) Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00088 
00089   int MyPID = Comm.MyPID();
00090   int NumProc = Comm.NumProc();
00091 
00092   if (verbose) {
00093     cout << EpetraExt::EpetraExt_Version() << endl << endl;
00094     cout << Comm << endl << flush;
00095   }
00096 
00097   Comm.Barrier();
00098   bool verbose_all = verbose;
00099 
00100   if (verbose) verbose = (MyPID==0);
00101 
00102   //Read in Matrix File and distribute
00103   int NumGlobalEqs;
00104   int NumLocalEqs;
00105   int NumNZs;
00106   int *NumRowNZs;
00107   double *Values;
00108   double *Xprime;
00109   double *B;
00110   double *X;
00111   int *Bindx;
00112   int NUpdate;
00113   int *Update;
00114 
00115   Trilinos_Util_read_hb(argv[1], Comm.MyPID(), &NumGlobalEqs, &NumNZs, &Values, &Bindx, &Xprime, &B, &X );
00116   Trilinos_Util_distrib_msr_matrix( Comm, &NumGlobalEqs, &NumNZs, &NUpdate, &Update, &Values, &Bindx, &Xprime, &B, &X );
00117 
00118   NumLocalEqs = NUpdate;
00119   NumRowNZs = new int[NumLocalEqs];
00120   for( int i = 0; i < NumLocalEqs; ++i ) NumRowNZs[i] = Bindx[i+1]-Bindx[i]+1;
00121 
00122   Epetra_Map Map(NumGlobalEqs,NumLocalEqs,Update,0,Comm);
00123 
00124   if( verbose ) cout << "Building Epetra_CrsMatrix" << endl;
00125 
00126   Epetra_CrsMatrix A(Copy, Map, NumRowNZs );
00127 
00128   //Add individual rows
00129   double *RowVals;
00130   int *ColInds;
00131   int NumEntries;
00132   for( int i = 0; i < NumLocalEqs; ++i )
00133   {
00134     RowVals = Values + Bindx[i];
00135     ColInds = Bindx + Bindx[i];
00136     NumEntries = Bindx[i+1] - Bindx[i];
00137     ierr = A.InsertGlobalValues( Update[i], NumEntries, RowVals, ColInds );
00138     if( ierr ) { printf("Row %d:", Update[i] ); perror1("Error Putting Row: ",ierr); }
00139     ierr = A.InsertGlobalValues( Update[i], 1, Values+i, Update+i);
00140     if( ierr ) { perror1("Error Putting Diag: ",ierr); }
00141   }
00142 
00143   ierr = A.FillComplete();
00144   if( ierr ) perror1("Error in FillComplete",ierr);
00145 
00146   Epetra_Vector XX(Copy,Map,X);
00147   Epetra_Vector BB(Copy,Map,B);
00148 
00149   Epetra_LinearProblem Prob(&A,&XX,&BB);
00150 
00151   //Generate Zoltan Load Balanced Version of Linear Problem
00152   if( verbose ) cout << "Creating Zoltan Partitioning Transform!\n";
00153 
00154   EpetraExt::Zoltan_CrsGraph * ZoltanTrans = new EpetraExt::Zoltan_CrsGraph();
00155   EpetraExt::LinearProblem_GraphTrans * ZoltanLPTrans =
00156     new EpetraExt::LinearProblem_GraphTrans( 
00157          *(dynamic_cast<EpetraExt::StructuralSameTypeTransform<Epetra_CrsGraph>*>(ZoltanTrans)) );
00158 
00159   if( verbose ) cout << "Creating Load Balanced Linear Problem\n";
00160   Epetra_LinearProblem &BalancedProb = (*ZoltanLPTrans)(Prob);
00161 
00162   EpetraExt::CrsGraph_SymmRCM * RCMTrans = new EpetraExt::CrsGraph_SymmRCM();
00163   EpetraExt::LinearProblem_GraphTrans * RCMLPTrans =
00164     new EpetraExt::LinearProblem_GraphTrans( 
00165          *(dynamic_cast<EpetraExt::StructuralSameTypeTransform<Epetra_CrsGraph>*>(RCMTrans)) );
00166 
00167   if( verbose ) cout << "Creating SymmRCMed Linear Problem\n";
00168   Epetra_LinearProblem &RCMProb = (*RCMLPTrans)(BalancedProb);
00169 
00170 #ifdef EPETRA_MPI
00171   MPI_Finalize();
00172 #endif
00173 
00174   return returnierr;
00175 }
00176 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines