DGtal  1.4.beta
Public Types | Public Member Functions | Static Public Member Functions
DGtal::DirichletConditions< TLinearAlgebraBackend > Class Template Reference

Aim: A helper class to solve a system with Dirichlet boundary conditions. More...

#include <DGtal/math/linalg/DirichletConditions.h>

Public Types

typedef TLinearAlgebraBackend LinearAlgebraBackend
 
typedef LinearAlgebraBackend::DenseVector::Index Index
 
typedef LinearAlgebraBackend::DenseVector::Scalar Scalar
 
typedef LinearAlgebraBackend::DenseVector DenseVector
 
typedef LinearAlgebraBackend::IntegerVector IntegerVector
 
typedef LinearAlgebraBackend::DenseMatrix DenseMatrix
 
typedef LinearAlgebraBackend::SparseMatrix SparseMatrix
 
typedef LinearAlgebraBackend::Triplet Triplet
 

Public Member Functions

 BOOST_CONCEPT_ASSERT ((concepts::CDynamicVector< DenseVector >))
 
 BOOST_CONCEPT_ASSERT ((concepts::CDynamicMatrix< DenseMatrix >))
 
 BOOST_CONCEPT_ASSERT ((concepts::CDynamicMatrix< SparseMatrix >))
 
 BOOST_CONCEPT_ASSERT ((concepts::CLinearAlgebra< DenseVector, SparseMatrix >))
 
 BOOST_CONCEPT_ASSERT ((concepts::CLinearAlgebra< DenseVector, DenseMatrix >))
 

Static Public Member Functions

static SparseMatrix dirichletOperator (const SparseMatrix &A, const IntegerVector &p)
 
static DenseVector dirichletVector (const SparseMatrix &A, const DenseVector &b, const IntegerVector &p, const DenseVector &u)
 
static DenseVector dirichletSolution (const DenseVector &xd, const IntegerVector &p, const DenseVector &u)
 
static IntegerVector nullBoundaryVector (const DenseVector &b)
 

Detailed Description

template<typename TLinearAlgebraBackend>
class DGtal::DirichletConditions< TLinearAlgebraBackend >

Aim: A helper class to solve a system with Dirichlet boundary conditions.

Description of template class 'DirichletConditions'

Typically you have a system of the form \( A x = b \), where you wish to set up Dirichlet boundary conditions \( u \) at some places \( p \) where \( p=1 \) ( \( p=0 \) within the domain). In this case, you solve a modified smaller system \( A_d x_d = b_d \).

A typical usage is:

typedef DirichletConditions< EigenLinearAlgebraBackend > DC;
DC::SparseMatrix A = ...; // the matrix of the linear system
DC::DenseVector b = ...; // contains right-hand side and boundary conditions
// set up boundary locations and values
DC::IntegerVector p = IntegerVector::Zero( b.rows() );
DC::DenseVector u = DenseVector::Zero( b.rows() );
std::vector<DC::Index> boundary_nodes = { 12, 17, 25, 38, ... };
std::vector<DC::Scalar> boundary_values = { 10.0, 20.0, 15.0, 7.0, ... };
for ( int i = 0; i < boundary_nodes.size(); i++ )
{
auto j = boundary_nodes[ i ];
p[ j ] = 1;
u[ j ] = boundary_values[ i ];
}
DC::SparseMatrix A_d = DC::dirichletOperator( A, p );
DC::DenseVector b_d = DC::dirichletVector ( A, b, p, u );
solver.compute( A_d ); // prefactorization
ASSERT(solver.info()==Eigen::Success);
DC::DenseVector x_d = solver.solve( b_d ); // solve Dirichlet system
DC::DenseVector x = DC::dirichletSolution( x_d, p, u ); // get whole solution

It is a model of boost::CopyConstructible, boost::DefaultConstructible, boost::Assignable.

Template Parameters
TLinearAlgebraBackendlinear algebra backend used (i.e. EigenLinearAlgebraBackend).

Definition at line 97 of file DirichletConditions.h.

Member Typedef Documentation

◆ DenseMatrix

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::DenseMatrix DGtal::DirichletConditions< TLinearAlgebraBackend >::DenseMatrix

Definition at line 106 of file DirichletConditions.h.

◆ DenseVector

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::DenseVector DGtal::DirichletConditions< TLinearAlgebraBackend >::DenseVector

Definition at line 104 of file DirichletConditions.h.

◆ Index

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::DenseVector::Index DGtal::DirichletConditions< TLinearAlgebraBackend >::Index

