1 #ifndef DUNE_FEM_OPERATOR_LINEAR_BLOCKDIAGONAL_HH 2 #define DUNE_FEM_OPERATOR_LINEAR_BLOCKDIAGONAL_HH 6 #include <dune/common/fmatrix.hh> 7 #include <dune/common/nullptr.hh> 26 class LocalBlock = Dune::FieldMatrix<
typename DiscreteFunctionSpace ::
27 RangeFieldType, DiscreteFunctionSpace::localBlockSize, DiscreteFunctionSpace::localBlockSize > >
42 typedef typename RangeFunctionType::DiscreteFunctionSpaceType
RangeSpaceType;
54 typedef typename LocalBlockType::row_type
DofType ;
57 template<
class Operation >
60 typedef typename DiscreteFunctionSpaceType
86 if( &domainSpace != &rangeSpace )
87 DUNE_THROW( InvalidStateException,
"BlockDiagonalLinearOperator must be created with identical spaces." );
90 void operator() (
const DomainFunctionType &u, RangeFunctionType &w )
const 95 template <
class DomainSpace,
class RangeSpace >
102 template <
class DomainSpace,
class RangeSpace >
106 typedef typename std::vector< LocalBlockType >::const_iterator Iterator;
107 const typename DomainSpace :: RangeFieldType *uit = u.
leakPointer();
108 typename RangeSpace :: RangeFieldType *wit = w.
leakPointer();
110 for( Iterator dit =
diagonal_.begin(); dit != dend; ++dit )
122 typedef typename std::vector< LocalBlockType >::iterator Iterator;
124 for( Iterator dit =
diagonal_.begin(); dit != dend; ++dit )
128 template<
class Functor >
131 typedef typename std::vector< LocalBlockType >::iterator Iterator;
133 for( Iterator dit =
diagonal_.begin(); dit != dend; ++dit )
139 typedef typename std::vector< LocalBlockType >::iterator Iterator;
141 for( Iterator dit =
diagonal_.begin(); dit != dend; ++dit )
148 typedef typename std::vector< LocalBlockType >::iterator Iterator;
149 typedef typename std::vector< LocalBlockType >::const_iterator ConstIterator;
150 ConstIterator otherIt = other.
diagonal_.begin();
152 for( Iterator dit =
diagonal_.begin(); dit != dend; ++dit, ++otherIt )
154 (*dit).rightmultiply( *otherIt );
161 typedef typename std::vector< LocalBlockType >::iterator Iterator;
162 typedef typename std::vector< LocalBlockType >::const_iterator ConstIterator;
163 ConstIterator otherIt = other.
diagonal_.begin();
165 for( Iterator dit =
diagonal_.begin(); dit != dend; ++dit, ++otherIt )
167 (*dit).leftmultiply( *otherIt );
197 template<
class Operation >
201 return space().createDataHandle( *
this, operation );
204 template<
class Stencil >
211 LocalColumnObjectType
localColumn (
const DomainEntityType &domainEntity )
const 216 LocalMatrixType
localMatrix (
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity )
const;
240 template<
class DiscreteFunctionSpace,
class LocalBlock >
261 template<
class DiscreteFunctionSpace,
class LocalBlock >
280 typedef DomainBasisFunctionSetType BasisFunctionSetType;
283 struct SetLocalBlockFunctor
285 SetLocalBlockFunctor ( std::vector< LocalBlockType > &diagonal,
286 LocalBlockType *&localBlock )
288 localBlock_( localBlock )
291 void operator() (
int localDoF, std::size_t globalDoF )
293 assert( localDoF == 0 );
299 std::vector< LocalBlockType > &
diagonal_;
300 LocalBlockType *&localBlock_;
306 localBlock_( nullptr )
309 void init (
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity )
311 basisFunctionSet_ =
domainSpace().basisFunctionSet( domainEntity );
312 SetLocalBlockFunctor f( op_->diagonal_, localBlock_ );
313 domainSpace().blockMapper().mapEach( domainEntity, f );
314 if( &domainEntity != &rangeEntity )
316 static LocalBlockType dummyBlock( 0 );
318 LocalBlockType *otherBlock = 0;
319 SetLocalBlockFunctor f( op_->diagonal_, otherBlock );
320 rangeSpace().blockMapper().mapEach( rangeEntity, f );
323 if( otherBlock != localBlock_ )
324 localBlock_ = &dummyBlock ;
329 void scale (
const RangeFieldType &a ) { localBlock() *= a; }
331 RangeFieldType
get (
int i,
int j )
const {
return localBlock()[ i ][ j ]; }
332 void add (
int i,
int j,
const RangeFieldType &value ) { localBlock()[ i ][ j ] += value; }
333 void set (
int i,
int j,
const RangeFieldType &value ) { localBlock()[ i ][ j ] = value; }
339 for(
int i = 0; i < rows(); ++i )
343 template<
class DomainLocalFunction,
class RangeLocalFunction >
344 void multiplyAdd (
const DomainLocalFunction &x, RangeLocalFunction &y )
const 346 localBlock().umv( x, y );
352 int rows ()
const {
return localBlock().N(); }
353 int columns ()
const {
return localBlock().M(); }
361 const DomainEntityType &
domainEntity ()
const {
return domainBasisFunctionSet().entity(); }
362 const RangeEntityType &
rangeEntity ()
const {
return rangeBasisFunctionSet().entity(); }
365 const LocalBlockType &localBlock ()
const { assert( localBlock_ );
return *localBlock_; }
366 LocalBlockType &localBlock () { assert( localBlock_ );
return *localBlock_; }
369 BasisFunctionSetType basisFunctionSet_;
370 LocalBlockType *localBlock_;
378 template<
class DiscreteFunctionSpace,
class LocalBlock >
388 ObjectType *
newObject ()
const {
return new ObjectType( *op_ ); }
399 template<
class DiscreteFunctionSpace,
class LocalBlock >
411 #endif // #ifndef DUNE_FEM_OPERATOR_LINEAR_BLOCKDIAGONAL_HH RangeFunctionType::DiscreteFunctionSpaceType RangeSpaceType
Definition: blockdiagonal.hh:42
Traits::DomainSpaceType DomainSpaceType
type of domain discrete function space
Definition: localmatrix.hh:51
AdaptiveDiscreteFunction< DiscreteFunctionSpace > RangeFunctionType
type of discrete function in the operator's range
Definition: operator.hh:30
Definition: blockdiagonal.hh:58
void clear()
Definition: blockdiagonal.hh:120
void clear()
Definition: blockdiagonal.hh:328
DomainSpaceType::EntityType DomainEntityType
Definition: blockdiagonal.hh:44
Definition: blockdiagonal.hh:28
Definition: localmatrixwrapper.hh:16
Definition: blockdiagonal.hh:241
DiscreteFunctionSpaceType::template CommDataHandle< ThisType, Operation >::Type Type
Definition: blockdiagonal.hh:62
default implementation for a general operator stencil
Definition: stencil.hh:31
void add(int i, int j, const RangeFieldType &value)
Definition: blockdiagonal.hh:332
LocalColumnObjectType localColumn(const DomainEntityType &domainEntity) const
Definition: blockdiagonal.hh:211
void communicate()
copy matrices to ghost cells to make this class work in parallel
Definition: blockdiagonal.hh:189
BaseType::RangeFieldType RangeFieldType
Definition: blockdiagonal.hh:271
const std::string & name() const
Definition: blockdiagonal.hh:224
Definition: adaptivefunction/adaptivefunction.hh:34
BaseType::RangeFunctionType RangeFunctionType
Definition: blockdiagonal.hh:36
const DomainSpaceType & domainSpace() const
Definition: blockdiagonal.hh:355
RangeSpaceType::EntityType RangeEntityType
Definition: blockdiagonal.hh:45
void clearRow(int i)
Definition: blockdiagonal.hh:335
void forEach(const Functor &functor)
Definition: blockdiagonal.hh:129
const RangeBasisFunctionSetType & rangeBasisFunctionSet() const
Definition: blockdiagonal.hh:359
LocalBlockType::row_type DofType
Definition: blockdiagonal.hh:54
static const int localBlockSize
Definition: blockdiagonal.hh:47
DomainSpaceType::BasisFunctionSetType DomainBasisFunctionSetType
type of base function sets within domain function space
Definition: localmatrix.hh:58
const RangeSpaceType & space_
Definition: blockdiagonal.hh:228
void leftmultiply(const ThisType &other)
Definition: blockdiagonal.hh:158
const DomainEntityType & domainEntity() const
Definition: blockdiagonal.hh:361
LocalMatrixStackType localMatrixStack_
Definition: blockdiagonal.hh:232
Interface for local matrix classes.
Definition: localmatrix.hh:28
OperatorType::RangeSpaceType RangeSpaceType
Definition: blockdiagonal.hh:251
const DomainSpaceType & domainSpace() const
Definition: blockdiagonal.hh:218
LocalMatrixFactory localMatrixFactory_
Definition: blockdiagonal.hh:231
LocalMatrixFactory(OperatorType &op)
Definition: blockdiagonal.hh:384
BlockDiagonalLinearOperator< DiscreteFunctionSpace, LocalBlock > OperatorType
Definition: blockdiagonal.hh:269
SizeType size() const
Return the number of blocks in the block vector.
Definition: discretefunction.hh:740
void rightmultiply(const ThisType &other)
Definition: blockdiagonal.hh:145
BaseType::DomainFieldType DomainFieldType
Definition: blockdiagonal.hh:38
Definition: coordinate.hh:4
Definition: blockdiagonal.hh:262
ColumnObject< ThisType > LocalColumnObjectType
Definition: blockdiagonal.hh:75
BlockDiagonalLinearOperator< DiscreteFunctionSpace, LocalBlock > OperatorType
Definition: blockdiagonal.hh:381
LocalBlock LocalBlockType
Definition: blockdiagonal.hh:49
BaseType::DomainBasisFunctionSetType DomainBasisFunctionSetType
Definition: blockdiagonal.hh:273
void clearCol(int j)
Definition: blockdiagonal.hh:337
DomainSpaceType DiscreteFunctionSpaceType
Definition: blockdiagonal.hh:55
BaseType::DomainEntityType DomainEntityType
Definition: blockdiagonal.hh:276
std::string name_
Definition: blockdiagonal.hh:227
BaseType::RangeFieldType RangeFieldType
Definition: blockdiagonal.hh:39
CommDataHandle< Operation >::Type dataHandle(const Operation *operation)
return reference to data handle object (needed to make this class work with CommunicationManager) ...
Definition: blockdiagonal.hh:199
void invert()
Definition: blockdiagonal.hh:137
DofType * leakPointer()
Definition: adaptivefunction/adaptivefunction.hh:95
void reserve(const Stencil &stencil, bool verbose=false)
Definition: blockdiagonal.hh:205
AdaptiveDiscreteFunction< DiscreteFunctionSpace >::RangeFieldType DomainFieldType
field type of the operator's domain
Definition: operator.hh:33
Traits::RangeSpaceType RangeSpaceType
type of range discrete function space
Definition: localmatrix.hh:54
AdaptiveDiscreteFunction< DiscreteFunctionSpace >::RangeFieldType RangeFieldType
field type of the operator's range
Definition: operator.hh:35
OperatorType::DomainSpaceType DomainSpaceType
Definition: blockdiagonal.hh:250
RangeFieldType LittleBlockType
Definition: blockdiagonal.hh:253
BlockDiagonalLinearOperator(const std::string &name, const DomainSpaceType &domainSpace, const RangeSpaceType &rangeSpace)
Definition: blockdiagonal.hh:78
LocalMatrix ObjectType
Definition: blockdiagonal.hh:382
ObjectStack< LocalMatrixFactory > LocalMatrixStackType
Definition: blockdiagonal.hh:69
LocalBlockType * DofBlockPtrType
Definition: blockdiagonal.hh:52
ConstDofBlockPtrType block(const size_t block) const
return block matrix for given block number (== entity number)
Definition: blockdiagonal.hh:179
OperatorType::RangeFieldType RangeFieldType
Definition: blockdiagonal.hh:248
const RangeEntityType & rangeEntity() const
Definition: blockdiagonal.hh:362
void resort()
Definition: blockdiagonal.hh:350
ObjectType * newObject() const
Definition: blockdiagonal.hh:388
OperatorType::LocalMatrix LocalMatrixType
Definition: blockdiagonal.hh:246
LocalMatrixWrapper< LocalMatrixStackType > LocalMatrixType
Definition: blockdiagonal.hh:73
void operator()(const DomainFunctionType &u, RangeFunctionType &w) const
Definition: blockdiagonal.hh:90
abstract matrix operator
Definition: operator.hh:106
DomainFunctionType::DiscreteFunctionSpaceType DomainSpaceType
Definition: blockdiagonal.hh:41
DofBlockPtrType block(const size_t block)
return block matrix for given block number (== entity number)
Definition: blockdiagonal.hh:172
void scale(const RangeFieldType &a)
Definition: blockdiagonal.hh:329
BaseType::RangeBasisFunctionSetType RangeBasisFunctionSetType
Definition: blockdiagonal.hh:274
int rows() const
Definition: blockdiagonal.hh:352
void init(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity)
Definition: blockdiagonal.hh:309
Definition: columnobject.hh:11
void multiply(const AdaptiveDiscreteFunction< DomainSpace > &u, AdaptiveDiscreteFunction< RangeSpace > &w) const
Definition: blockdiagonal.hh:103
BaseType::RangeEntityType RangeEntityType
Definition: blockdiagonal.hh:277
const RangeSpaceType & rangeSpace() const
Definition: blockdiagonal.hh:219
LocalMatrix(OperatorType &op)
Definition: blockdiagonal.hh:304
void multiplyAdd(const DomainLocalFunction &x, RangeLocalFunction &y) const
Definition: blockdiagonal.hh:344
const DomainBasisFunctionSetType & domainBasisFunctionSet() const
Definition: blockdiagonal.hh:358
const RangeSpaceType & rangeSpace() const
Definition: blockdiagonal.hh:356
BaseType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
Definition: adaptivefunction/adaptivefunction.hh:53
const DomainSpaceType & space() const
return reference to space (needed to make this class work with CommunicationManager) ...
Definition: blockdiagonal.hh:222
LocalMatrixType localMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity) const
Definition: blockdiagonal.hh:402
BaseType::DomainFunctionType DomainFunctionType
Definition: blockdiagonal.hh:35
int columns() const
Definition: blockdiagonal.hh:353
Definition: blockdiagonal.hh:379
Definition: bartonnackmaninterface.hh:15
void finalize()
Definition: blockdiagonal.hh:349
std::vector< LocalBlockType > diagonal_
Definition: blockdiagonal.hh:229
const LocalBlockType * ConstDofBlockPtrType
Definition: blockdiagonal.hh:53