Ifpack Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/UseTranspose/cxx_main.cpp
Go to the documentation of this file.
1/*@HEADER
2// ***********************************************************************
3//
4// Ifpack: Object-Oriented Algebraic Preconditioner Package
5// Copyright (2002) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ***********************************************************************
40//@HEADER
41*/
42
43#include "Ifpack_ConfigDefs.h"
44#ifdef HAVE_MPI
45#include "mpi.h"
46#include "Epetra_MpiComm.h"
47#else
48#include "Epetra_SerialComm.h"
49#endif
50#include "Epetra_Comm.h"
51#include "Epetra_Map.h"
52#include "Epetra_MultiVector.h"
53#include "Epetra_CrsMatrix.h"
54#include "Ifpack_ILU.h"
55#include "Ifpack_ILUT.h"
57
58template<class T>
59void Test(const std::string what, Epetra_RowMatrix& A)
60{
61 using std::cout;
62 using std::endl;
63
64 T Prec(&A);
65
66 bool UseTranspose = true;
67
68 IFPACK_CHK_ERRV(Prec.Initialize());
69 IFPACK_CHK_ERRV(Prec.Compute());
70 IFPACK_CHK_ERRV(Prec.SetUseTranspose(UseTranspose));
71
72 Epetra_MultiVector LHS_exact(A.OperatorDomainMap(), 2);
75
76 LHS_exact.Random(); LHS.PutScalar(0.0);
77
78 A.Multiply(UseTranspose, LHS_exact, RHS);
79
80 Prec.ApplyInverse(RHS, LHS);
81
82 LHS.Update(1.0, LHS_exact, -1.0);
83 double norm[2];
84
85 LHS.Norm2(norm);
86 norm[0] += norm[1];
87
88 if (norm[0] > 1e-5)
89 {
90 cout << what << ": Test failed: norm = " << norm[0] << endl;
91 exit(EXIT_FAILURE);
92 }
93
94 cout << what << ": Test passed: norm = " << norm[0] << endl;
95}
96
97// =========== //
98// main driver //
99// =========== //
100
101int main(int argc, char *argv[])
102{
103#ifdef HAVE_MPI
104 MPI_Init(&argc,&argv);
105 Epetra_MpiComm Comm (MPI_COMM_WORLD);
106#else
108#endif
109
110 if (Comm.NumProc() != 1)
111 {
112 cerr << "To be run with one processor only" << endl;
113#ifdef HAVE_MPI
114 MPI_Finalize();
115#endif
116 exit(EXIT_SUCCESS);
117 }
118
119 Epetra_Map Map(8, 0, Comm);
120
121 Epetra_CrsMatrix A(Copy, Map, 0);
122
123 // for this matrix the incomplete factorization
124 // is the exact one, so ILU and ILUT must be exact solvers.
125 for (int row = 0; row < 8; ++row)
126 {
127 double value = 2.0 + row;
128 A.InsertGlobalValues(row, 1, &value, &row);
129 if (row)
130 {
131 int col = row - 1;
132 value = 1.0 + row;
133 A.InsertGlobalValues(row, 1, &value, &col);
134 }
135#if 0
136 if (row != Map.NumGlobalElements() - 1)
137 {
138 int col = row + 1;
139 value = 0.0;
140 A.InsertGlobalValues(row, 1, &value, &col);
141 }
142#endif
143 }
144
145 A.FillComplete();
146
147 Test<Ifpack_ILU>("Ifpack_ILU", A);
148 Test<Ifpack_ILUT>("Ifpack_ILUT", A);
149 Test<Ifpack_AdditiveSchwarz<Ifpack_ILU> >("AS, Ifpack_ILU", A);
150 Test<Ifpack_AdditiveSchwarz<Ifpack_ILUT> >("AS, Ifpack_ILUT", A);
151
152#ifdef HAVE_MPI
153 MPI_Finalize() ;
154#endif
155
156 return(EXIT_SUCCESS);
157}
Copy
#define IFPACK_CHK_ERRV(ifpack_err)
#define RHS(a)
Definition MatGenFD.c:60
int NumGlobalElements() const
int FillComplete(bool OptimizeDataStorage=true)
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
int NumProc() const
int Update(double ScalarA, const Epetra_MultiVector &A, double ScalarThis)
int Norm2(double *Result) const
int PutScalar(double ScalarConstant)
virtual const Epetra_Map & OperatorDomainMap() const=0
virtual const Epetra_Map & OperatorRangeMap() const=0
virtual int Multiply(bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const=0
int main(int argc, char *argv[])
void Test(const std::string what, Epetra_RowMatrix &A)