KASKADE 7 development version
Classes | Public Types | Public Member Functions | Static Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule > Class Template Reference

A class for assembling Galerkin representation matrices and right hand sides for variational functionals depending on multiple variables. More...

#include <assemblerCore.hh>

Detailed Description

template<class F, class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
class Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >

A class for assembling Galerkin representation matrices and right hand sides for variational functionals depending on multiple variables.

Template parameters:

Template Parameters
FThe variational functional. This has to be a model of the VariationalFunctional concept.
BoundaryDetectorA policy class type for detecting cells incident to the boundary. Currently CachingBoundaryDetector, ForwardingBoundaryDetector, and TrivialBoundaryDetector may be used.
QuadRuleallows choice of quadrature formula. Default: provided by DUNE, some special formulas: special_quadrature.hh

After a grid refinement an assembler is not valid anymore and has to be constructed again. This is, because the assembler copies the index set at construction, which changes, if the grid changes

Definition at line 1430 of file assemblerCore.hh.

Inheritance diagram for Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >:

Classes

struct  BlockToSequence
 
struct  BlockToTriplet
 
struct  CountNonzeros
 

Public Types

typedef VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule > Self
 
typedef F Functional
 functional More...
 
typedef Functional::AnsatzVars AnsatzVariableSetDescription
 ansatz variables More...
 
typedef Functional::TestVars TestVariableSetDescription
 test variables More...
 
typedef AnsatzVariableSetDescription::Grid Grid
 grid More...
 
typedef AnsatzVariableSetDescription::GridView GridView
 The grid view on which the variables live. More...
 
typedef AnsatzVariableSetDescription::Spaces Spaces
 spaces More...
 
typedef AnsatzVariableSetDescription::IndexSet IndexSet
 index set More...
 
typedef Functional::Scalar field_type
 Underlying field type. More...
 
using Scalar = field_type
 
using MatrixBlockArray = typename AssemblyDetail::BlockArray< Policy, AnsatzVariables, TestVariables, SparseIndex >::type
 A boost::fusion sequence of AssemblyDetail::MatrixBlock elements for present matrix blocks. More...
 
typedef AssemblyDetail::RhsBlockArray< Policy, TestVariables >::type RhsBlockArray
 A boost::fusion sequence of AssemblyDetail::RhsBock elements for present rhs blocks. More...
 
typedef TestVariableSetDescription::template CoefficientVectorRepresentation ::type RhsArray
 A LinearProductSpace type of right hand side coefficient vectors. More...
 

Public Member Functions

 VariationalFunctionalAssembler (Spaces const &spaces)
 Construct an empty assembler, gridManager is implicitly obtained from the first space. More...
 
void assemble (F const &f, unsigned int flags=Assembler::EVERYTHING, int nThreads=0, bool verbose=false)
 Assembly without block filter or cell filter. More...
 
template<class BlockFilter , class CellFilter >
void assemble (F const &f, CellFilter const &cellFilter, unsigned int flags=Assembler::EVERYTHING, int nTasks=0, bool verbose=false)
 Create data in assembler. More...
 
int valid ()
 Returns a bitfield specifyign which of the parts have been assembled since construction or flush() according to the format (VALUE|RHS|MATRX). More...
 
void flush (int flags=(VALUE|RHS|MATRIX))
 Destructs parts of the assembled quatities according to the format (VALUE|RHS|MATRX) More...
 
std::pair< size_t, size_t > size (int row, int col) const
 the size of a matrix block (in terms of scalar rows/columns) More...
 
template<class DataOutIter >
void toSequence (int rbegin, int rend, DataOutIter xi) const
 Writes the subvector defined by the half-open block range [rbegin,rend) to the output iterator xi. More...
 
MatrixBlockArraygetMatrix () const
 Returns a mutable reference to the sequence of matrix blocks. More...
 
template<int first = 0, int last = TestVariableSetDescription::noOfVariables>
TestVariableSetDescription::template CoefficientVectorRepresentation< first, last >::type rhs () const
 Returns a contiguous subrange of the rhs coefficient vectors. More...
 
Tuning parameters
SelfsetNSimultaneousBlocks (int n)
 Defines how many cells are assembled locally before scattering them together into the global data structures. More...
 
SelfsetRowBlockFactor (double a)
 Defines how many more row blocks in each matrix are used compared to the number of threads. More...
 
SelfsetLocalStorageSize (size_t s)
 Defines how many memory the locally assembled matrices may occupy before they are scattered. More...
 
SelfsetGatherTimings (bool gt)
 Whether to gather timing statistics using Kaskade::Timings. More...
 
