EpetraExt Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
MatrixMarket_IO.cpp
Go to the documentation of this file.
1/*
2//@HEADER
3// ***********************************************************************
4//
5// EpetraExt: Epetra Extended - Linear Algebra Services Package
6// Copyright (2011) Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39//
40// ***********************************************************************
41//@HEADER
42*/
43
44#include "Epetra_ConfigDefs.h"
45#include "EpetraExt_Version.h"
46
47#ifdef EPETRA_MPI
48#include "mpi.h"
49#include "Epetra_MpiComm.h"
50#else
51#include "Epetra_SerialComm.h"
52#endif
53
54#include "Trilinos_Util.h"
55#include "Epetra_Comm.h"
56#include "Epetra_Map.h"
57#include "Epetra_Time.h"
58#include "Epetra_BlockMap.h"
59#include "Epetra_MultiVector.h"
60#include "Epetra_Vector.h"
61#include "Epetra_Export.h"
62
63#include "Epetra_VbrMatrix.h"
64#include "Epetra_CrsMatrix.h"
65
67
68int main(int argc, char *argv[]) {
69
70#ifdef EPETRA_MPI
71 MPI_Init(&argc,&argv);
72 Epetra_MpiComm Comm (MPI_COMM_WORLD);
73#else
75#endif
76
77 int MyPID = Comm.MyPID();
78
79 bool verbose = true;
80 if (MyPID==0) verbose = true;
81
82 if (verbose)
83 cout << EpetraExt::EpetraExt_Version() << endl << endl;
84
85 cout << Comm << endl;
86
87 if(argc < 2 && verbose) {
88 cerr << "Usage: " << argv[0]
89 << " HB_filename" << endl;
90 return(1);
91
92 }
93
94 // Uncomment the next three lines to debug in mpi mode
95 //int tmp;
96 //if (MyPID==0) cin >> tmp;
97 //Comm.Barrier();
98
99 Epetra_Map * readMap;
100 Epetra_CrsMatrix * readA;
101 Epetra_Vector * readx;
102 Epetra_Vector * readb;
103 Epetra_Vector * readxexact;
104
105 // Call routine to read in HB problem
106 Trilinos_Util_ReadHb2Epetra(argv[1], Comm, readMap, readA, readx, readb, readxexact);
107
108 // Create uniform distributed map
109 Epetra_Map map(readMap->NumGlobalElements(), 0, Comm);
110
111 // Create Exporter to distribute read-in matrix and vectors
112
113 Epetra_Export exporter(*readMap, map);
114 Epetra_CrsMatrix A(Copy, map, 0);
115 Epetra_Vector x(map);
116 Epetra_Vector b(map);
117 Epetra_Vector xexact(map);
118
119 Epetra_Time FillTimer(Comm);
120 x.Export(*readx, exporter, Add);
121 b.Export(*readb, exporter, Add);
122 xexact.Export(*readxexact, exporter, Add);
123 Comm.Barrier();
124 double vectorRedistributeTime = FillTimer.ElapsedTime();
125 A.Export(*readA, exporter, Add);
126 Comm.Barrier();
127 double matrixRedistributeTime = FillTimer.ElapsedTime() - vectorRedistributeTime;
128 assert(A.FillComplete()==0);
129 Comm.Barrier();
130 double fillCompleteTime = FillTimer.ElapsedTime() - matrixRedistributeTime;
131 if (Comm.MyPID()==0) {
132 cout << "\n\n****************************************************" << endl;
133 cout << "\n Vector redistribute time (sec) = " << vectorRedistributeTime<< endl;
134 cout << " Matrix redistribute time (sec) = " << matrixRedistributeTime << endl;
135 cout << " Transform to Local time (sec) = " << fillCompleteTime << endl<< endl;
136 }
137 Epetra_Vector tmp1(*readMap);
138 Epetra_Vector tmp2(map);
139 readA->Multiply(false, *readxexact, tmp1);
140
141 A.Multiply(false, xexact, tmp2);
142 double residual;
143 tmp1.Norm2(&residual);
144 if (verbose) cout << "Norm of Ax from file = " << residual << endl;
145 tmp2.Norm2(&residual);
146 if (verbose) cout << "Norm of Ax after redistribution = " << residual << endl << endl << endl;
147
148 //cout << "A from file = " << *readA << endl << endl << endl;
149
150 //cout << "A after dist = " << A << endl << endl << endl;
151
152 delete readA;
153 delete readx;
154 delete readb;
155 delete readxexact;
156 delete readMap;
157
158 Comm.Barrier();
159
160 EpetraExt::RowMatrixToMatrixMarketFile("test.mm", A, "test matrix", "This is a test matrix");
161
162#ifdef EPETRA_MPI
163 MPI_Finalize() ;
164#endif
165
166return 0 ;
167}
Add
Copy
int main(int argc, char *argv[])
int NumGlobalElements() const
int FillComplete(bool OptimizeDataStorage=true)
int Multiply(bool TransA, const Epetra_Vector &x, Epetra_Vector &y) const
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
void Barrier() const
int MyPID() const
int Norm2(double *Result) const
double ElapsedTime(void) const
int RowMatrixToMatrixMarketFile(const char *filename, const Epetra_RowMatrix &A, const char *matrixName, const char *matrixDescription, bool writeHeader)
Writes an Epetra_RowMatrix object to a Matrix Market format file.
std::string EpetraExt_Version()