1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
|
///
/// This file is part of Rheolef.
///
/// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
///
/// Rheolef is free software; you can redistribute it and/or modify
/// it under the terms of the GNU General Public License as published by
/// the Free Software Foundation; either version 2 of the License, or
/// (at your option) any later version.
///
/// Rheolef is distributed in the hope that it will be useful,
/// but WITHOUT ANY WARRANTY; without even the implied warranty of
/// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
/// GNU General Public License for more details.
///
/// You should have received a copy of the GNU General Public License
/// along with Rheolef; if not, write to the Free Software
/// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
///
/// =========================================================================
#include "rheolef/solver.h"
#include "rheolef/pretty_name.h"
#include "solver_gmres_cg.h"
#include "solver_eigen.h"
#include "solver_cholmod.h"
#include "solver_umfpack.h"
#include "solver_mumps.h"
// no more supported:
#include "solver_pastix.h" // portage pb on new version
namespace rheolef {
// ----------------------------------------------------------------------------
// utility
// ----------------------------------------------------------------------------
template <class T, class M>
std::string
solver_abstract_rep<T,M>::name() const
{
return typeid_name(typeid(*this).name(),false);
}
// ----------------------------------------------------------------------------
// choose the best solver
// ----------------------------------------------------------------------------
template<class T, class M>
solver_abstract_rep<T,M>*
solver_abstract_rep<T,M>::make_solver_ptr (const csr<T,M>& a, const solver_option& opt)
{
using namespace std;
// --------------------------------------------------------------------------
// 1d, 2d => direct ; 3d => iterative
// --------------------------------------------------------------------------
bool is_direct = opt.iterative == 0 ||
(opt.iterative == solver_option::decide && a.pattern_dimension() <= 2);
// --------------------------------------------------------------------------
// iterative: no choice here, only the cg/gmres pair (choice on preconditioner)
// --------------------------------------------------------------------------
if (!is_direct) {
using rep = solver_gmres_cg_rep<T,M>;
return new_macro (rep(a,opt));
}
string lib = solver_option::used_library (a,opt);
// --------------------------------------------------------------------------
// direct & sequential
// * when symmetry: eigen, mumps, cholmod(bof)
// * when no-symmetry: umfpack, mumps, eigen
// direct & distributed: only mumps
// --------------------------------------------------------------------------
#ifdef _RHEOLEF_HAVE_MUMPS
if (lib == "mumps") {
using rep = solver_mumps_rep<T,M>;
return new_macro (rep(a,opt));
}
#endif // _RHEOLEF_HAVE_MUMPS
#if defined(_RHEOLEF_HAVE_CHOLMOD)
if (lib == "suitesparse" && a.is_symmetric() && a.is_definite_positive()) {
using rep = solver_cholmod_rep<T,M>;
return new_macro (rep(a,opt));
}
#endif // _RHEOLEF_HAVE_CHOLMOD
#if defined(_RHEOLEF_HAVE_UMFPACK)
if (lib == "suitesparse") {
using rep = solver_umfpack_rep<T,M>;
return new_macro (rep(a,opt));
}
#endif // _RHEOLEF_HAVE_UMFPACK
#ifdef _RHEOLEF_HAVE_EIGEN
if (lib == "eigen") {
using rep = solver_eigen_rep<T,M>;
return new_macro (rep(a,opt));
}
#endif // _RHEOLEF_HAVE_EIGEN
error_macro ("no direct solver available"); // should not apppend
return 0;
}
// ----------------------------------------------------------------------------
// instanciation in library
// ----------------------------------------------------------------------------
#define _RHEOLEF_instanciation(T,M) \
template solver_abstract_rep<T,M>* solver_abstract_rep<T,M>::make_solver_ptr (const csr<T,M>&, const solver_option&); \
template std::string solver_abstract_rep<T,M>::name() const;
_RHEOLEF_instanciation(Float,sequential)
#ifdef _RHEOLEF_HAVE_MPI
_RHEOLEF_instanciation(Float,distributed)
#endif // _RHEOLEF_HAVE_MPI
#undef _RHEOLEF_instanciation
} // namespace rheolef
|
