Solvers¶
Solvers for linear, nonlinear, and time-dependent problems.
-
pyiga.solvers.
GaussSeidelSmoother
(iterations=1, sweep='forward')¶ Gauss-Seidel smoother.
By default, iterations is 1. The direction to be used is specified by sweep and may be either ‘forward’, ‘backward’, or ‘symmetric’.
-
exception
pyiga.solvers.
NoConvergenceError
(method, num_iter, last_iterate)¶
-
pyiga.solvers.
OperatorSmoother
(S)¶ A smoother which applies an arbitrary operator S to the residual and uses the result as an update, i.e.,
\[u \leftarrow u + S(f - Au).\]
-
pyiga.solvers.
SequentialSmoother
(smoothers)¶ Smoother which applies several smoothers in sequence.
-
pyiga.solvers.
crank_nicolson
(M, F, J, x, tau, t_end, *, t0=0.0, progress=False)¶ Solve a time-dependent problem using the Crank-Nicolson method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- t0 (float) – the initial time; 0 by default
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
dirk34
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the DIRK34 Runge-Kutta method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
esdirk23
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the ESDIRK23 Runge-Kutta method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
esdirk34
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the ESDIRK34 Runge-Kutta method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
fastdiag_solver
(KM)¶ The fast diagonalization solver as described in [Sangalli, Tani 2016].
Parameters: KM – a sequence of length d (dimension of the problem) containing pairs of symmetric matrices (K_i, M_i) Returns: A LinearOperator which realizes the inverse of the generalized Laplacian matrix described by the input matrices.
-
pyiga.solvers.
gauss_seidel
(A, x, b, iterations=1, indices=None, sweep='forward')¶ Perform Gauss-Seidel relaxation on the linear system Ax=b, updating x in place.
Parameters: - A – the matrix; either sparse or dense, but should be in CSR format if sparse
- x – the current guess for the solution
- b – the right-hand side
- iterations – the number of iterations to perform
- indices – if given, relaxation is only performed on this list of indices
- sweep – the direction; either ‘forward’, ‘backward’, or ‘symmetric’
-
pyiga.solvers.
iterative_solve
(step, A, f, x0=None, active_dofs=None, tol=1e-08, maxiter=5000)¶ Solve the linear system Ax=f using a basic iterative method.
Parameters: - step (callable) – a function which performs the update x_old -> x_new for the iterative method
- A – matrix or linear operator describing the linear system of equations
- f (ndarray) – the right-hand side
- x0 – the starting vector; 0 is used if not specified
- active_dofs (list or ndarray) – list of active dofs on which the residual is computed. Useful for eliminating Dirichlet dofs without changing the matrix. If not specified, all dofs are active.
- tol (float) – the desired reduction in the Euclidean norm of the residual relative to the starting residual
- maxiter (int) – the maximum number of iterations
Returns: a pair (x, iterations) containing the solution and the number of iterations performed. If maxiter was reached without convergence, the returned number of iterations is infinite.
-
pyiga.solvers.
newton
(F, J, x0, atol=1e-06, rtol=1e-06, maxiter=100, freeze_jac=1)¶ Solve the nonlinear problem F(x) == 0 using Newton iteration.
Parameters: - F (function) – function computing the residual of the nonlinear equation
- J (function) – function computing the Jacobian matrix of F
- x0 (ndarray) – the initial guess as a vector
- atol (float) – absolute tolerance for the norm of the residual
- rtol (float) – relative tolerance with respect to the initial residual
- maxiter (int) – the maximum number of iterations
- freeze_jac (int) – if >1, the Jacobian is only updated every freeze_jac steps
Returns: ndarray – a vector x which approximately satisfies F(x) == 0
-
pyiga.solvers.
rodasp
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the RODASP Rosenbrock method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
ros3p
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the ROS3P Rosenbrock method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
ros3pw
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the ROS3PW Rosenbrock method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
rosi2p1
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the ROSI2P1 Rosenbrock method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
rowdaind2
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the ROWDAIND2 Rosenbrock method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
sdirk21
(M, F, J, x, tau0, t_end, tol, *, t0=0.0, step_factor=0.9, progress=False)¶ Solve a time-dependent problem using the SDIRK21 (Ellsiepen) Runge-Kutta method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- tol (float) – error tolerance for choosing the adaptive time step; if None, use constant time steps
- t0 (float) – the initial time; 0 by default
- step_factor (float) – the safety factor for choosing the step size
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
sdirk3
(M, F, J, x, tau, t_end, *, t0=0.0, progress=False)¶ Solve a time-dependent problem using the SDIRK3 Runge-Kutta method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- t0 (float) – the initial time; 0 by default
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
sdirk3_b
(M, F, J, x, tau, t_end, *, t0=0.0, progress=False)¶ Solve a time-dependent problem using the SDIRK3 (alternate) Runge-Kutta method.
Parameters: - M (matrix) – the mass matrix
- F (function) – the right-hand side
- J (function) – function computing the Jacobian of F
- x (vector) – the initial value
- tau0 (float) – the initial time step
- t_end (float) – the final time up to which to integrate
- t0 (float) – the initial time; 0 by default
- progress (bool) – whether to show a progress bar
Returns: A tuple (times, solutions), where times is a list of increasing times in the interval (t0, t_end), and solutions is a list of vectors which contains the computed solutions at these times.
-
pyiga.solvers.
solve_hmultigrid
(hs, A, f, strategy='cell_supp', smoother='gs', smooth_steps=2, tol=1e-08, maxiter=5000)¶ Solve a linear scalar problem in a hierarchical spline space using local multigrid.
Parameters: - hs – the
HSpace
which describes the HB- or THB-spline space - A – the matrix describing the discretization of the problem
- f – the right-hand side vector
- strategy (string) –
how to choose the smoothing sets. Valid options are
"new"
: only the new dofs per level"trunc"
: all dofs which interact via the truncation operator"cell_supp"
: all dofs whose support intersects that of the new ones (support extension)"func_supp"
- smoother (string) –
the multigrid smoother to use. Valid options are
"gs"
: forward Gauss-Seidel for pre-smoothing, backward Gauss-Seidel for post-smoothing"forward_gs"
: always use forward Gauss-Seidel"backward_gs"
: always use backward Gauss-Seidel"symmetric_gs"
: use complete symmetric Gauss-Seidel sweep for both pre- and post-smoothing"exact"
: use an exact direct solver as a pre-smoother (no post-smoothing)
- smooth_steps (int) – the number of pre- and post-smoothing steps
- tol (float) – the desired reduction in the residual
- maxiter (int) – the maximum number of iterations
Returns: a pair (x, iterations) containing the solution and the number of iterations performed. If maxiter was reached without convergence, the returned number of iterations is infinite.
- hs – the
-
pyiga.solvers.
twogrid
(A, f, P, smoother, u0=None, tol=1e-08, smooth_steps=2, maxiter=1000)¶ Generic two-grid method with arbitrary smoother.
Parameters: - A – stiffness matrix on fine grid
- f – right-hand side
- P – prolongation matrix from coarse to fine grid
- smoother – a function with arguments (A,u,f) which applies one smoothing iteration in-place to u
- u0 – starting value; 0 if not given
- tol – desired reduction relative to initial residual
- smooth_steps – number of smoothing steps
- maxiter – maximum number of iterations
Returns: ndarray – the computed solution to the equation Au = f