Fork me on GitHub
 All Classes Files Functions Variables Groups Pages
Functions
Krylov Methods

Iterative Krylov methods for hermitian and non-hermitian linear systems. More...

Functions

template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::bicg (LinearOperator &A, LinearOperator &At, Vector &x, Vector &b, Monitor &monitor, Preconditioner &M, Preconditioner &Mt)
 Biconjugate Gradient method. More...
 
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::bicgstab (LinearOperator &A, Vector &x, Vector &b, Monitor &monitor, Preconditioner &M)
 Biconjugate Gradient Stabilized method. More...
 
template<class LinearOperator , class VectorType1 , class VectorType2 , class VectorType3 , class Monitor >
thrust::detail::enable_if_convertible
< typename
LinearOperator::format,
cusp::known_format >::type 
cusp::krylov::bicgstab_m (LinearOperator &A, VectorType1 &x, VectorType2 &b, VectorType3 &sigma, Monitor &monitor)
 Multi-mass Biconjugate Gradient stabilized method. More...
 
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::cg (LinearOperator &A, Vector &x, Vector &b, Monitor &monitor, Preconditioner &M)
 Conjugate Gradient method. More...
 
template<class LinearOperator , class VectorType1 , class VectorType2 , class VectorType3 , class Monitor >
thrust::detail::enable_if_convertible
< typename
LinearOperator::format,
cusp::known_format >::type 
cusp::krylov::cg_m (LinearOperator &A, VectorType1 &x, VectorType2 &b, VectorType3 &sigma, Monitor &monitor)
 Multi-mass Conjugate Gradient method. More...
 
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::cr (LinearOperator &A, Vector &x, Vector &b, Monitor &monitor, Preconditioner &M)
 Conjugate Residual method. More...
 
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::gmres (LinearOperator &A, Vector &x, Vector &b, const size_t restart, Monitor &monitor, Preconditioner &M)
 GMRES method. More...
 

Function Documentation

template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::bicg ( LinearOperator &  A,
LinearOperator &  At,
Vector &  x,
Vector &  b,
Monitor &  monitor,
Preconditioner &  M,
Preconditioner &  Mt 
)

Biconjugate Gradient method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
Atconjugate tranpose of the matrix of the linear system
xapproximate solution of the linear system
bright-hand side of the linear system
monitormontiors iteration and determines stopping conditions
Mpreconditioner for A
Mtconjugate tranpose of the preconditioner for A
Overview
Solves the linear system A x = b with preconditioner M.
Example

The following code snippet demonstrates how to use bicg to solve a 10x10 Poisson problem.

#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
// set preconditioner (identity)
// solve the linear system A x = b
// because both A and M are hermitian we can use
// them for their own conjugate transpose
cusp::krylov::bicg(A, A, x, b, monitor, M, M);
return 0;
}
See Also
monitor
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::bicgstab ( LinearOperator &  A,
Vector &  x,
Vector &  b,
Monitor &  monitor,
Preconditioner &  M 
)

Biconjugate Gradient Stabilized method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
xapproximate solution of the linear system
bright-hand side of the linear system
monitormontiors iteration and determines stopping conditions
Mpreconditioner for A
Overview

Solves the linear system A x = b with preconditioner M.

Example

The following code snippet demonstrates how to use bicgstab to solve a 10x10 Poisson problem.

#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
// set preconditioner (identity)
// solve the linear system A x = b
cusp::krylov::bicgstab(A, x, b, monitor, M);
return 0;
}
See Also
monitor
Examples:
bicgstab.cu.
template<class LinearOperator , class VectorType1 , class VectorType2 , class VectorType3 , class Monitor >
thrust::detail::enable_if_convertible<typename LinearOperator::format,cusp::known_format>::type cusp::krylov::bicgstab_m ( LinearOperator &  A,
VectorType1 &  x,
VectorType2 &  b,
VectorType3 &  sigma,
Monitor &  monitor 
)

Multi-mass Biconjugate Gradient stabilized method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
xapproximate solution of the linear system
bright-hand side of the linear system
sigmaarray of shifts
monitormontiors iteration and determines stopping conditions
Overview

This routine solves systems of the type

(A+sigma Id)x = b

for a number of different sigma, iteratively, for sparse A, without additional matrix-vector multiplication.

See Also
http://arxiv.org/abs/hep-lat/9612014
Example

