Amesos Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
example_AmesosFactory.cpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Amesos: An Interface to Direct Solvers
5// Copyright (2004) 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// This library is free software; you can redistribute it and/or modify
11// it under the terms of the GNU Lesser General Public License as
12// published by the Free Software Foundation; either version 2.1 of the
13// License, or (at your option) any later version.
14//
15// This library is distributed in the hope that it will be useful, but
16// WITHOUT ANY WARRANTY; without even the implied warranty of
17// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18// Lesser General Public License for more details.
19//
20// You should have received a copy of the GNU Lesser General Public
21// License along with this library; if not, write to the Free Software
22// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23// USA
24// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
25//
26// ***********************************************************************
27// @HEADER
28
29#include "Amesos_ConfigDefs.h"
30
31// This example needs Galeri to generate the linear system.
32// You must have configured Trilinos with --enable-galeri
33// in order to compile this example
34
35#ifdef HAVE_MPI
36#include "mpi.h"
37#include "Epetra_MpiComm.h"
38#else
39#include "Epetra_SerialComm.h"
40#endif
41#include "Amesos.h"
42#include "Epetra_RowMatrix.h"
43#include "Epetra_MultiVector.h"
44#include "Epetra_LinearProblem.h"
45// following header file and namespace declaration
46// are required by this example to generate the linear system,
47// not by Amesos itself.
48#include "Galeri_Maps.h"
49#include "Galeri_CrsMatrices.h"
50using namespace Teuchos;
51using namespace Galeri;
52
53// ==================== //
54// M A I N D R I V E R //
55// ==================== //
56//
57// This example will:
58// 1.- create a linear system, stored as an
59// Epetra_LinearProblem. The matrix corresponds
60// to a 5pt Laplacian (2D on Cartesian grid).
61// The user can change the global size of the problem
62// by modifying variables nx and ny.
63// 2.- The linear system matrix, solution and rhs
64// are distributed among the available processors,
65// using a linear distribution. This is for
66// simplicity only! Amesos can support any Epetra_Map.
67// 3.- Once the linear problem is created, we
68// create an Amesos Factory object.
69// 4.- Using the Factory, we create the required Amesos_BaseSolver
70// solver. Any supported (and compiled) Amesos
71// solver can be used. If the selected solver
72// is not available (that is, if Amesos has *not*
73// been configured with support for this solver),
74// the factory returns 0. Usually, Amesos_Klu
75// is always available.
76// 5.- At this point we can factorize the matrix,
77// and solve the linear system. Only three methods
78// should be used for any Amesos_BaseSolver object:
79// 1.- NumericFactorization();
80// 2.- SymbolicFactorization();
81// 3.- Solve();
82// 6.- We note that the header files of Amesos-supported
83// libraries are *not* required in this file. They are
84// actually needed to compile the Amesos library only.
85//
86// NOTE: this example can be run with any number of processors.
87//
88// Author: Marzio Sala, SNL 9214
89// Last modified: Oct-05.
90
91int main(int argc, char *argv[])
92{
93
94#ifdef HAVE_MPI
95 MPI_Init(&argc, &argv);
96 Epetra_MpiComm Comm(MPI_COMM_WORLD);
97#else
98 Epetra_SerialComm Comm;
99#endif
100
101 int nx = 100; // number of grid points in the x direction
102 int ny = 100 * Comm.NumProc(); // number of grid points in the y direction
103 int NumVectors = 1; // number of rhs's. Amesos
104 // supports single or
105 // multiple RHS.
106
107 // Initializes an Gallery object.
108 // NOTE: this example uses the Trilinos package Galeri
109 // to define in an easy way the linear system matrix.
110 // The user can easily change the matrix type; consult the
111 // Galeri documentation for mode details.
112 //
113 // Here the problem has size nx x ny, and the 2D Cartesian
114 // grid is divided into mx x my subdomains.
115 ParameterList GaleriList;
116 GaleriList.set("nx", nx);
117 GaleriList.set("ny", ny);
118 GaleriList.set("mx", 1);
119 GaleriList.set("my", Comm.NumProc());
120
121 Epetra_Map* Map = CreateMap("Cartesian2D", Comm, GaleriList);
122 Epetra_CrsMatrix* Matrix = CreateCrsMatrix("Laplace2D", Map, GaleriList);
123
124 // Creates vectors for right-hand side and solution, and the
125 // linear problem container.
126
127 Epetra_MultiVector LHS(*Map, NumVectors); LHS.PutScalar(0.0); // zero solution
128 Epetra_MultiVector RHS(*Map, NumVectors); RHS.Random(); // random rhs
129 Epetra_LinearProblem Problem(Matrix, &LHS, &RHS);
130
131 // ===================================================== //
132 // B E G I N N I N G O F T H E AM E S O S P A R T //
133 // ===================================================== //
134
135 // Initializes the Amesos solver. This is the base class for
136 // Amesos. It is a pure virtual class (hence objects of this
137 // class cannot be allocated, and can exist only as pointers
138 // or references).
139 //
140 Amesos_BaseSolver* Solver;
141
142 // Initializes the Factory. Factory is a function class (a
143 // class that contains methods only, no data). Factory
144 // will be used to create Amesos_BaseSolver derived objects.
145 //
146 Amesos Factory;
147
148 // Specifies the solver. String ``SolverType'' can assume one
149 // of the following values:
150 // - Lapack
151 // - Klu
152 // - Umfpack
153 // - Pardiso
154 // - Taucs
155 // - Superlu
156 // - Superludist
157 // - Mumps
158 // - Dscpack
159 //
160 std::string SolverType = "Klu";
161 Solver = Factory.Create(SolverType, Problem);
162
163 // Factory.Create() returns 0 if the requested solver
164 // is not available
165 //
166 if (Solver == 0) {
167 std::cerr << "Specified solver is not available" << std::endl;
168 // return ok not to break test harness even if
169 // the solver is not available
170#ifdef HAVE_MPI
171 MPI_Finalize();
172#endif
173 return(EXIT_SUCCESS);
174 }
175
176 // Parameters for all Amesos solvers are set through
177 // a call to SetParameters(List). List is a Teuchos
178 // parameter list (Amesos requires Teuchos to compile).
179 // In most cases, users can proceed without calling
180 // SetParameters(). Please refer to the Amesos guide
181 // for more details.
182 // NOTE: you can skip this call; then the solver will
183 // use default parameters.
184 //
185 // Parameters in the list are set using
186 // List.set("parameter-name", ParameterValue);
187 // In this example, we specify that we want more output.
188 //
189 Teuchos::ParameterList List;
190 List.set("PrintTiming", true);
191 List.set("PrintStatus", true);
192
193 Solver->SetParameters(List);
194
195 // Now we are ready to solve. Generally, users will
196 // call SymbolicFactorization(), then NumericFactorization(),
197 // and finally Solve(). Note that:
198 // - the numerical values of the linear system matrix
199 // are *not* required before NumericFactorization();
200 // - solution and rhs are *not* required before calling
201 // Solve().
202 if (Comm.MyPID() == 0)
203 std::cout << "Starting symbolic factorization..." << std::endl;
204 Solver->SymbolicFactorization();
205
206 // you can change the matrix values here
207 if (Comm.MyPID() == 0)
208 std::cout << "Starting numeric factorization..." << std::endl;
209 Solver->NumericFactorization();
210
211 // you can change LHS and RHS here
212 if (Comm.MyPID() == 0)
213 std::cout << "Starting solution phase..." << std::endl;
214 Solver->Solve();
215
216
217 // Nothing done with timing infomation below, so commenting out.
218
219 // you can get the timings here
220 //Teuchos::ParameterList TimingsList;
221 //Solver->GetTiming( TimingsList );
222
223 // you can find out how much time was spent in ...
224 //double sfact_time, nfact_time, solve_time;
225 //double mtx_conv_time, mtx_redist_time, vec_redist_time;
226
227 // 1) The symbolic factorization
228 // (parameter doesn't always exist)
229 //sfact_time = TimingsList.get( "Total symbolic factorization time", 0.0 );
230
231 // 2) The numeric factorization
232 // (always exists if NumericFactorization() is called)
233 //nfact_time = Teuchos::getParameter<double>( TimingsList, "Total numeric factorization time" );
234
235 // 3) Solving the linear system
236 // (always exists if Solve() is called)
237 //solve_time = Teuchos::getParameter<double>( TimingsList, "Total solve time" );
238
239 // 4) Converting the matrix to the accepted format for the solver
240 // (always exists if SymbolicFactorization() is called)
241 //mtx_conv_time = Teuchos::getParameter<double>( TimingsList, "Total solve time" );
242
243 // 5) Redistributing the matrix for each solve to the accepted format for the solver
244 //mtx_redist_time = TimingsList.get( "Total matrix redistribution time", 0.0 );
245
246 // 6) Redistributing the vector for each solve to the accepted format for the solver
247 //vec_redist_time = TimingsList.get( "Total vector redistribution time", 0.0 );
248
249 // =========================================== //
250 // E N D O F T H E A M E S O S P A R T //
251 // =========================================== //
252
253 // delete Solver. MPI calls can occur.
254 delete Solver;
255
256 // delete the objects created by Galeri
257 delete Matrix;
258 delete Map;
259
260#ifdef HAVE_MPI
261 MPI_Finalize();
262#endif
263
264 return(EXIT_SUCCESS);
265
266} // end of main()
int CreateCrsMatrix(const char *in_filename, const Epetra_Comm &Comm, Epetra_Map *&readMap, const bool transpose, const bool distribute, bool &symmetric, Epetra_CrsMatrix *&Matrix)
Definition: TestOptions.cpp:81
Amesos_BaseSolver: A pure virtual class for direct solution of real-valued double-precision operators...
virtual int Solve()=0
Solves A X = B (or AT x = B)
virtual int NumericFactorization()=0
Performs NumericFactorization on the matrix A.
virtual int SymbolicFactorization()=0
Performs SymbolicFactorization on the matrix A.
virtual int SetParameters(Teuchos::ParameterList &ParameterList)=0
Updates internal variables.
Factory for binding a third party direct solver to an Epetra_LinearProblem.
Definition: Amesos.h:44
Amesos_BaseSolver * Create(const char *ClassType, const Epetra_LinearProblem &LinearProblem)
Amesos Create method.
Definition: Amesos.cpp:69
int main(int argc, char *argv[])