Submatrix access
template<class Matrix >
Matrix get (bool extractOnlyLowerTriangle=false, int rbegin=0, int rend=TestVariableSetDescription::noOfVariables, int cbegin=0, int cend=AnsatzVariableSetDescription::noOfVariables) const
 Extracts the submatrix defined by the half-open block ranges [rbegin,rend), [cbegin,cend). More...
 
template<int row, int col>
auto const & get () const
 Extracts a raw single submatrix block indexed by row and column. More...
 
template<class MatrixType , class BlockInformation >
MatrixType get (bool extractOnlyLowerTriangle) const
 Extracts the submatrix defined by the half-open block ranges given as template parameter. More...
 
General information
size_t nnz (size_t rbegin=0, size_t rend=TestVariableSetDescription::noOfVariables, size_t cbegin=0, size_t cend=AnsatzVariableSetDescription::noOfVariables, bool extractOnlyLowerTriangle=false) const
 Returns the number of structurally nonzero entries in the submatrix defined by the half-open block ranges [rbegin,rend), [cbegin,cend). More...
 
size_t nrows (int firstBlock=0, int lastBlock=TestVariableSetDescription::noOfVariables) const
 Returns the number of scalar rows in the row block range [firstBlock,lastBlock[. More...
 
size_t ncols (int firstBlock=0, int lastBlock=AnsatzVariableSetDescription::noOfVariables) const
 Returns the number of scalar cols in the column block range [firstBlock,lastBlock[. More...
 
Spaces const & spaces () const
 Returns the list of spaces used. More...
 
GridView const & gridView () const
 The grid view used. More...
 

Static Public Attributes

static unsigned int const VALUE = Assembler::VALUE
 DEPRECATED, enums in the Assembler namespace. More...
 
static unsigned int const RHS = Assembler::RHS
 
static unsigned int const MATRIX = Assembler::MATRIX
 
static unsigned int const EVERYTHING = Assembler::EVERYTHING
 

Protected Types

typedef AssemblyDetail::FormulationPolicy< F > Policy
 
typedef AnsatzVariableSetDescription::Variables AnsatzVariables
 
typedef TestVariableSetDescription::Variables TestVariables
 
typedef GridView::template Codim< 0 >::Iterator CellIterator
 
typedef CellIterator::Entity Entity
 

Protected Member Functions

std::pair< RhsArray &, RhsBlockArray & > getRhs () const
 
void reactToRefinement (GridSignals::Status const status)
 

Protected Attributes

boost::signals2::scoped_connection refConnection
 
Spaces const & spaces_
 
GridManagerBase< Grid > const & gridManager
 
GridView const & gv
 
IndexSet const & indexSet
 
std::shared_ptr< MatrixBlockArraymatrixBlocks
 
std::shared_ptr< RhsArrayrhss
 
std::shared_ptr< RhsBlockArrayrhsBlocks
 
BoundaryDetector boundaryDetector
 
int validparts
 
int nSimultaneousBlocks
 
double rowBlockFactor
 
size_t localStorageSize
 
bool gatherTimings = false
 

Static Protected Attributes

static constexpr bool innerBoundaries = hasInnerBoundaryCache<Functional>
 

Member Typedef Documentation

◆ AnsatzVariables

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AnsatzVariableSetDescription::Variables Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::AnsatzVariables
protected

Definition at line 1468 of file assemblerCore.hh.

◆ AnsatzVariableSetDescription

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef Functional::AnsatzVars Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::AnsatzVariableSetDescription

ansatz variables

Definition at line 1442 of file assemblerCore.hh.

◆ CellIterator

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef GridView::template Codim<0>::Iterator Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::CellIterator
protected

Definition at line 1607 of file assemblerCore.hh.

◆ Entity

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef CellIterator::Entity Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Entity
protected

Definition at line 1608 of file assemblerCore.hh.

◆ field_type

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef Functional::Scalar Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::field_type

Underlying field type.

Definition at line 1463 of file assemblerCore.hh.

◆ Functional

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef F Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Functional

functional

Definition at line 1440 of file assemblerCore.hh.

◆ Grid

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AnsatzVariableSetDescription::Grid Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Grid

grid

Definition at line 1448 of file assemblerCore.hh.

◆ GridView

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AnsatzVariableSetDescription::GridView Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::GridView

The grid view on which the variables live.

Definition at line 1453 of file assemblerCore.hh.

◆ IndexSet

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AnsatzVariableSetDescription::IndexSet Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::IndexSet

index set

Definition at line 1458 of file assemblerCore.hh.

◆ MatrixBlockArray

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
using Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::MatrixBlockArray = typename AssemblyDetail::BlockArray<Policy,AnsatzVariables,TestVariables,SparseIndex>::type

A boost::fusion sequence of AssemblyDetail::MatrixBlock elements for present matrix blocks.

The elements of this sequence type include both block info and the global matrices.

Definition at line 1486 of file assemblerCore.hh.

◆ Policy

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AssemblyDetail::FormulationPolicy<F> Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Policy
protected

Definition at line 1434 of file assemblerCore.hh.

◆ RhsArray

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef TestVariableSetDescription::template CoefficientVectorRepresentation ::type Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::RhsArray

A LinearProductSpace type of right hand side coefficient vectors.

Definition at line 1498 of file assemblerCore.hh.

◆ RhsBlockArray

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AssemblyDetail::RhsBlockArray<Policy,TestVariables>::type Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::RhsBlockArray

A boost::fusion sequence of AssemblyDetail::RhsBock elements for present rhs blocks.

The elements of this sequence type include the rhs block infos.

Definition at line 1493 of file assemblerCore.hh.

◆ Scalar

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
using Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Scalar = field_type

Definition at line 1465 of file assemblerCore.hh.

◆ Self

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef VariationalFunctionalAssembler<F,SparseIndex,BoundaryDetector,QuadRule> Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Self

Definition at line 1437 of file assemblerCore.hh.

◆ Spaces

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef AnsatzVariableSetDescription::Spaces Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::Spaces

spaces

Definition at line 1456 of file assemblerCore.hh.

◆ TestVariables

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef TestVariableSetDescription::Variables Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::TestVariables
protected

Definition at line 1469 of file assemblerCore.hh.

◆ TestVariableSetDescription

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
typedef Functional::TestVars Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::TestVariableSetDescription

test variables

Definition at line 1445 of file assemblerCore.hh.

Constructor & Destructor Documentation

◆ VariationalFunctionalAssembler()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::VariationalFunctionalAssembler ( Spaces const &  spaces)
inlineexplicit

Construct an empty assembler, gridManager is implicitly obtained from the first space.

Definition at line 1521 of file assemblerCore.hh.

Member Function Documentation

◆ assemble() [1/2]

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
template<class BlockFilter , class CellFilter >
void Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::assemble ( F const &  f,
CellFilter const &  cellFilter,
unsigned int  flags = Assembler::EVERYTHING,
int  nTasks = 0,
bool  verbose = false 
)
inline

Create data in assembler.

Assembles the value and/or derivative and/or Hessian (according to the given flags, default is to assemble all three) of the FE-discretized variational functional given by f. The block filter can be used for partial assembly of specified blocks in the matrix and rhs.

The data that is not concerned by the given flags and block filter is left unmodified. Be careful: This allows to drive the Galerkin operator representation into a semantically inconsistent state. This need not be a problem (instead it can be useful and more efficient in certain situations), but one has to be aware of that fact. The data that is not concerned by the cell filter is cleared before assembly (must be, because new values are partially added up to existing data).

Template Parameters
BlockFiltera traits class implementing the BlockFilter concept
CellFiltera callable type for deciding on which cells to assemble

The assembled data can afterwards be accessed via the functional, rhs, and matrix methods.

Parameters
fthe variational functional or weak formulation to be assembled
cellFiltera callable mapping a grid cell to a boolean. Assembly is performed only on those cells for which the cell filter returns true. This can be used, e.g., for subassembly in non-overlapping domain decomposition such as BDDC, or for partial update of stiffness matrices/rhs in case of spatially local changes in Newton's method.
flagsa bitset determining what to assemble
nTasksnumber of tasks to submit to the thread pool for parallel execution. 0 (default) means to use the default hardware concurrency reported by the system. Sequential computation is performed if the grid is reported not to be thread-safe. You can call GridManagerBase::enforceConcurrentReads to lie about thread-safety in case you know (or believe) the implementation to be thread safe.
verbosewhether to output status messages to standard output
See also
NumaThreadPool

Definition at line 1652 of file assemblerCore.hh.

◆ assemble() [2/2]

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
void Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::assemble ( F const &  f,
unsigned int  flags = Assembler::EVERYTHING,
int  nThreads = 0,
bool  verbose = false 
)
inline

◆ flush()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
void Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::flush ( int  flags = (VALUE|RHS|MATRIX))
inline

Destructs parts of the assembled quatities according to the format (VALUE|RHS|MATRX)

The default flag leads to destruction of all three values.

Definition at line 2109 of file assemblerCore.hh.

Referenced by Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::reactToRefinement().

◆ get() [1/3]

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
template<int row, int col>
auto const & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get ( ) const
inline

Extracts a raw single submatrix block indexed by row and column.

Note that for symmetric matrices, only the lower triangular part is stored.

Compilation will fail if the referenced matrix block is not present.

Definition at line 2188 of file assemblerCore.hh.

◆ get() [2/3]

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
template<class MatrixType , class BlockInformation >
MatrixType Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get ( bool  extractOnlyLowerTriangle) const
inline

Extracts the submatrix defined by the half-open block ranges given as template parameter.

Provide the block information via IstlInterface_Detail::BlockInfo.

Definition at line 2199 of file assemblerCore.hh.

◆ get() [3/3]

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
template<class Matrix >
Matrix Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get ( bool  extractOnlyLowerTriangle = false,
int  rbegin = 0,
int  rend = TestVariableSetDescription::noOfVariables,
int  cbegin = 0,
int  cend = AnsatzVariableSetDescription::noOfVariables 
) const
inline

Extracts the submatrix defined by the half-open block ranges [rbegin,rend), [cbegin,cend).

If extractOnlyLowerTriangle is true, the Galerkin operator is treated as if its upper triangle was structurally zero.

In order to be able to read the assembled matrix into the smart pointer a specialization of AssemblyDetail::Fill, i.e. AssemblyDetail::Fill<MatrixType> must be available!!!

Kaskade provides spezializations of AssemblyDetail::Fill for

Parameters
extractOnlyLowerTriangleif true, and the matrix is symmetric, only the lower triangular part is extracted
rbeginnumber first block row to store
rendnumber of one past the last block row to store
cbeginnumber first block column to store
cendnumber of one past the last block column to store

Definition at line 2160 of file assemblerCore.hh.

◆ getMatrix()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
MatrixBlockArray & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::getMatrix ( ) const
inline

Returns a mutable reference to the sequence of matrix blocks.

The data structure for the matrix is created on the fly if needed.

Todo:
check const correctness - why is this method const and returns a mutable matrix?

Definition at line 2285 of file assemblerCore.hh.

Referenced by Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::assemble(), Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get(), and Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::nnz().

◆ getRhs()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
std::pair< RhsArray &, RhsBlockArray & > Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::getRhs ( ) const
inlineprotected

◆ gridView()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
GridView const & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::gridView ( ) const
inline

The grid view used.

Definition at line 2269 of file assemblerCore.hh.

◆ ncols()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
size_t Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::ncols ( int  firstBlock = 0,
int  lastBlock = AnsatzVariableSetDescription::noOfVariables 
) const
inline

Returns the number of scalar cols in the column block range [firstBlock,lastBlock[.

Definition at line 2253 of file assemblerCore.hh.

Referenced by Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get().

◆ nnz()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
size_t Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::nnz ( size_t  rbegin = 0,
size_t  rend = TestVariableSetDescription::noOfVariables,
size_t  cbegin = 0,
size_t  cend = AnsatzVariableSetDescription::noOfVariables,
bool  extractOnlyLowerTriangle = false 
) const
inline

Returns the number of structurally nonzero entries in the submatrix defined by the half-open block ranges [rbegin,rend), [cbegin,cend).

This is the number of elements written by the corresponding call of toTriplet().

Definition at line 2234 of file assemblerCore.hh.

Referenced by Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get(), Kaskade::HigherOrderRecovery< Grid, Space >::getErrorFunction(), and Kaskade::Limex< Eq >::step().

◆ nrows()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
size_t Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::nrows ( int  firstBlock = 0,
int  lastBlock = TestVariableSetDescription::noOfVariables 
) const
inline

Returns the number of scalar rows in the row block range [firstBlock,lastBlock[.

Definition at line 2245 of file assemblerCore.hh.

Referenced by Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::get().

◆ reactToRefinement()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
void Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::reactToRefinement ( GridSignals::Status const  status)
inlineprotected

Definition at line 2362 of file assemblerCore.hh.

◆ rhs()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
template<int first = 0, int last = TestVariableSetDescription::noOfVariables>
TestVariableSetDescription::template CoefficientVectorRepresentation< first, last >::type Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::rhs ( ) const
inline

◆ setGatherTimings()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Self & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::setGatherTimings ( bool  gt)
inline

Whether to gather timing statistics using Kaskade::Timings.

Gathering timings can be switched on and off. Note that if several assembler calls are done in parallel, it must be switched off, since the timings functionality is not thread-safe. Defaults to off.

Definition at line 1587 of file assemblerCore.hh.

◆ setLocalStorageSize()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Self & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::setLocalStorageSize ( size_t  s)
inline

Defines how many memory the locally assembled matrices may occupy before they are scattered.

The amount of memory (in bytes) should scale with CPU cache, roughly a quarter of the second level cache size should be a reasonable value. Note that this is merely a hint, not a hard limit. The actual local matrix storage can be slightly larger.

Definition at line 1574 of file assemblerCore.hh.

◆ setNSimultaneousBlocks()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Self & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::setNSimultaneousBlocks ( int  n)
inline

Defines how many cells are assembled locally before scattering them together into the global data structures.

Higher numbers reduce the overhead of locking for mutually exclusive write access to global matrix/rhs and improve memory access locality, but increase the amount of thread-local computation and storage. The default value is 42, but the performance appears to be not very sensitive to this parameter in the range [20,60].

Definition at line 1545 of file assemblerCore.hh.

◆ setRowBlockFactor()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Self & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::setRowBlockFactor ( double  a)
inline

Defines how many more row blocks in each matrix are used compared to the number of threads.

When scattering the local matrices/rhs into the global data structures, blocks of sequential rows are processed in turn. Each block of rows is equipped with a lock in order to guarantee exclusive write access.

A high number of row blocks allows fine granularity and improves cache locality and simultaneous scattering into different parts of the global data structures, but incurs some inefficiency due to more frequent locking and more thread-local processing. The default value is to use two times the number of threads.

Definition at line 1561 of file assemblerCore.hh.

◆ size()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
std::pair< size_t, size_t > Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::size ( int  row,
int  col 
) const
inline

the size of a matrix block (in terms of scalar rows/columns)

Returns
a pair (nrows,ncols)

Definition at line 2125 of file assemblerCore.hh.

◆ spaces()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Spaces const & Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::spaces ( ) const
inline

◆ toSequence()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
template<class DataOutIter >
void Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::toSequence ( int  rbegin,
int  rend,
DataOutIter  xi 
) const
inline

Writes the subvector defined by the half-open block range [rbegin,rend) to the output iterator xi.

Definition at line 2213 of file assemblerCore.hh.

Referenced by Kaskade::HigherOrderRecovery< Grid, Space >::getErrorFunction(), and Kaskade::Limex< Eq >::step().

◆ valid()

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
int Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::valid ( )
inline

Returns a bitfield specifyign which of the parts have been assembled since construction or flush() according to the format (VALUE|RHS|MATRX).

Definition at line 2102 of file assemblerCore.hh.

Member Data Documentation

◆ boundaryDetector

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
BoundaryDetector Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::boundaryDetector
protected

Definition at line 2380 of file assemblerCore.hh.

◆ EVERYTHING

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
unsigned int const Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::EVERYTHING = Assembler::EVERYTHING
static

Definition at line 1531 of file assemblerCore.hh.

◆ gatherTimings

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
bool Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::gatherTimings = false
protected

◆ gridManager

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
GridManagerBase<Grid> const& Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::gridManager
protected

◆ gv

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
GridView const& Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::gv
protected

◆ indexSet

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
IndexSet const& Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::indexSet
protected

Definition at line 2375 of file assemblerCore.hh.

◆ innerBoundaries

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
constexpr bool Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::innerBoundaries = hasInnerBoundaryCache<Functional>
staticconstexprprotected

Definition at line 1477 of file assemblerCore.hh.

◆ localStorageSize

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
size_t Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::localStorageSize
protected

◆ MATRIX

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
unsigned int const Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::MATRIX = Assembler::MATRIX
static

◆ matrixBlocks

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
std::shared_ptr<MatrixBlockArray> Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::matrixBlocks
mutableprotected

◆ nSimultaneousBlocks

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
int Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::nSimultaneousBlocks
protected

◆ refConnection

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
boost::signals2::scoped_connection Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::refConnection
protected

Definition at line 2370 of file assemblerCore.hh.

◆ RHS

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
unsigned int const Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::RHS = Assembler::RHS
static

◆ rhsBlocks

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
std::shared_ptr<RhsBlockArray> Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::rhsBlocks
mutableprotected

◆ rhss

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
std::shared_ptr<RhsArray> Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::rhss
mutableprotected

◆ rowBlockFactor

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
double Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::rowBlockFactor
protected

◆ spaces_

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
Spaces const& Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::spaces_
protected

◆ validparts

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
int Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::validparts
protected

◆ VALUE

template<class F , class SparseIndex = size_t, class BoundaryDetector = CachingBoundaryDetector<typename F::AnsatzVars::GridView>, class QuadRule = Dune::QuadratureRule<typename F::AnsatzVars::Grid::ctype, F::AnsatzVars::Grid::dimension>>
unsigned int const Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::VALUE = Assembler::VALUE
static

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