Definition at line 102 of file DirichletConditions.h.

◆ IntegerVector

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::IntegerVector DGtal::DirichletConditions< TLinearAlgebraBackend >::IntegerVector

Definition at line 105 of file DirichletConditions.h.

◆ LinearAlgebraBackend

template<typename TLinearAlgebraBackend >
typedef TLinearAlgebraBackend DGtal::DirichletConditions< TLinearAlgebraBackend >::LinearAlgebraBackend

Definition at line 100 of file DirichletConditions.h.

◆ Scalar

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::DenseVector::Scalar DGtal::DirichletConditions< TLinearAlgebraBackend >::Scalar

Definition at line 103 of file DirichletConditions.h.

◆ SparseMatrix

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::SparseMatrix DGtal::DirichletConditions< TLinearAlgebraBackend >::SparseMatrix

Definition at line 107 of file DirichletConditions.h.

◆ Triplet

template<typename TLinearAlgebraBackend >
typedef LinearAlgebraBackend::Triplet DGtal::DirichletConditions< TLinearAlgebraBackend >::Triplet

Definition at line 108 of file DirichletConditions.h.

Member Function Documentation

◆ BOOST_CONCEPT_ASSERT() [1/5]

template<typename TLinearAlgebraBackend >
DGtal::DirichletConditions< TLinearAlgebraBackend >::BOOST_CONCEPT_ASSERT ( (concepts::CDynamicMatrix< DenseMatrix >)  )

◆ BOOST_CONCEPT_ASSERT() [2/5]

template<typename TLinearAlgebraBackend >
DGtal::DirichletConditions< TLinearAlgebraBackend >::BOOST_CONCEPT_ASSERT ( (concepts::CDynamicMatrix< SparseMatrix >)  )

◆ BOOST_CONCEPT_ASSERT() [3/5]

template<typename TLinearAlgebraBackend >
DGtal::DirichletConditions< TLinearAlgebraBackend >::BOOST_CONCEPT_ASSERT ( (concepts::CDynamicVector< DenseVector >)  )

◆ BOOST_CONCEPT_ASSERT() [4/5]

template<typename TLinearAlgebraBackend >
DGtal::DirichletConditions< TLinearAlgebraBackend >::BOOST_CONCEPT_ASSERT ( (concepts::CLinearAlgebra< DenseVector, DenseMatrix >)  )

◆ BOOST_CONCEPT_ASSERT() [5/5]

template<typename TLinearAlgebraBackend >
DGtal::DirichletConditions< TLinearAlgebraBackend >::BOOST_CONCEPT_ASSERT ( (concepts::CLinearAlgebra< DenseVector, SparseMatrix >)  )

◆ dirichletOperator()

template<typename TLinearAlgebraBackend >
static SparseMatrix DGtal::DirichletConditions< TLinearAlgebraBackend >::dirichletOperator ( const SparseMatrix A,
const IntegerVector p 
)
inlinestatic

Typically you have a system of the form \( A x = b \), where you wish to set up Dirichlet boundary conditions \( u \) at some places \( p \) where \( p=1 \) ( \( p=0 \) within the domain). In this case, you solve a modified smaller system \( A_d x_d = b_d \).

Parameters
Athe matrix representing the linear operator in Ax=b.
pthe vector such that p[i]=1 whenever the i-th node is a boundary node, p[i]=0 otherwise.
Precondition
#row(A) = #col(A) = #col(p)
Returns
the (smaller) linear matrix \( A_d \) to prefactor.

Definition at line 130 of file DirichletConditions.h.

132  {
133  ASSERT( A.cols() == A.rows() );
134  ASSERT( p.rows() == A.rows() );
135  const auto n = p.rows();
136  std::vector< Index > relabeling( n );
137  Index j = 0;
138  for ( Index i = 0; i < p.rows(); i++ )
139  relabeling[ i ] = ( p[ i ] == 0 ) ? j++ : n;
140  // Building matrix
141  SparseMatrix Ap( j, j );
142  std::vector< Triplet > triplets;
143  for ( int k = 0; k < A.outerSize(); ++k )
144  for ( typename SparseMatrix::InnerIterator it( A, k ); it; ++it )
145  {
146  if ( ( relabeling[ it.row() ] != n ) && ( relabeling[ it.col() ] != n ) )
147  triplets.push_back( { relabeling[ it.row() ], relabeling[ it.col() ],
148  it.value() } );
149  }
150  Ap.setFromTriplets( triplets.cbegin(), triplets.cend() );
151  return Ap;
152  }