The following code snippet demonstrates how to use bicgstab to solve a 10x10 Poisson problem.

#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
size_t N_s = 4;
// set sigma values
sigma[0] = 0.1;
sigma[1] = 0.5;
sigma[2] = 1.0;
sigma[3] = 5.0;
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
// solve the linear system A x = b
cusp::krylov::bicgstab_m(A, x, b, sigma, monitor);
return 0;
}
See Also
monitor
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::cg ( LinearOperator &  A,
Vector &  x,
Vector &  b,
Monitor &  monitor,
Preconditioner &  M 
)

Conjugate Gradient method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
xapproximate solution of the linear system
bright-hand side of the linear system
monitormonitors iteration and determines stopping conditions
Mpreconditioner for A
Overview
Solves the symmetric, positive-definite linear system A x = b with preconditioner M.
Note
A and M must be symmetric and positive-definite.
Example
The following code snippet demonstrates how to use cg to solve a 10x10 Poisson problem.
#include <cusp/monitor.h>
#include <cusp/krylov/cg.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
// set preconditioner (identity)
// solve the linear system A x = b
cusp::krylov::cg(A, x, b, monitor, M);
return 0;
}
See Also
monitor
Examples:
ainv.cu, cg.cu, cg_raw.cu, custom_amg.cu, diagonal.cu, monitor.cu, smoothed_aggregation.cu, stencil.cu, and verbose_monitor.cu.
template<class LinearOperator , class VectorType1 , class VectorType2 , class VectorType3 , class Monitor >
thrust::detail::enable_if_convertible<typename LinearOperator::format,cusp::known_format>::type cusp::krylov::cg_m ( LinearOperator &  A,
VectorType1 &  x,
VectorType2 &  b,
VectorType3 &  sigma,
Monitor &  monitor 
)

Multi-mass Conjugate Gradient method.

Parameters
Amatrix of the linear system
xsolutions of the system
bright-hand side of the linear system
sigmaarray of shifts
monitormonitors interation and determines stoppoing conditions
Overview

Solves the symmetric, positive-definited linear system (A+sigma) x = b for some set of constant shifts sigma for the price of the smallest shift iteratively, for sparse A, without additional matrix-vector multiplication.

See Also
http://arxiv.org/abs/hep-lat/9612014
Example

The following code snippet demonstrates how to use bicgstab to solve a 10x10 Poisson problem.

#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
size_t N_s = 4;
// set sigma values
sigma[0] = 0.1;
sigma[1] = 0.5;
sigma[2] = 1.0;
sigma[3] = 5.0;
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
// solve the linear system A x = b
cusp::krylov::cg_m(A, x, b, sigma, monitor);
return 0;
}
See Also
monitor
Examples:
cg_m.cu.
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::cr ( LinearOperator &  A,
Vector &  x,
Vector &  b,
Monitor &  monitor,
Preconditioner &  M 
)

Conjugate Residual method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
xapproximate solution of the linear system
bright-hand side of the linear system
monitormontiors iteration and determines stopping conditions
Mpreconditioner for A
Overview
Solves a linear system using the conjugate residual method
Note
A and M must be symmetric and semi-definite.
Example
The following code snippet demonstrates how to use cr to solve a 10x10 Poisson problem.
#include <cusp/monitor.h>
#include <cusp/krylov/cr.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
// set preconditioner (identity)
// solve the linear system A x = b
cusp::krylov::cr(A, x, b, monitor, M);
return 0;
}
See Also
monitor
Examples:
cr.cu.
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::krylov::gmres ( LinearOperator &  A,
Vector &  x,
Vector &  b,
const size_t  restart,
Monitor &  monitor,
Preconditioner &  M 
)

GMRES method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor such as default_monitor or verbose_monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
xapproximate solution of the linear system
bright-hand side of the linear system
restartthe method every restart inner iterations
monitormontiors iteration and determines stopping conditions
Mpreconditioner for A
Overview
Solves the nonsymmetric, linear system A x = b with preconditioner M.
Example

The following code snippet demonstrates how to use gmres to solve a 10x10 Poisson problem.

#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage for solution (x) and right hand side (b)
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<float> monitor(b, 100, 1e-6, 0, true);
int restart = 50;
// set preconditioner (identity)
// solve the linear system A x = b
cusp::krylov::gmres(A, x, b,restart, monitor, M);
return 0;
}
See Also
monitor
Examples:
gmres.cu.