sfepy.solvers.ls module

class sfepy.solvers.ls.MultiProblem(conf, problem, **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.

presolve : bool (default: False)

If True, pre-factorize the matrix.

warn : bool (default: True)

If True, allow warnings.

others : list

The list of auxiliary problem definition files.

coupling_variables : list

The list of coupling variables.

name = 'ls.cm_pb'
sparse_submat(Ad, Ar, Ac, gr, gc, S)[source]

A[gr,gc] = S

class sfepy.solvers.ls.PETScKrylovSolver(conf, comm=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:

method : str (default: ‘cg’)

The actual solver to use.

precond : str (default: ‘icc’)

The preconditioner.

sub_precond : str

The preconditioner for matrix blocks (in parallel runs).

precond_side : {‘left’, ‘right’, ‘symmetric’, None}

The preconditioner side.

i_max : int (default: 100)

The maximum number of iterations.

eps_a : float (default: 1e-08)

The absolute tolerance for the residual.

eps_r : float (default: 1e-08)

The relative tolerance for the residual.

eps_d : float (default: 100000.0)

The divergence tolerance for the residual.

create_ksp(comm=None)[source]
create_petsc_matrix(mtx, comm=None)[source]
name = 'ls.petsc'
set_field_split(field_ranges, comm=None)[source]

Setup local PETSc ranges for fields to be used with ‘fieldsplit’ preconditioner.

This function must be called before solving the linear system.

class sfepy.solvers.ls.PETScParallelKrylovSolver(conf, comm=None, **kwargs)[source]

PETSc Krylov subspace solver able to run in parallel by storing the system to disk and running a separate script via mpiexec.

The solver and preconditioner types are set upon the solver object creation. Tolerances can be overriden 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_parallel’

For common configuration parameters, see Solver.

Specific configuration parameters:

Parameters:

method : str (default: ‘cg’)

The actual solver to use.

precond : str (default: ‘icc’)

The preconditioner.

sub_precond : str

The preconditioner for matrix blocks (in parallel runs).

precond_side : {‘left’, ‘right’, ‘symmetric’, None}

The preconditioner side.

i_max : int (default: 100)

The maximum number of iterations.

eps_a : float (default: 1e-08)

The absolute tolerance for the residual.

eps_r : float (default: 1e-08)

The relative tolerance for the residual.

eps_d : float (default: 100000.0)

The divergence tolerance for the residual.

log_dir : str (default: ‘.’)

The directory for storing logs.

n_proc : int (default: 1)

The number of processes.

sub_precond : str (default: ‘icc’)

The preconditioner for matrix blocks.

name = 'ls.petsc_parallel'
class sfepy.solvers.ls.PyAMGSolver(conf, **kwargs)[source]

Interface to PyAMG solvers.

Kind: ‘ls.pyamg’

For common configuration parameters, see Solver.

Specific configuration parameters:

Parameters:

method : str (default: ‘smoothed_aggregation_solver’)

The actual solver to use.

accel : str

The accelerator.

i_max : int (default: 100)

The maximum number of iterations.

eps_r : float (default: 1e-08)

The relative tolerance for the residual.

name = 'ls.pyamg'
class sfepy.solvers.ls.SchurComplement(conf, **kwargs)[source]

Schur complement.

Solution of the linear system

\left[ \begin{array}{cc} A & B \\ C & D \end{array} \right] \cdot \left[ \begin{array}{c} u \\ v \end{array} \right] = \left[ \begin{array}{c} f \\ g \end{array} \right]

is obtained by solving the following equation:

(D - C A^{-1} B) \cdot v = g - C A^{-1} f

variable(s) u are specified in “eliminate” list, variable(s) v are specified in “keep” list,

See: http://en.wikipedia.org/wiki/Schur_complement

Kind: ‘ls.schur_complement’

For common configuration parameters, see Solver.

Specific configuration parameters:

Parameters:

method : {‘auto’, ‘umfpack’, ‘superlu’} (default: ‘auto’)

The actual solver to use.

presolve : bool (default: False)

If True, pre-factorize the matrix.

warn : bool (default: True)

If True, allow warnings.

eliminate : list

The list of variables to eliminate.

keep : list

The list of variables to keep.

name = 'ls.schur_complement'
static schur_fun(res, mtx, rhs, nn)[source]
class sfepy.solvers.ls.SchurGeneralized(conf, problem, **kwargs)[source]

Generalized Schur complement.

Defines the matrix blocks and calls user defined function.

Kind: ‘ls.schur_generalized’

For common configuration parameters, see Solver.

Specific configuration parameters:

Parameters:

method : {‘auto’, ‘umfpack’, ‘superlu’} (default: ‘auto’)

The actual solver to use.

presolve : bool (default: False)

If True, pre-factorize the matrix.

warn : bool (default: True)

If True, allow warnings.

blocks : dict

The description of blocks: {block_name1 : [variable_name1, ...], ...}.

function : callable

The user defined function.

name = 'ls.schur_generalized'
class sfepy.solvers.ls.ScipyDirect(conf, **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.

presolve : bool (default: False)

If True, pre-factorize the matrix.

warn : bool (default: True)

If True, allow warnings.

name = 'ls.scipy_direct'
presolve(mtx)[source]
class sfepy.solvers.ls.ScipyIterative(conf, **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.

A preconditioner can be anything that the SciPy solvers accept (sparse matrix, dense matrix, LinearOperator).

Kind: ‘ls.scipy_iterative’

For common configuration parameters, see Solver.

Specific configuration parameters:

Parameters:

method : str (default: ‘cg’)

The actual solver to use.

precond : {sparse matrix, dense matrix, LinearOperator}

The preconditioner.

callback : function

User-supplied function to call after each iteration. It is called as callback(xk), where xk is the current solution vector.

i_max : int (default: 100)

The maximum number of iterations.

eps_r : float (default: 1e-08)

The relative or absolute tolerance for the residual.

name = 'ls.scipy_iterative'
sfepy.solvers.ls.petsc_call(call)[source]

Decorator handling argument preparation and timing for PETSc-based linear solvers.

sfepy.solvers.ls.solve(mtx, rhs, solver_class=None, solver_conf=None)[source]

Solve the linear system with the matrix mtx and the right-hand side rhs.

Convenience wrapper around the linear solver classes below.

sfepy.solvers.ls.standard_call(call)[source]

Decorator handling argument preparation and timing for linear solvers.