sfepy.solvers.ls module¶
- class sfepy.solvers.ls.CholeskySolver(conf, **kwargs)[source]¶
Interface to scikit-sparse.cholesky solver.
Kind: ‘ls.cholesky’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- name = 'ls.cholesky'¶
- class sfepy.solvers.ls.MUMPSParallelSolver(conf, **kwargs)[source]¶
Interface to MUMPS parallel solver.
Kind: ‘ls.mumps_par’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- memory_relaxationint (default: 20)
The percentage increase in the estimated working space.
- name = 'ls.mumps_par'¶
- class sfepy.solvers.ls.MUMPSSolver(conf, **kwargs)[source]¶
Interface to MUMPS solver.
Kind: ‘ls.mumps’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- memory_relaxationint (default: 20)
The percentage increase in the estimated working space.
- name = 'ls.mumps'¶
- class sfepy.solvers.ls.MultiProblem(conf, context=None, **kwargs)[source]¶
Conjugate multiple problems.
Allows to define conjugate multiple problems.
Kind: ‘ls.cm_pb’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- method{‘auto’, ‘umfpack’, ‘superlu’} (default: ‘auto’)
The actual solver to use.
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- otherslist
The list of auxiliary problem definition files.
- coupling_variableslist
The list of coupling variables.
- name = 'ls.cm_pb'¶
- class sfepy.solvers.ls.PETScKrylovSolver(conf, comm=None, context=None, **kwargs)[source]¶
PETSc Krylov subspace solver.
The solver supports parallel use with a given MPI communicator (see comm argument of
PETScKrylovSolver.__init__()
) and allows passing in PETSc matrices and vectors. Returns a (global) PETSc solution vector instead of a (local) numpy array, when given a PETSc right-hand side vector.The solver and preconditioner types are set upon the solver object creation. Tolerances can be overridden when called by passing a conf object.
Convergence is reached when rnorm < max(eps_r * rnorm_0, eps_a), where, in PETSc, rnorm is by default the norm of preconditioned residual.
Kind: ‘ls.petsc’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- methodstr (default: ‘cg’)
The actual solver to use.
- setup_precondcallable
User-supplied function for the preconditioner initialization/setup. It is called as setup_precond(mtx, context), where mtx is the matrix, context is a user-supplied context, and should return an object with setUp(self, pc) and apply(self, pc, x, y) methods. Has precedence over the precond/sub_precond parameters.
- precondstr (default: ‘icc’)
The preconditioner.
- sub_precondstr (default: ‘none’)
The preconditioner for matrix blocks (in parallel runs).
- precond_side{‘left’, ‘right’, ‘symmetric’, None}
The preconditioner side.
- i_maxint (default: 100)
The maximum number of iterations.
- eps_afloat (default: 1e-08)
The absolute tolerance for the residual.
- eps_rfloat (default: 1e-08)
The relative tolerance for the residual.
- eps_dfloat (default: 100000.0)
The divergence tolerance for the residual.
- force_reusebool (default: False)
If True, skip the check whether the KSP solver object corresponds to the mtx argument: it is always reused.
- **
Additional parameters supported by the method. Can be used to pass all PETSc options supported by
petsc.Options()
.
- name = 'ls.petsc'¶
- class sfepy.solvers.ls.PyAMGKrylovSolver(conf, context=None, **kwargs)[source]¶
Interface to PyAMG Krylov solvers.
Kind: ‘ls.pyamg_krylov’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- methodstr (default: ‘cg’)
The actual solver to use.
- setup_precondcallable (default: <function PyAMGKrylovSolver.<lambda> at 0x75f65056a3b0>)
User-supplied function for the preconditioner initialization/setup. It is called as setup_precond(mtx, context), where mtx is the matrix, context is a user-supplied context, and should return one of {sparse matrix, dense matrix, LinearOperator}.
- callbackcallable
User-supplied function to call after each iteration. It is called as callback(xk), where xk is the current solution vector, except the gmres method, where the argument is the residual norm.
- i_maxint (default: 100)
The maximum number of iterations.
- eps_rfloat (default: 1e-08)
The relative tolerance for the residual.
- **
Additional parameters supported by the method.
- name = 'ls.pyamg_krylov'¶
- class sfepy.solvers.ls.PyAMGSolver(conf, **kwargs)[source]¶
Interface to PyAMG solvers.
The method parameter can be one of: ‘smoothed_aggregation_solver’, ‘ruge_stuben_solver’. The accel parameter specifies the Krylov solver name, that is used as an accelerator for the multigrid solver.
Kind: ‘ls.pyamg’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- methodstr (default: ‘smoothed_aggregation_solver’)
The actual solver to use.
- accelstr
The accelerator.
- callbackcallable
User-supplied function to call after each iteration. It is called as callback(xk), where xk is the current solution vector, except the gmres accelerator, where the argument is the residual norm.
- i_maxint (default: 100)
The maximum number of iterations.
- eps_rfloat (default: 1e-08)
The relative tolerance for the residual.
- force_reusebool (default: False)
If True, skip the check whether the MG solver object corresponds to the mtx argument: it is always reused.
- **
Additional parameters supported by the method. Use the ‘method:’ prefix for arguments of the method construction function (e.g. ‘method:max_levels’ : 5), and the ‘solve:’ prefix for the subsequent solver call.
- name = 'ls.pyamg'¶
- class sfepy.solvers.ls.RMMSolver(conf, context=None, **kwargs)[source]¶
Special solver for explicit transient elastodynamics.
The solver uses the reciprocal mass matrix algorithm [1], [2] to directly construct a sparse inverse mass matrix. Instead of solving a linear system, calling the solver simply performs a sparse matrix multiplication.
Limitations:
Assumes that the density is constant in time.
Uses the direct EBC application, i.e., no EBC projection matrix.
[1]González, J.A., Kolman, R., Cho, S.S., Felippa, C.A., Park, K.C., 2018. Inverse mass matrix via the method of localized Lagrange multipliers. International Journal for Numerical Methods in Engineering 113, 277–295. https://doi.org/10.1002/nme.5613
[2]González, J.A., Kopačka, J., Kolman, R., Cho, S.S., Park, K.C., 2019. Inverse mass matrix for isogeometric explicit transient analysis via the method of localized Lagrange multipliers. International Journal for Numerical Methods in Engineering 117, 939–966. https://doi.org/10.1002/nme.5986
Kind: ‘ls.rmm’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- rmm_termstr
The RMM term definition, see
MassTerm
.- debugbool (default: False)
If True, run in debug mode.
- name = 'ls.rmm'¶
- class sfepy.solvers.ls.SchurMumps(conf, **kwargs)[source]¶
Mumps Schur complement solver.
Kind: ‘ls.schur_mumps’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- memory_relaxationint (default: 20)
The percentage increase in the estimated working space.
- schur_variableslist
The list of Schur variables.
- name = 'ls.schur_mumps'¶
- class sfepy.solvers.ls.ScipyDirect(conf, method=None, **kwargs)[source]¶
Direct sparse solver from SciPy.
Kind: ‘ls.scipy_direct’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- method{‘auto’, ‘umfpack’, ‘superlu’} (default: ‘auto’)
The actual solver to use.
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- name = 'ls.scipy_direct'¶
- class sfepy.solvers.ls.ScipyIterative(conf, context=None, **kwargs)[source]¶
Interface to SciPy iterative solvers.
The eps_r tolerance is both absolute and relative - the solvers stop when either the relative or the absolute residual is below it.
Kind: ‘ls.scipy_iterative’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- methodstr (default: ‘cg’)
The actual solver to use.
- setup_precondcallable (default: <function ScipyIterative.<lambda> at 0x75f650569fc0>)
User-supplied function for the preconditioner initialization/setup. It is called as setup_precond(mtx, context), where mtx is the matrix, context is a user-supplied context, and should return one of {sparse matrix, dense matrix, LinearOperator}.
- callbackcallable
User-supplied function to call after each iteration. It is called as callback(xk), where xk is the current solution vector, except the gmres method, where the argument is the residual.
- i_maxint (default: 100)
The maximum number of iterations.
- eps_afloat (default: 1e-08)
The absolute tolerance for the residual.
- eps_rfloat (default: 1e-08)
The relative tolerance for the residual.
- **
Additional parameters supported by the method.
- name = 'ls.scipy_iterative'¶
- class sfepy.solvers.ls.ScipySuperLU(conf, **kwargs)[source]¶
SuperLU - direct sparse solver from SciPy.
Kind: ‘ls.scipy_superlu’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- name = 'ls.scipy_superlu'¶
- class sfepy.solvers.ls.ScipyUmfpack(conf, **kwargs)[source]¶
UMFPACK - direct sparse solver from SciPy.
Kind: ‘ls.scipy_umfpack’
For common configuration parameters, see
Solver
.Specific configuration parameters:
- Parameters:
- use_presolvebool (default: False)
If True, pre-factorize the matrix.
- use_mtx_digestbool (default: True)
If True, determine automatically a reused matrix using its SHA1 digest. If False, .clear() has to be called manually whenever the matrix changes - expert use only!
- name = 'ls.scipy_umfpack'¶
- sfepy.solvers.ls.petsc_call(call)[source]¶
Decorator handling argument preparation and timing for PETSc-based linear solvers.