dune-fem 2.8-git
discretefunction_inline.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_DISCRETEFUNCTION_INLINE_HH
2#define DUNE_FEM_DISCRETEFUNCTION_INLINE_HH
3
4#include <fstream>
5
6#include <dune/geometry/referenceelements.hh>
7
13
14#include "discretefunction.hh"
15
16namespace Dune
17{
18
19 namespace Fem
20 {
21
22 // DiscreteFunctionDefault
23 // -----------------------
24
25 template< class Impl >
27 :: DiscreteFunctionDefault ( const std::string &name,
28 const DiscreteFunctionSpaceType &dfSpace )
29 : dfSpace_( referenceToSharedPtr( dfSpace ) ),
30 ldvStack_( std::max( std::max( sizeof( DofType ), sizeof( DofType* ) ),
31 sizeof(typename LocalDofVectorType::value_type) ) // for PetscDiscreteFunction
32 * space().blockMapper().maxNumDofs() * DiscreteFunctionSpaceType::localBlockSize ),
33 ldvAllocator_( &ldvStack_ ),
34 localFunction_( space() ),
35 name_( name ),
36 scalarProduct_( dfSpace )
37 {
38 }
39
40
41 template< class Impl >
42 inline DiscreteFunctionDefault< Impl >::DiscreteFunctionDefault ( std::string name, std::shared_ptr< const DiscreteFunctionSpaceType > dfSpace )
43 : dfSpace_( std::move( dfSpace ) ),
44 ldvStack_( std::max( std::max( sizeof( DofType ), sizeof( DofType* ) ),
45 sizeof(typename LocalDofVectorType::value_type) ) // for PetscDiscreteFunction
46 * space().blockMapper().maxNumDofs() * DiscreteFunctionSpaceType::localBlockSize ),
47 ldvAllocator_( &ldvStack_ ),
48 localFunction_( space() ),
49 name_( std::move( name ) ),
50 scalarProduct_( dfSpace )
51 {}
52
53
54 template< class Impl >
56 : BaseType( static_cast< const BaseType & >( other ) ),
57 dfSpace_( other.dfSpace_ ),
58 ldvStack_( std::max( std::max( sizeof( DofType ), sizeof( DofType* ) ),
59 sizeof(typename LocalDofVectorType::value_type) ) // for PetscDiscreteFunction
60 * space().blockMapper().maxNumDofs() * DiscreteFunctionSpaceType::localBlockSize ),
61 ldvAllocator_( &ldvStack_ ),
62 localFunction_( space() ),
63 name_( other.name_ ),
64 scalarProduct_( other.scalarProduct_ )
65 {
66 if( other.assembleOperation_ != std::type_index( typeid( void ) ) )
67 DUNE_THROW( InvalidStateException, "Cannot copy discrete function during assembly" );
68 assert( other.assembleCount_ == 0 );
69 }
70
71
72 template< class Impl >
74 :: DiscreteFunctionDefault ( DiscreteFunctionDefault && other )
75 : BaseType( static_cast< BaseType&& >( other ) ),
76 dfSpace_( std::move( other.dfSpace_ ) ),
77 ldvStack_( std::move( other.ldvStack_ ) ),
78 ldvAllocator_( &ldvStack_ ),
79 localFunction_( space() ),
80 name_( std::move( other.name_ ) ),
81 scalarProduct_( std::move( other.scalarProduct_ ) )
82 {
83 if( other.assembleOperation_ != std::type_index( typeid( void ) ) )
84 DUNE_THROW( InvalidStateException, "Cannot move discrete function during assembly" );
85 assert( other.assembleCount_ == 0 );
86 }
87
88
89 template< class Impl >
91 :: print ( std::ostream &out ) const
92 {
93 const auto end = BaseType :: dend();
94 for( auto dit = BaseType :: dbegin(); dit != end; ++dit )
95 out << (*dit) << std::endl;
96 }
97
98
99 template< class Impl >
102 {
103 const auto end = BaseType :: dend();
104 for( auto it = BaseType :: dbegin(); it != end; ++it )
105 {
106 if( ! std::isfinite( *it ) )
107 return false ;
108 }
109
110 return true;
111 }
112
113
114 template< class Impl >
115 template< class DFType >
118 {
119 if( BaseType::size() != g.size() )
120 DUNE_THROW(InvalidStateException,"DiscreteFunctionDefault: sizes do not match in axpy");
121
122 // apply axpy to all dofs from g
123 const auto end = BaseType::dend();
124 auto git = g.dbegin();
125 for( auto it = BaseType::dbegin(); it != end; ++it, ++git )
126 *it += s * (*git );
127 }
128
129
130 template< class Impl >
131 template< class DFType >
134 {
135 if( BaseType::size() != g.size() )
136 {
137 std::cout << BaseType::size() << " vs " << g.size() << std::endl;
138 DUNE_THROW(InvalidStateException,"DiscreteFunctionDefault: sizes do not match in assign");
139 }
140
141 // copy all dofs from g to this
142 const auto end = BaseType::dend();
143 auto git = g.dbegin();
144 for( auto it = BaseType::dbegin(); it != end; ++it, ++git )
145 *it = *git;
146 }
147
148
149 template< class Impl >
150 template< class Operation >
151 inline typename DiscreteFunctionDefault< Impl >
152 :: template CommDataHandle< Operation > :: Type
154 {
155 return BaseType :: space().createDataHandle( asImp(), operation );
156 }
157
158
159 template< class Impl >
160 template< class Functor >
162 ::evaluateGlobal ( const DomainType &x, Functor functor ) const
163 {
164 typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
165 EntitySearch< GridPartType, EntityType::codimension > entitySearch( BaseType::space().gridPart() );
166
167 const auto& entity(entitySearch( x ));
168 const auto geometry = entity.geometry();
169
170 // bind local function to entity
171 auto guard = bindGuard( localFunction_, entity );
172 // obtain dofs (since it's a temporary local function)
173 getLocalDofs( entity, localFunction_.localDofVector());
174 // evaluate functor
175 functor( geometry.local( x ), localFunction_ );
176 }
177
178
179 template< class Impl >
180 template< class DFType >
181 inline typename DiscreteFunctionDefault< Impl > :: DiscreteFunctionType &
184 {
185 if( BaseType::size() != g.size() )
186 DUNE_THROW(InvalidStateException,"DiscreteFunctionDefault: sizes do not match in operator +=");
187
188 const auto end = BaseType::dend();
189 auto git = g.dbegin();
190 for( auto it = BaseType::dbegin(); it != end; ++it, ++git )
191 *it += *git;
192 return asImp();
193 }
194
195
196 template< class Impl >
197 template< class DFType >
198 inline typename DiscreteFunctionDefault< Impl > :: DiscreteFunctionType &
201 {
202 if( BaseType::size() != g.size() )
203 DUNE_THROW(InvalidStateException,"DiscreteFunctionDefault: sizes do not match in operator -=");
204
205 const auto end = BaseType :: dend();
206 auto git = g.dbegin();
207 for( auto it = BaseType :: dbegin(); it != end; ++it, ++git )
208 *it -= *git;
209 return asImp();
210 }
211
212
213 template< class Impl >
214 template< class StreamTraits >
217 {
218 auto versionId = in.readUnsignedInt();
219 if( versionId < DUNE_VERSION_ID(0,9,1) )
220 DUNE_THROW( IOError, "Trying to read outdated file." );
221 else if( versionId > DUNE_MODULE_VERSION_ID(DUNE_FEM) )
222 std :: cerr << "Warning: Reading discrete function from newer version: "
223 << versionId << std :: endl;
224
225 // verify space id for files written with dune-fem version 1.5 or newer
226 if( versionId >= DUNE_VERSION_ID(1,5,0) )
227 {
228 // make sure that space of discrete function matches the space
229 // of the data that was written
230 const auto spaceId = space().type();
231 int mySpaceIdInt;
232 in >> mySpaceIdInt;
233 const auto mySpaceId = static_cast<DFSpaceIdentifier>(mySpaceIdInt);
234
235 if( spaceId != mySpaceId )
236 DUNE_THROW( IOError, "Trying to read discrete function from different space: DFSpace (" << spaceName( spaceId ) << ") != DataSpace (" << spaceName( mySpaceId ) << ")" );
237 }
238
239 // read name
240 in >> name_;
241
242 // read size as integer
243 int32_t mysize;
244 in >> mysize;
245
246 // check size
247 if( static_cast< int >( mysize ) != BaseType::size() )
248 DUNE_THROW( IOError, "Trying to read discrete function of different size." );
249 if( BaseType::size() != static_cast< int >( this->space().size() ) )
250 DUNE_THROW( InvalidStateException, "Trying to read discrete function in uncompressed state." );
251
252 // read all dofs
253 typedef typename DiscreteFunctionSpaceType::BlockMapperType BlockMapperType;
254 typedef typename DiscreteFunctionSpaceType::LocalBlockIndices LocalBlockIndices;
255 typedef typename BlockMapperType::SizeType SizeType;
256
257 const BlockMapperType &blockMapper = this->space().blockMapper();
258 const SizeType nBlocks = blockMapper.size();
259 for( SizeType i = 0; i < nBlocks; ++i )
260 {
261 auto &&block = this->dofVector()[ i ];
262 Hybrid::forEach( LocalBlockIndices(), [ &in, &block ] ( auto j ) { in >> block[ j ]; } );
263 }
264
265 // a discrete function that is part of backup/restore
266 // most likely needs to keep the dofs consistent
267 asImp().enableDofCompression();
268 }
269
270
271 template< class Impl >
272 template< class StreamTraits >
275 {
276 unsigned int versionId = DUNE_MODULE_VERSION_ID(DUNE_FEM);
277 out << versionId ;
278
279 // write space id to for testing when function is read
280 auto spaceId = space().type();
281 out << spaceId ;
282
283 // write name
284 out << name_;
285
286 // only allow write when vector is compressed
287 if( BaseType::size() != static_cast< int >( this->space().size() ) )
288 DUNE_THROW( InvalidStateException, "Trying to write DiscreteFunction in uncompressed state." );
289
290 // write size as 32-bit integer
291 const int32_t mysize = BaseType::size();
292 out << mysize;
293
294 // write all dofs
295 typedef typename DiscreteFunctionSpaceType::BlockMapperType BlockMapperType;
296 typedef typename DiscreteFunctionSpaceType::LocalBlockIndices LocalBlockIndices;
297 typedef typename BlockMapperType::SizeType SizeType;
298
299 const BlockMapperType &blockMapper = this->space().blockMapper();
300 const SizeType nBlocks = blockMapper.size();
301 for( SizeType i = 0; i < nBlocks; ++i )
302 {
303 const auto block = this->dofVector()[ i ];
304 Hybrid::forEach( LocalBlockIndices(), [ &out, &block ] ( auto j ) { out << block[ j ]; } );
305 }
306 }
307
308
309 template< class Impl >
312 {
313 typedef typename DiscreteFunctionSpaceType::IndexSetType IndexSetType;
314 IndexSetType& indexSet = (IndexSetType&)space().indexSet();
316 {
317 auto persistentIndexSet = Dune::Fem::Capabilities::isPersistentIndexSet< IndexSetType >::map( indexSet );
318
319 // this marks the index set in the DofManager's list of index set as persistent
320 if( persistentIndexSet )
321 persistentIndexSet->addBackupRestore();
322 }
323 }
324
325 template< class Impl >
328 {
329 typedef typename DiscreteFunctionSpaceType::IndexSetType IndexSetType;
330 IndexSetType& indexSet = (IndexSetType&)space().indexSet();
332 {
333 auto persistentIndexSet = Dune::Fem::Capabilities::isPersistentIndexSet< IndexSetType >::map( indexSet );
334
335 // this unmarks the index set in the DofManager's list of index set as persistent
336 if( persistentIndexSet )
337 persistentIndexSet->removeBackupRestore();
338 }
339 }
340
341
342 template< class Impl >
343 template< class DFType >
346 {
347 if( BaseType :: size() != g.size() )
348 return false;
349
350 const auto end = BaseType :: dend();
351
352 auto fit = BaseType :: dbegin();
353 auto git = g.dbegin();
354 for( ; fit != end; ++fit, ++git )
355 {
356 if( std::abs( *fit - *git ) > 1e-15 )
357 return false;
358 }
359
360 return true;
361 }
362
363
364 // Stream Operators
365 // ----------------
366
375 template< class Impl >
376 inline std :: ostream &
377 operator<< ( std :: ostream &out,
379 {
380 df.print( out );
381 return out;
382 }
383
384
385
395 template< class StreamTraits, class Impl >
399 {
400 df.write( out );
401 return out;
402 }
403
404
405
415 template< class StreamTraits, class Impl >
416 inline InStreamInterface< StreamTraits > &
419 {
420 df.read( in );
421 return in;
422 }
423
424 } // end namespace Fem
425
426} // end namespace Dune
427#endif // #ifndef DUNE_FEM_DISCRETEFUNCTION_INLINE_HH
DFSpaceIdentifier
enumerator for identification of spaces
Definition: discretefunctionspace.hh:94
std::string spaceName(const DFSpaceIdentifier id)
Definition: discretefunctionspace.hh:109
STL namespace.
Dune::Fem::Double abs(const Dune::Fem::Double &a)
Definition: double.hh:942
Definition: bindguard.hh:11
OutStreamInterface< StreamTraits > & operator<<(OutStreamInterface< StreamTraits > &out, const DiscreteFunctionInterface< Impl > &df)
write a discrete function into an output stream
Definition: discretefunction_inline.hh:397
static double max(const Double &v, const double p)
Definition: double.hh:398
static auto bindGuard(Object &object, Args &&... args) -> std::enable_if_t< isBindable< Object, Args... >::value, BindGuard< Object > >
Definition: bindguard.hh:67
static std::shared_ptr< T > referenceToSharedPtr(T &t)
Definition: memory.hh:19
InStreamInterface< StreamTraits > & operator>>(InStreamInterface< StreamTraits > &in, DiscreteFunctionInterface< Impl > &df)
read a discrete function from an input stream
Definition: discretefunction_inline.hh:417
static void forEach(IndexRange< T, sz > range, F &&f)
Definition: hybrid.hh:129
Definition: common/discretefunction.hh:578
void print(std ::ostream &out) const
print all DoFs to a stream (for debugging purposes)
Definition: discretefunction_inline.hh:91
DiscreteFunctionDefault(const std::string &name, const DiscreteFunctionSpaceType &dfSpace)
Constructor storing discrete function space and local function factory.
Definition: discretefunction_inline.hh:27
DiscreteFunctionType & operator-=(const DiscreteFunctionInterface< DFType > &g)
substract all degrees of freedom from given discrete function using the dof iterators
Definition: discretefunction_inline.hh:200
void write(OutStreamInterface< StreamTraits > &out) const
write the discrete function into a stream
Definition: discretefunction_inline.hh:274
DiscreteFunctionType & operator+=(const DiscreteFunctionInterface< DFType > &g)
add another discrete function to this one
Definition: discretefunction_inline.hh:183
bool dofsValid() const
check for NaNs
Definition: discretefunction_inline.hh:101
DofVectorType::SizeType SizeType
size type of the block vector
Definition: common/discretefunction.hh:646
Traits::LocalDofVectorType LocalDofVectorType
type of LocalDofVector
Definition: common/discretefunction.hh:628
std::type_index assembleOperation_
Definition: common/discretefunction.hh:1042
void evaluateGlobal(const DomainType &x, Functor functor) const
evaluate functor in global coordinate
Definition: discretefunction_inline.hh:162
void read(InStreamInterface< StreamTraits > &in)
read the discrete function from a stream
Definition: discretefunction_inline.hh:216
virtual void insertSubData()
Definition: discretefunction_inline.hh:311
BaseType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
type of discrete function space
Definition: common/discretefunction.hh:600
void assign(const DiscreteFunctionInterface< DFType > &g)
Definition: discretefunction_inline.hh:133
virtual void removeSubData()
Definition: discretefunction_inline.hh:327
CommDataHandle< Operation >::Type dataHandle(const Operation &operation)
return reference to data handle object
void axpy(const RangeFieldType &s, const DiscreteFunctionInterface< DFType > &g)
axpy operation
Definition: discretefunction_inline.hh:117
std::size_t assembleCount_
Definition: common/discretefunction.hh:1043
bool operator==(const DiscreteFunctionInterface< DFType > &g) const
Definition: discretefunction_inline.hh:345
Definition: common/discretefunction.hh:86
void write(OutStreamInterface< StreamTraits > &out) const
write the discrete function into a stream
Definition: common/discretefunction.hh:537
ConstDofIteratorType dbegin() const
obtain an iterator pointing to the first DoF (read-only)
Definition: common/discretefunction.hh:354
void print(std ::ostream &out) const
print all DoFs to a stream (for debugging purposes)
Definition: common/discretefunction.hh:437
void read(InStreamInterface< StreamTraits > &in)
read the discrete function from a stream
Definition: common/discretefunction.hh:527
DiscreteFunctionSpaceType::RangeFieldType RangeFieldType
type of range field, i.e. dof type
Definition: common/discretefunction.hh:109
DiscreteFunctionSpaceType::DomainType DomainType
type of domain, i.e. type of coordinates
Definition: common/discretefunction.hh:111
Traits::DofType DofType
Definition: common/discretefunction.hh:136
DiscreteFunctionSpaceType::GridPartType GridPartType
type of the underlying grid part
Definition: common/discretefunction.hh:118
int size() const
obtain total number of DoFs
Definition: common/discretefunction.hh:333
Definition: entitysearch.hh:131
capability for persistent index sets
Definition: persistentindexset.hh:92
static constexpr PersistentIndexSetInterface * map(IndexSet &indexSet) noexcept
please doc me
Definition: persistentindexset.hh:101
abstract interface for an output stream
Definition: streams.hh:46
abstract interface for an input stream
Definition: streams.hh:179