dune-istl  2.7.1
Public Types | Public Member Functions | Protected Types | Protected Member Functions | List of all members
Dune::Cholmod< BlockVector< FieldVector< T, k >, A > > Class Template Reference

Dune wrapper for SuiteSparse/CHOLMOD solver. More...

#include <dune/istl/cholmod.hh>

Inheritance diagram for Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >:
Inheritance graph

Public Types

using X = BlockVector< FieldVector< T, k >, A >
 
using B = BlockVector< FieldVector< T, k >, A >
 
typedef BlockVector< FieldVector< T, k >, A > domain_type
 Type of the domain of the operator to be inverted. More...
 
typedef BlockVector< FieldVector< T, k >, A > range_type
 Type of the range of the operator to be inverted. More...
 
typedef X::field_type field_type
 The field type of the operator. More...
 
typedef FieldTraits< field_type >::real_type real_type
 The real type of the field type (is the same if using real numbers, but differs for std::complex) More...
 
typedef Simd::Scalar< real_typescalar_real_type
 scalar type underlying the field_type More...
 

Public Member Functions

 Cholmod ()
 Default constructor. More...
 
 ~Cholmod ()
 Destructor. More...
 
 Cholmod (const Cholmod &)=delete
 
Cholmodoperator= (const Cholmod &)=delete
 
void apply (X &x, B &b, double reduction, InverseOperatorResult &res)
 simple forward to apply(X&, Y&, InverseOperatorResult&) More...
 
void apply (X &x, B &b, InverseOperatorResult &res)
 solve the linear system Ax=b (possibly with respect to some ignore field) More...
 
void setMatrix (const BCRSMatrix< FieldMatrix< T, k, k >> &matrix)
 Set matrix without ignore nodes. More...
 
template<class Ignore >
void setMatrix (const BCRSMatrix< FieldMatrix< T, k, k >> &matrix, const Ignore *ignore)
 Set matrix and ignore nodes. More...
 
virtual SolverCategory::Category category () const
 Category of the solver (see SolverCategory::Category) More...
 

Protected Types

enum  
 

Protected Member Functions

void printHeader (std::ostream &s) const
 helper function for printing header of solver output More...
 
void printOutput (std::ostream &s, const CountType &iter, const DataType &norm, const DataType &norm_old) const
 helper function for printing solver output More...
 
void printOutput (std::ostream &s, const CountType &iter, const DataType &norm) const
 helper function for printing solver output More...
 

Detailed Description

template<class T, class A, int k>
class Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >

Dune wrapper for SuiteSparse/CHOLMOD solver.

This class implements an InverseOperator between BlockVector types

Member Typedef Documentation

◆ B

template<class T , class A , int k>
using Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::B = BlockVector<FieldVector<T,k>, A>

◆ domain_type

typedef BlockVector< FieldVector< T, k >, A > Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::domain_type
inherited

Type of the domain of the operator to be inverted.

◆ field_type

typedef X::field_type Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::field_type
inherited

The field type of the operator.

◆ range_type

typedef BlockVector< FieldVector< T, k >, A > Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::range_type
inherited

Type of the range of the operator to be inverted.

◆ real_type

typedef FieldTraits<field_type>::real_type Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::real_type
inherited

The real type of the field type (is the same if using real numbers, but differs for std::complex)

◆ scalar_real_type

typedef Simd::Scalar<real_type> Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::scalar_real_type
inherited

scalar type underlying the field_type

◆ X

template<class T , class A , int k>
using Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::X = BlockVector<FieldVector<T,k>, A>

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
protectedinherited

Constructor & Destructor Documentation

◆ Cholmod() [1/2]

template<class T , class A , int k>
Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::Cholmod ( )
inline

Default constructor.

Calls the the cholmod runtime, see CHOLMOD doc

◆ ~Cholmod()

template<class T , class A , int k>
Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::~Cholmod ( )
inline

Destructor.

Free space and calls cholmod_finish, see CHOLMOD doc

◆ Cholmod() [2/2]

template<class T , class A , int k>
Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::Cholmod ( const Cholmod< BlockVector< FieldVector< T, k >, A > > &  )
delete

Member Function Documentation

◆ apply() [1/2]

template<class T , class A , int k>
void Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::apply ( X x,
B b,
double  reduction,
InverseOperatorResult res 
)
inlinevirtual

simple forward to apply(X&, Y&, InverseOperatorResult&)

Implements Dune::InverseOperator< BlockVector< FieldVector< T, k >, A >, BlockVector< FieldVector< T, k >, A > >.

◆ apply() [2/2]

template<class T , class A , int k>
void Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::apply ( X x,
B b,
InverseOperatorResult res 
)
inlinevirtual

solve the linear system Ax=b (possibly with respect to some ignore field)

The method assumes that setMatrix() was called before In the case of a given ignore field the corresponding entries of both in x and b will stay untouched in this method.

Implements Dune::InverseOperator< BlockVector< FieldVector< T, k >, A >, BlockVector< FieldVector< T, k >, A > >.

◆ category()

template<class T , class A , int k>
virtual SolverCategory::Category Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::category ( ) const
inlinevirtual

◆ operator=()

template<class T , class A , int k>
Cholmod& Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::operator= ( const Cholmod< BlockVector< FieldVector< T, k >, A > > &  )
delete

◆ printHeader()

void Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::printHeader ( std::ostream &  s) const
inlineprotectedinherited

helper function for printing header of solver output

◆ printOutput() [1/2]

void Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::printOutput ( std::ostream &  s,
const CountType &  iter,
const DataType &  norm 
) const
inlineprotectedinherited

helper function for printing solver output

◆ printOutput() [2/2]

void Dune::InverseOperator< BlockVector< FieldVector< T, k >, A > , BlockVector< FieldVector< T, k >, A > >::printOutput ( std::ostream &  s,
const CountType &  iter,
const DataType &  norm,
const DataType &  norm_old 
) const
inlineprotectedinherited

helper function for printing solver output

◆ setMatrix() [1/2]

template<class T , class A , int k>
void Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::setMatrix ( const BCRSMatrix< FieldMatrix< T, k, k >> &  matrix)
inline

Set matrix without ignore nodes.

This method forwards a nullptr to the setMatrix method with ignore nodes

◆ setMatrix() [2/2]

template<class T , class A , int k>
template<class Ignore >
void Dune::Cholmod< BlockVector< FieldVector< T, k >, A > >::setMatrix ( const BCRSMatrix< FieldMatrix< T, k, k >> &  matrix,
const Ignore *  ignore 
)
inline

Set matrix and ignore nodes.

The input matrix is copied to CHOLMOD compatible form. It is possible to ignore some degrees of freedom, provided an ignore field is given with same block structure like the BlockVector template of the class.

The ignore field causes the method to ignore both rows and cols of the matrix and therefore operates only on the reduced quadratic matrix. In case of, e.g., Dirichlet values the user has to take care of proper adjusting of the rhs before calling apply().

Decomposing the matrix at the end takes a lot of time

Parameters
[in]matrixMatrix to be decomposed. In BCRS compatible form
[in]ignorePointer to a compatible BitVector

The documentation for this class was generated from the following file: