1from __future__ import absolute_import, division, print_function, unicode_literals
2
3import sys
4import logging
5logger = logging.getLogger(__name__)
6
7def solvers(includes, storage, operator):
8    _, dfIncludes, dfTypeName, linearOperatorType, _, _ = storage
9    includes += dfIncludes + ["dune/fempy/parameter.hh"]
10    typeName = operator(dfTypeName,linearOperatorType)
11    return includes, typeName
12
13
14def femsolver(storage,solverType=None):
15    includes = ["dune/fem/solver/krylovinverseoperators.hh"]
16
17    operator = lambda df,_: "Dune::Fem::KrylovInverseOperator< " + df + " >"
18
19    includes, typeName = solvers(includes,storage,operator)
20    parameter = {"newton.linear.method":solverType} if solverType is not None else {}
21    return "fem",includes,typeName, parameter
22
23def istlsolver(storage,solverType=None):
24    includes = ["dune/fem/solver/istlinverseoperators.hh"]
25
26    operator = lambda df,_: "Dune::Fem::ISTLInverseOperator< " + df + " >"
27
28    includes, typeName = solvers(includes,storage,operator)
29    parameter = {"newton.linear.method":solverType} if solverType is not None else {}
30    return "istl",includes,typeName, parameter
31
32def suitesparsesolver(storage,solverType="umfpack"):
33    includes = ["dune/fem/solver/ldlsolver.hh", "dune/fem/solver/spqrsolver.hh", "dune/fem/solver/umfpacksolver.hh"]
34
35    if solverType == "ldl":
36        operator = lambda df,linop: "Dune::Fem::LDLInverseOperator< "    +df+", typename "+linop+"::MatrixType" + " >"
37    elif solverType == "spqr_symmetric":
38        operator = lambda df,linop: "Dune::Fem::SPQRInverseOperator< "   +df+", true, typename "+linop+"::MatrixType" + " >"
39    elif solverType == "spqr_nonsymmetric":
40        operator = lambda df,linop: "Dune::Fem::SPQRInverseOperator< "   +df+", false, typename "+linop+"::MatrixType" + " >"
41    elif solverType == "umfpack":
42        operator = lambda df,linop: "Dune::Fem::UMFPACKInverseOperator< "+df+", typename "+linop+"::MatrixType" + " >"
43    else:
44        raise ValueError("wrong krylov solver - only ldl,spqr_symmetric,spqr_nonsymmetric,umfpack available")
45
46    includes, typeName = solvers(includes,storage,operator)
47    return "suitesparse",includes,typeName,{}
48
49def eigensolver(storage,solverType="bicgstab"):
50    includes = ["dune/fem/solver/eigen.hh"]
51
52    if solverType == "cg":
53        operator = lambda df,_: "Dune::Fem::EigenCGInverseOperator< " + df + " >"
54    elif solverType == "bicgstab":
55        operator = lambda df,_: "Dune::Fem::EigenBiCGStabInverseOperator< " + df + " >"
56    else:
57        raise ValueError("wrong krylov solver - only cg,bicgstab available")
58
59    includes, typeName = solvers(includes,storage,operator)
60    return "eigen",includes,typeName, {}
61
62def viennaclsolver(storage,solverType="gmres"):
63    includes = ["dune/fem/solver/viennacl.hh"]
64
65    if solverType == "cg":
66        operator = lambda df,_: "Dune::Fem::ViennaCLCGInverseOperator< " + df + " >"
67    elif solverType == "gmres":
68        operator = lambda df,_: "Dune::Fem::ViennalCLGMResInverseOperator< " + df + " >"
69    elif solverType == "bicgstab":
70        operator = lambda df,_: "Dune::Fem::ViennalCLBiCGStabInverseOperator< " + df + " >"
71    else:
72        raise ValueError("wrong krylov solver - only cg,gmres,bicgstab available")
73
74    includes, typeName = solvers(includes,storage,operator)
75    return "viennacl",includes,typeName, {}
76
77def petscsolver(storage,solverType=None):
78    includes = ["dune/fem/solver/petscinverseoperators.hh"]
79
80    operator = lambda df,_: "Dune::Fem::PetscInverseOperator< " + df + " >"
81
82    includes, typeName = solvers(includes,storage,operator)
83    parameter = {"newton.linear.method":solverType} if solverType is not None else {}
84    return "petsc",includes,typeName, parameter
85
86def amgxsolver(storage,solverType="gmres"):
87    print("AMGX sovler used ")
88    includes = ["dune/fem/solver/amgxsolver.hh"]
89
90    operator = lambda df,_: "Dune::Fem::AMGXInverseOperator< " + df + " >"
91
92    includes, typeName = solvers(includes,storage,operator)
93    parameter = {"newton.linear.method":solverType} if solverType is not None else {}
94    return "amgx",includes,typeName, parameter
95