Amesos Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Amesos_TestMultiSolver.cpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Amesos: Direct Sparse Solver Package
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#include "Teuchos_ParameterList.hpp"
31//#include "Trilinos_Util_ReadTriples2Epetra.h"
32//#include "Trilinos_Util_ReadMatrixMarket2Epetra.h"
33#include "Trilinos_Util.h"
34#include "Epetra_LocalMap.h"
35#include "Epetra_Map.h"
36#include "Epetra_Vector.h"
37#include "Epetra_Export.h"
38#include "Epetra_CrsMatrix.h"
39#include "Epetra_LinearProblem.h"
40#include "Epetra_Time.h"
41#ifdef HAVE_AMESOS_SLUD
42#include "SuperludistOO.h"
43#endif
44#ifdef HAVE_AMESOS_SLUS
45#include "Epetra_SLU.h"
46#endif
47#ifdef HAVE_AMESOS_SLUD2
48#include "Superludist2_OO.h"
49#endif
50#ifdef HAVE_AMESOS_DSCPACK
51#include "Amesos_Dscpack.h"
52#endif
53#ifdef HAVE_AMESOS_UMFPACK
54#include "Amesos_Umfpack.h"
55#endif
56#ifdef HAVE_AMESOS_KLU
57#include "Amesos_Klu.h"
58#endif
59#ifdef HAVE_AMESOS_LAPACK
60#include "Amesos_Lapack.h"
61#endif
62#ifdef HAVE_AMESOS_TAUCS
63#include "Amesos_Taucs.h"
64#endif
65#ifdef HAVE_AMESOS_PARDISO
66#include "Amesos_Pardiso.h"
67#endif
68#ifdef HAVE_AMESOS_PARAKLETE
69#include "Amesos_Paraklete.h"
70#endif
71#if defined(HAVE_AMESOS_MUMPS) && defined(HAVE_MPI)
72#include "Amesos_Mumps.h"
73#endif
74#ifdef HAVE_AMESOS_SCALAPACK
75#include "Amesos_Scalapack.h"
76#endif
77#ifdef HAVE_AMESOS_SUPERLUDIST
78#include "Amesos_Superludist.h"
79#endif
80#ifdef HAVE_AMESOS_SUPERLU
81#include "Amesos_Superlu.h"
82#endif
83#ifdef TEST_SPOOLES
84#include "SpoolesOO.h"
85#endif
86
87#include "Amesos_TestSolver.h"
88#include "CrsMatrixTranspose.h"
90
91//
92// Amesos_TestMultiSolver.cpp reads in a matrix in Harwell-Boeing format,
93// calls one of the sparse direct solvers, using blocked right hand sides
94// and computes the error and residual.
95//
96// TestSolver ignores the Harwell-Boeing right hand sides, creating
97// random right hand sides instead.
98//
99// Amesos_TestMultiSolver can test either A x = b or A^T x = b.
100// This can be a bit confusing because sparse direct solvers
101// use compressed column storage - the transpose of Trilinos'
102// sparse row storage.
103//
104// Matrices:
105// readA - Serial. As read from the file.
106// transposeA - Serial. The transpose of readA.
107// serialA - if (transpose) then transposeA else readA
108// distributedA - readA distributed to all processes
109// passA - if ( distributed ) then distributedA else serialA
110//
111//
112int Amesos_TestMultiSolver( Epetra_Comm &Comm, char *matrix_file, int numsolves,
113 SparseSolverType SparseSolver, bool transpose,
114 int special, AMESOS_MatrixType matrix_type ) {
115
116
117 int iam = Comm.MyPID() ;
118
119
120 // int hatever;
121 // if ( iam == 0 ) std::cin >> hatever ;
122 Comm.Barrier();
123
124
125 Epetra_Map * readMap;
126 Epetra_CrsMatrix * readA;
127 Epetra_Vector * readx;
128 Epetra_Vector * readb;
129 Epetra_Vector * readxexact;
130
131 std::string FileName = matrix_file ;
132 int FN_Size = FileName.size() ;
133 std::string LastFiveBytes = FileName.substr( EPETRA_MAX(0,FN_Size-5), FN_Size );
134 std::string LastFourBytes = FileName.substr( EPETRA_MAX(0,FN_Size-4), FN_Size );
135 bool NonContiguousMap = false;
136
137 if ( LastFiveBytes == ".triU" ) {
138 NonContiguousMap = true;
139 // Call routine to read in unsymmetric Triplet matrix
140 EPETRA_CHK_ERR( Trilinos_Util_ReadTriples2Epetra( matrix_file, false, Comm, readMap, readA, readx,
141 readb, readxexact, NonContiguousMap ) );
142 } else {
143 if ( LastFiveBytes == ".triS" ) {
144 NonContiguousMap = true;
145 // Call routine to read in symmetric Triplet matrix
146 EPETRA_CHK_ERR( Trilinos_Util_ReadTriples2Epetra( matrix_file, true, Comm,
147 readMap, readA, readx,
148 readb, readxexact, NonContiguousMap ) );
149 } else {
150 if ( LastFourBytes == ".mtx" ) {
151 EPETRA_CHK_ERR( Trilinos_Util_ReadMatrixMarket2Epetra( matrix_file, Comm, readMap,
152 readA, readx, readb, readxexact) );
153 } else {
154 // Call routine to read in HB problem
155 Trilinos_Util_ReadHb2Epetra( matrix_file, Comm, readMap, readA, readx,
156 readb, readxexact) ;
157 }
158 }
159 }
160
161 Epetra_CrsMatrix transposeA(Copy, *readMap, 0);
162 Epetra_CrsMatrix *serialA ;
163
164 if ( transpose ) {
165 assert( CrsMatrixTranspose( readA, &transposeA ) == 0 );
166 serialA = &transposeA ;
167 } else {
168 serialA = readA ;
169 }
170
171 // Create uniform distributed map
172 Epetra_Map map(readMap->NumGlobalElements(), 0, Comm);
173 Epetra_Map* map_;
174
175 if( NonContiguousMap ) {
176 //
177 // map gives us NumMyElements and MyFirstElement;
178 //
179 int NumGlobalElements = readMap->NumGlobalElements();
180 int NumMyElements = map.NumMyElements();
181 int MyFirstElement = map.MinMyGID();
182 std::vector<int> MapMap_( NumGlobalElements );
183 readMap->MyGlobalElements( &MapMap_[0] ) ;
184 Comm.Broadcast( &MapMap_[0], NumGlobalElements, 0 ) ;
185 map_ = new Epetra_Map( NumGlobalElements, NumMyElements, &MapMap_[MyFirstElement], 0, Comm);
186 } else {
187 map_ = new Epetra_Map( map ) ;
188 }
189
190
191 // Create Exporter to distribute read-in matrix and vectors
192 Epetra_Export exporter(*readMap, *map_);
193 Epetra_CrsMatrix A(Copy, *map_, 0);
194
195 Epetra_RowMatrix * passA = 0;
196 Epetra_MultiVector * passx = 0;
197 Epetra_MultiVector * passb = 0;
198 Epetra_MultiVector * passxexact = 0;
199 Epetra_MultiVector * passresid = 0;
200 Epetra_MultiVector * passtmp = 0;
201
202 Epetra_MultiVector x(*map_,numsolves);
203 Epetra_MultiVector b(*map_,numsolves);
204 Epetra_MultiVector xexact(*map_,numsolves);
205 Epetra_MultiVector resid(*map_,numsolves);
206 Epetra_MultiVector tmp(*map_,numsolves);
207
208 Epetra_MultiVector serialx(*readMap,numsolves);
209 Epetra_MultiVector serialb(*readMap,numsolves);
210 Epetra_MultiVector serialxexact(*readMap,numsolves);
211 Epetra_MultiVector serialresid(*readMap,numsolves);
212 Epetra_MultiVector serialtmp(*readMap,numsolves);
213
214 bool distribute_matrix = ( matrix_type == AMESOS_Distributed ) ;
215 if ( distribute_matrix ) {
216 //
217 // Initialize x, b and xexact to the values read in from the file
218 //
219
220 A.Export(*serialA, exporter, Add);
221 Comm.Barrier();
222
223 assert(A.FillComplete()==0);
224 Comm.Barrier();
225
226 passA = &A;
227 passx = &x;
228 passb = &b;
229 passxexact = &xexact;
230 passresid = &resid;
231 passtmp = &tmp;
232 } else {
233 passA = serialA;
234 passx = &serialx;
235 passb = &serialb;
236 passxexact = &serialxexact;
237 passresid = &serialresid;
238 passtmp = &serialtmp;
239 }
240
241 passxexact->SetSeed(131) ;
242 passxexact->Random();
243 passx->SetSeed(11231) ;
244 passx->Random();
245
246 passb->PutScalar( 0.0 );
247 passA->Multiply( transpose, *passxexact, *passb ) ;
248
249 Epetra_MultiVector CopyB( *passb ) ;
250
251 double Anorm = passA->NormInf() ;
252 SparseDirectTimingVars::SS_Result.Set_Anorm(Anorm) ;
253
254 Epetra_LinearProblem Problem( (Epetra_RowMatrix *) passA,
255 (Epetra_MultiVector *) passx,
256 (Epetra_MultiVector *) passb );
257
258 double max_resid = 0.0;
259 for ( int j = 0 ; j < special+1 ; j++ ) {
260
261 Epetra_Time TotalTime( Comm ) ;
262 if ( false ) {
263#ifdef TEST_UMFPACK
264
265 unused code
266
267 } else if ( SparseSolver == UMFPACK ) {
268 UmfpackOO umfpack( (Epetra_RowMatrix *) passA,
269 (Epetra_MultiVector *) passx,
270 (Epetra_MultiVector *) passb ) ;
271
272 umfpack.SetTrans( transpose ) ;
273 umfpack.Solve() ;
274#endif
275#ifdef TEST_SUPERLU
276 } else if ( SparseSolver == SuperLU ) {
277 SuperluserialOO superluserial( (Epetra_RowMatrix *) passA,
278 (Epetra_MultiVector *) passx,
279 (Epetra_MultiVector *) passb ) ;
280
281 superluserial.SetPermc( SuperLU_permc ) ;
282 superluserial.SetTrans( transpose ) ;
283 superluserial.SetUseDGSSV( special == 0 ) ;
284 superluserial.Solve() ;
285#endif
286#ifdef HAVE_AMESOS_SLUD
287 } else if ( SparseSolver == SuperLUdist ) {
288 SuperludistOO superludist( Problem ) ;
289 superludist.SetTrans( transpose ) ;
290 EPETRA_CHK_ERR( superludist.Solve( true ) ) ;
291#endif
292#ifdef HAVE_AMESOS_SLUD2
293 } else if ( SparseSolver == SuperLUdist2 ) {
294 Superludist2_OO superludist2( Problem ) ;
295 superludist2.SetTrans( transpose ) ;
296 EPETRA_CHK_ERR( superludist2.Solve( true ) ) ;
297#endif
298#ifdef TEST_SPOOLES
299 } else if ( SparseSolver == SPOOLES ) {
300 SpoolesOO spooles( (Epetra_RowMatrix *) passA,
301 (Epetra_MultiVector *) passx,
302 (Epetra_MultiVector *) passb ) ;
303
304 spooles.SetTrans( transpose ) ;
305 spooles.Solve() ;
306#endif
307#ifdef HAVE_AMESOS_DSCPACK
308 } else if ( SparseSolver == DSCPACK ) {
309 Teuchos::ParameterList ParamList ;
310 Amesos_Dscpack dscpack( Problem ) ;
311 ParamList.set( "MaxProcs", -3 );
312 EPETRA_CHK_ERR( dscpack.SetParameters( ParamList ) );
313
314 EPETRA_CHK_ERR( dscpack.Solve( ) );
315#endif
316#ifdef HAVE_AMESOS_UMFPACK
317 } else if ( SparseSolver == UMFPACK ) {
318 Teuchos::ParameterList ParamList ;
319 Amesos_Umfpack umfpack( Problem ) ;
320 ParamList.set( "MaxProcs", -3 );
321 EPETRA_CHK_ERR( umfpack.SetParameters( ParamList ) );
322 EPETRA_CHK_ERR( umfpack.SetUseTranspose( transpose ) );
323
324 EPETRA_CHK_ERR( umfpack.Solve( ) );
325#endif
326#ifdef HAVE_AMESOS_KLU
327 } else if ( SparseSolver == KLU ) {
328 Teuchos::ParameterList ParamList ;
329 Amesos_Klu klu( Problem ) ;
330 ParamList.set( "MaxProcs", -3 );
331 EPETRA_CHK_ERR( klu.SetParameters( ParamList ) );
332 EPETRA_CHK_ERR( klu.SetUseTranspose( transpose ) );
333
334 EPETRA_CHK_ERR( klu.SymbolicFactorization( ) );
335 EPETRA_CHK_ERR( klu.NumericFactorization( ) );
336 EPETRA_CHK_ERR( klu.Solve( ) );
337#endif
338#ifdef HAVE_AMESOS_PARAKLETE
339 } else if ( SparseSolver == PARAKLETE ) {
340 Teuchos::ParameterList ParamList ;
341 Amesos_Paraklete paraklete( Problem ) ;
342 ParamList.set( "MaxProcs", -3 );
343 EPETRA_CHK_ERR( paraklete.SetParameters( ParamList ) );
344 EPETRA_CHK_ERR( paraklete.SetUseTranspose( transpose ) );
345
346 EPETRA_CHK_ERR( paraklete.SymbolicFactorization( ) );
347 EPETRA_CHK_ERR( paraklete.NumericFactorization( ) );
348 EPETRA_CHK_ERR( paraklete.Solve( ) );
349#endif
350#ifdef HAVE_AMESOS_SLUS
351 } else if ( SparseSolver == SuperLU ) {
352 Epetra_SLU superluserial( &Problem ) ;
353 EPETRA_CHK_ERR( superluserial.SetUseTranspose( transpose ) );
354
355 EPETRA_CHK_ERR( superluserial.SymbolicFactorization( ) );
356 EPETRA_CHK_ERR( superluserial.NumericFactorization( ) );
357
358 EPETRA_CHK_ERR( superluserial.Solve( ) );
359#endif
360#ifdef HAVE_AMESOS_LAPACK
361 } else if ( SparseSolver == LAPACK ) {
362 Teuchos::ParameterList ParamList ;
363 ParamList.set( "MaxProcs", -3 );
364 Amesos_Lapack lapack( Problem ) ;
365 EPETRA_CHK_ERR( lapack.SetUseTranspose( transpose ) );
366
367 EPETRA_CHK_ERR( lapack.SymbolicFactorization( ) );
368 EPETRA_CHK_ERR( lapack.NumericFactorization( ) );
369 EPETRA_CHK_ERR( lapack.Solve( ) );
370#endif
371#ifdef HAVE_AMESOS_TAUCS
372 } else if ( SparseSolver == TAUCS ) {
373 Teuchos::ParameterList ParamList ;
374 Amesos_Taucs taucs( Problem ) ;
375 ParamList.set( "MaxProcs", -3 );
376 EPETRA_CHK_ERR( taucs.SetParameters( ParamList ) );
377 EPETRA_CHK_ERR( taucs.SetUseTranspose( transpose ) );
378
379 EPETRA_CHK_ERR( taucs.SymbolicFactorization( ) );
380 EPETRA_CHK_ERR( taucs.NumericFactorization( ) );
381 EPETRA_CHK_ERR( taucs.Solve( ) );
382#endif
383#ifdef HAVE_AMESOS_PARDISO
384 } else if ( SparseSolver == PARDISO ) {
385 Teuchos::ParameterList ParamList ;
386 Amesos_Pardiso pardiso( Problem ) ;
387 ParamList.set( "MaxProcs", -3 );
388 EPETRA_CHK_ERR( pardiso.SetParameters( ParamList ) );
389 EPETRA_CHK_ERR( pardiso.SetUseTranspose( transpose ) );
390
391 EPETRA_CHK_ERR( pardiso.SymbolicFactorization( ) );
392 EPETRA_CHK_ERR( pardiso.NumericFactorization( ) );
393 EPETRA_CHK_ERR( pardiso.Solve( ) );
394#endif
395#ifdef HAVE_AMESOS_PARKLETE
396 } else if ( SparseSolver == PARKLETE ) {
397 Teuchos::ParameterList ParamList ;
398 Amesos_Parklete parklete( Problem ) ;
399 ParamList.set( "MaxProcs", -3 );
400 EPETRA_CHK_ERR( parklete.SetParameters( ParamList ) );
401 EPETRA_CHK_ERR( parklete.SetUseTranspose( transpose ) );
402
403 EPETRA_CHK_ERR( parklete.SymbolicFactorization( ) );
404 EPETRA_CHK_ERR( parklete.NumericFactorization( ) );
405 EPETRA_CHK_ERR( parklete.Solve( ) );
406#endif
407#if defined(HAVE_AMESOS_MUMPS) && defined(HAVE_MPI)
408 } else if ( SparseSolver == MUMPS ) {
409 Teuchos::ParameterList ParamList ;
410 Amesos_Mumps mumps( Problem ) ;
411 ParamList.set( "MaxProcs", -3 );
412 EPETRA_CHK_ERR( mumps.SetParameters( ParamList ) );
413 EPETRA_CHK_ERR( mumps.SetUseTranspose( transpose ) );
414
415 EPETRA_CHK_ERR( mumps.SymbolicFactorization( ) );
416 EPETRA_CHK_ERR( mumps.NumericFactorization( ) );
417 EPETRA_CHK_ERR( mumps.Solve( ) );
418#endif
419#ifdef HAVE_AMESOS_SCALAPACK
420 } else if ( SparseSolver == SCALAPACK ) {
421 Teuchos::ParameterList ParamList ;
422 Amesos_Scalapack scalapack( Problem ) ;
423 ParamList.set( "MaxProcs", -3 );
424 EPETRA_CHK_ERR( scalapack.SetParameters( ParamList ) );
425 EPETRA_CHK_ERR( scalapack.SetUseTranspose( transpose ) );
426
427 EPETRA_CHK_ERR( scalapack.SymbolicFactorization( ) );
428 EPETRA_CHK_ERR( scalapack.NumericFactorization( ) );
429 EPETRA_CHK_ERR( scalapack.Solve( ) );
430#endif
431#ifdef HAVE_AMESOS_SUPERLUDIST
432 } else if ( SparseSolver == SUPERLUDIST ) {
433 Teuchos::ParameterList ParamList ;
434 Amesos_Superludist superludist( Problem ) ;
435 ParamList.set( "MaxProcs", -3 );
436 EPETRA_CHK_ERR( superludist.SetParameters( ParamList ) );
437
438 EPETRA_CHK_ERR( superludist.SetUseTranspose( transpose ) );
439
440 EPETRA_CHK_ERR( superludist.SymbolicFactorization( ) );
441 EPETRA_CHK_ERR( superludist.NumericFactorization( ) );
442 EPETRA_CHK_ERR( superludist.Solve( ) );
443#endif
444#ifdef HAVE_AMESOS_SUPERLU
445 } else if ( SparseSolver == SUPERLU ) {
446 Teuchos::ParameterList ParamList ;
447 Amesos_Superlu superlu( Problem ) ;
448 ParamList.set( "MaxProcs", -3 );
449 EPETRA_CHK_ERR( superlu.SetParameters( ParamList ) );
450 EPETRA_CHK_ERR( superlu.SetUseTranspose( transpose ) );
451
452 EPETRA_CHK_ERR( superlu.SymbolicFactorization( ) );
453 EPETRA_CHK_ERR( superlu.NumericFactorization( ) );
454 EPETRA_CHK_ERR( superlu.Solve( ) );
455#endif
456#ifdef TEST_SPOOLESSERIAL
457 } else if ( SparseSolver == SPOOLESSERIAL ) {
458 SpoolesserialOO spoolesserial( (Epetra_RowMatrix *) passA,
459 (Epetra_MultiVector *) passx,
460 (Epetra_MultiVector *) passb ) ;
461
462 spoolesserial.Solve() ;
463#endif
464 } else {
465 SparseDirectTimingVars::log_file << "Solver not implemented yet" << std::endl ;
466 std::cerr << "\n\n#################### Requested solver not available (Or not tested with blocked RHS) on this platform #####################\n" << std::endl ;
467 }
468
469 SparseDirectTimingVars::SS_Result.Set_Total_Time( TotalTime.ElapsedTime() );
470 // SparseDirectTimingVars::SS_Result.Set_First_Time( 0.0 );
471 // SparseDirectTimingVars::SS_Result.Set_Middle_Time( 0.0 );
472 // SparseDirectTimingVars::SS_Result.Set_Last_Time( 0.0 );
473
474 //
475 // Compute the error = norm(xcomp - xexact )
476 //
477 std::vector <double> error(numsolves) ;
478 double max_error = 0.0;
479
480 passresid->Update(1.0, *passx, -1.0, *passxexact, 0.0);
481
482 passresid->Norm2(&error[0]);
483 for ( int i = 0 ; i< numsolves; i++ )
484 if ( error[i] > max_error ) max_error = error[i] ;
485 SparseDirectTimingVars::SS_Result.Set_Error(max_error) ;
486
487 // passxexact->Norm2(&error[0] ) ;
488 // passx->Norm2(&error ) ;
489
490 //
491 // Compute the residual = norm(Ax - b)
492 //
493 std::vector <double> residual(numsolves) ;
494
495 passtmp->PutScalar(0.0);
496 passA->Multiply( transpose, *passx, *passtmp);
497 passresid->Update(1.0, *passtmp, -1.0, *passb, 0.0);
498 // passresid->Update(1.0, *passtmp, -1.0, CopyB, 0.0);
499 passresid->Norm2(&residual[0]);
500
501 for ( int i = 0 ; i< numsolves; i++ )
502 if ( residual[i] > max_resid ) max_resid = residual[i] ;
503
504
505 SparseDirectTimingVars::SS_Result.Set_Residual(max_resid) ;
506
507 std::vector <double> bnorm(numsolves);
508 passb->Norm2( &bnorm[0] ) ;
509 SparseDirectTimingVars::SS_Result.Set_Bnorm(bnorm[0]) ;
510
511 std::vector <double> xnorm(numsolves);
512 passx->Norm2( &xnorm[0] ) ;
513 SparseDirectTimingVars::SS_Result.Set_Xnorm(xnorm[0]) ;
514
515
516 if ( false && iam == 0 ) {
517
518 std::cout << " Amesos_TestMutliSolver.cpp " << std::endl ;
519 for ( int i = 0 ; i< numsolves && i < 10 ; i++ ) {
520 std::cout << "i=" << i
521 << " error = " << error[i]
522 << " xnorm = " << xnorm[i]
523 << " residual = " << residual[i]
524 << " bnorm = " << bnorm[i]
525 << std::endl ;
526
527 }
528
529 std::cout << std::endl << " max_resid = " << max_resid ;
530 std::cout << " max_error = " << max_error << std::endl ;
531 std::cout << " Get_residual() again = " << SparseDirectTimingVars::SS_Result.Get_Residual() << std::endl ;
532
533 }
534 }
535 delete readA;
536 delete readx;
537 delete readb;
538 delete readxexact;
539 delete readMap;
540 delete map_;
541
542 Comm.Barrier();
543
544return 0 ;
545}
int SuperLU_permc
int Amesos_TestMultiSolver(Epetra_Comm &Comm, char *matrix_file, int numsolves, SparseSolverType SparseSolver, bool transpose, int special, AMESOS_MatrixType matrix_type)
AMESOS_MatrixType
@ AMESOS_Distributed
SparseSolverType
@ TAUCS
@ PARAKLETE
@ SCALAPACK
@ LAPACK
@ UMFPACK
@ SUPERLUDIST
@ MUMPS
@ DSCPACK
@ SUPERLU
@ PARDISO
int CrsMatrixTranspose(Epetra_CrsMatrix *In, Epetra_CrsMatrix *Out)
#define EPETRA_CHK_ERR(xxx)
Amesos_Dscpack: An object-oriented wrapper for Dscpack.
Amesos_Klu: A serial, unblocked code ideal for getting started and for very sparse matrices,...
Definition Amesos_Klu.h:116
Amesos_Lapack: an interface to LAPACK.
Amesos_Mumps: An object-oriented wrapper for the double precision version of MUMPS.
Amesos_Paraklete: A serial, unblocked code ideal for getting started and for very sparse matrices,...
Amesos_Pardiso: Interface to the PARDISO package.
Amesos_Scalapack: A serial and parallel dense solver. For now, we implement only the unsymmetric ScaL...
Amesos_Superlu: Amesos interface to Xioye Li's SuperLU 3.0 serial code.
Amesos_Superludist: An object-oriented wrapper for Superludist.
Amesos_Taucs: An interface to the TAUCS package.
Class Amesos_Umfpack: An object-oriented wrapper for UMFPACK.
Epetra_SLU: An object-oriented wrapper for Xiaoye Li's serial sparse solver package: Superlu.
Definition Epetra_SLU.h:52
static std::ofstream log_file
static SparseSolverResult SS_Result
SpoolesOO: An object-oriented wrapper for Spooles.
Definition SpoolesOO.h:51
Superludist2_OO: An object-oriented wrapper for Superludist.
SuperludistOO: An object-oriented wrapper for Xiaoye Li's Superludist.