Referenced by DGtal::VectorsInHeat< TPolygonalCalculus >::init(), and DGtal::GeodesicsInHeat< TPolygonalCalculus >::init().

◆ dirichletSolution()

template<typename TLinearAlgebraBackend >
static DenseVector DGtal::DirichletConditions< TLinearAlgebraBackend >::dirichletSolution ( const DenseVector xd,
const IntegerVector p,
const DenseVector u 
)
inlinestatic

Typically you have a system of the form \( A x = b \), where you wish to set up Dirichlet boundary conditions \( u \) at some places \( p \) where \( p=1 \) ( \( p=0 \) within the domain). In this case, you solve a modified smaller system \( A_d x_d = b_d \).

Parameters
xdthe solution to the smaller system \( A_d x_d = b_d \).
pthe vector such that p[i]=1 whenever the i-th node is a boundary node, p[i]=0 otherwise.
uthe vector giving the constrained Dirichlet values on places such that p[i]=1.
Returns
the solution \( x \) to the original system \( A x = b \),

Definition at line 203 of file DirichletConditions.h.

206  {
207  DenseVector x( p.rows() );
208  Index j = 0;
209  for ( Index i = 0; i < p.rows(); i++ )
210  x[ i ] = ( p[ i ] == 0 ) ? xd[ j++ ] : u[ i ];
211  return x;
212  }

Referenced by DGtal::GeodesicsInHeat< TPolygonalCalculus >::compute(), and DGtal::VectorsInHeat< TPolygonalCalculus >::compute().

◆ dirichletVector()

template<typename TLinearAlgebraBackend >
static DenseVector DGtal::DirichletConditions< TLinearAlgebraBackend >::dirichletVector ( const SparseMatrix A,
const DenseVector b,
const IntegerVector p,
const DenseVector u 
)
inlinestatic

Typically you have a system of the form \( A x = b \), where you wish to set up Dirichlet boundary conditions \( u \) at some places \( p \) where \( p=1 \) ( \( p=0 \) within the domain). In this case, you solve a modified smaller system \( A_d x_d = b_d \).

Parameters
Athe matrix representing the linear operator in \( Ax=b \).
bthe vector representing the vector to solve for in \( Ax=b \).
pthe vector such that p[i]=1 whenever the i-th node is a boundary node, p[i]=0 otherwise.
uthe vector giving the constrained Dirichlet values on places such that p[i]=1.
Precondition
#row(A) = #col(A) = #col(p) = #col(b) = #col(u)
Returns
the (smaller) form \( b_d \) to solve for

Definition at line 170 of file DirichletConditions.h.

174  {
175  ASSERT( A.cols() == A.rows() );
176  ASSERT( p.rows() == A.rows() );
177  const auto n = p.rows();
178  DenseVector up = p.array().template cast<double>() * u.array();
179  DenseVector tmp = b.array() - (A * up).array();
180  std::vector< Index > relabeling( n );
181  Index j = 0;
182  for ( Index i = 0; i < p.rows(); i++ )
183  relabeling[ i ] = ( p[ i ] == 0 ) ? j++ : n;
184  DenseVector bp( j );
185  for ( Index i = 0; i < p.rows(); i++ )
186  if ( p[ i ] == 0 ) bp[ relabeling[ i ] ] = tmp[ i ];
187  return bp;
188  }

Referenced by DGtal::GeodesicsInHeat< TPolygonalCalculus >::compute(), and DGtal::VectorsInHeat< TPolygonalCalculus >::compute().

◆ nullBoundaryVector()

template<typename TLinearAlgebraBackend >
static IntegerVector DGtal::DirichletConditions< TLinearAlgebraBackend >::nullBoundaryVector ( const DenseVector b)
inlinestatic

Utility method to build a null vector that will be useful to define the boundary characteristic set for Dirichlet conditions.

Parameters
bany vector of the size of your problem (e.g., \( b \) in your problem \( Ax = b \))
Returns
an integer vector of same size as b with zero everywhere.

Definition at line 223 of file DirichletConditions.h.

224  {
225  return IntegerVector::Zero( b.rows() );
226  }

The documentation for this class was generated from the following file:
Index
SMesh::Index Index
Definition: fullConvexitySphereGeodesics.cpp:117
DGtal::DirichletConditions::DenseVector
LinearAlgebraBackend::DenseVector DenseVector
Definition: DirichletConditions.h:104
SparseMatrix
EigenLinearAlgebraBackend::SparseMatrix SparseMatrix
Definition: testHeatLaplace.cpp:50