dune-fem 2.8.0
Loading...
Searching...
No Matches
elliptic.hh
Go to the documentation of this file.
1/**************************************************************************
2
3 The dune-fem module is a module of DUNE (see www.dune-project.org).
4 It is based on the dune-grid interface library
5 extending the grid interface by a number of discretization algorithms
6 for solving non-linear systems of partial differential equations.
7
8 Copyright (C) 2003 - 2015 Robert Kloefkorn
9 Copyright (C) 2003 - 2010 Mario Ohlberger
10 Copyright (C) 2004 - 2015 Andreas Dedner
11 Copyright (C) 2005 Adrian Burri
12 Copyright (C) 2005 - 2015 Mirko Kraenkel
13 Copyright (C) 2006 - 2015 Christoph Gersbacher
14 Copyright (C) 2006 - 2015 Martin Nolte
15 Copyright (C) 2011 - 2015 Tobias Malkmus
16 Copyright (C) 2012 - 2015 Stefan Girke
17 Copyright (C) 2013 - 2015 Claus-Justus Heine
18 Copyright (C) 2013 - 2014 Janick Gerstenberger
19 Copyright (C) 2013 Sven Kaulman
20 Copyright (C) 2013 Tom Ranner
21 Copyright (C) 2015 Marco Agnese
22 Copyright (C) 2015 Martin Alkaemper
23
24
25 The dune-fem module is free software; you can redistribute it and/or
26 modify it under the terms of the GNU General Public License as
27 published by the Free Software Foundation; either version 2 of
28 the License, or (at your option) any later version.
29
30 The dune-fem module is distributed in the hope that it will be useful,
31 but WITHOUT ANY WARRANTY; without even the implied warranty of
32 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
33 GNU General Public License for more details.
34
35 You should have received a copy of the GNU General Public License along
36 with this program; if not, write to the Free Software Foundation, Inc.,
37 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
38
39**************************************************************************/
40#ifndef DUNE_FEM_SCHEMES_ELLIPTIC_HH
41#define DUNE_FEM_SCHEMES_ELLIPTIC_HH
42
43#include <cstddef>
44
45#include <dune/common/timer.hh>
46#include <dune/common/fmatrix.hh>
47
51
55
57
58
59// include parameter handling
62
63// fempy includes
64#include <dune/fempy/quadrature/fempyquadratures.hh>
65
66// EllipticOperator
67// ----------------
68
70template< class DomainDiscreteFunction, class RangeDiscreteFunction, class Model>
72: public virtual Dune::Fem::Operator< DomainDiscreteFunction, RangeDiscreteFunction >
74{
75 typedef DomainDiscreteFunction DomainFunctionType;
76 typedef RangeDiscreteFunction RangeFunctionType;
77 typedef Model ModelType;
78 typedef Model DirichletModelType;
79
80 typedef typename DomainFunctionType::DiscreteFunctionSpaceType DomainDiscreteFunctionSpaceType;
81 typedef typename DomainFunctionType::LocalFunctionType DomainLocalFunctionType;
82 typedef typename DomainLocalFunctionType::RangeType DomainRangeType;
83 typedef typename DomainLocalFunctionType::JacobianRangeType DomainJacobianRangeType;
84 typedef typename RangeFunctionType::DiscreteFunctionSpaceType RangeDiscreteFunctionSpaceType;
85 typedef typename RangeFunctionType::LocalFunctionType RangeLocalFunctionType;
86 typedef typename RangeLocalFunctionType::RangeType RangeRangeType;
87 typedef typename RangeLocalFunctionType::JacobianRangeType RangeJacobianRangeType;
88
89 // the following types must be identical for domain and range
90 typedef typename RangeDiscreteFunctionSpaceType::IteratorType IteratorType;
91 typedef typename IteratorType::Entity EntityType;
92 typedef typename EntityType::Geometry GeometryType;
93 typedef typename RangeDiscreteFunctionSpaceType::DomainType DomainType;
94 typedef typename RangeDiscreteFunctionSpaceType::GridPartType GridPartType;
95 typedef typename GridPartType::IntersectionIteratorType IntersectionIteratorType;
96 typedef typename IntersectionIteratorType::Intersection IntersectionType;
97
100
104 : model_( model ),
105 dSpace_(rangeSpace), rSpace_(rangeSpace)
106 {}
108 const RangeDiscreteFunctionSpaceType &rSpace,
111 : model_( model ),
112 dSpace_(dSpace), rSpace_(rSpace),
114 {}
115
117 virtual void operator() ( const DomainFunctionType &u, RangeFunctionType &w ) const
118 { apply(u,w); }
119 template <class GF>
120 void operator()( const GF &u, RangeFunctionType &w ) const
121 { apply(u,w); }
122 template <class GF>
123 void apply( const GF &u, RangeFunctionType &w ) const;
124
126 {
127 return dSpace_;
128 }
130 {
131 return rSpace_;
132 }
133 const ModelType &model () const { return model_; }
134 ModelType &model () { return model_; }
135
136 void setQuadratureOrders(unsigned int interior, unsigned int surface)
137 {
138 interiorOrder_ = interior;
139 surfaceOrder_ = surface;
140 }
141
142protected:
144private:
145 ModelType &model_;
146 const DomainDiscreteFunctionSpaceType &dSpace_;
147 const RangeDiscreteFunctionSpaceType &rSpace_;
148};
149
150
151
152// DifferentiableEllipticOperator
153// ------------------------------
154
156template< class JacobianOperator, class Model >
158: public EllipticOperator< typename JacobianOperator::DomainFunctionType, typename JacobianOperator::RangeFunctionType, Model >,
159 public Dune::Fem::DifferentiableOperator< JacobianOperator >
161{
163
164 typedef JacobianOperator JacobianOperatorType;
165
169
170 typedef typename DomainFunctionType::DiscreteFunctionSpaceType DomainDiscreteFunctionSpaceType;
171 typedef typename DomainFunctionType::LocalFunctionType DomainLocalFunctionType;
172 typedef typename DomainLocalFunctionType::RangeType DomainRangeType;
173 typedef typename DomainLocalFunctionType::JacobianRangeType DomainJacobianRangeType;
174 typedef typename RangeFunctionType::DiscreteFunctionSpaceType RangeDiscreteFunctionSpaceType;
175 typedef typename RangeFunctionType::LocalFunctionType RangeLocalFunctionType;
176 typedef typename RangeLocalFunctionType::RangeType RangeRangeType;
177 typedef typename RangeLocalFunctionType::JacobianRangeType RangeJacobianRangeType;
178
179 // the following types must be identical for domain and range
180 typedef typename RangeDiscreteFunctionSpaceType::IteratorType IteratorType;
181 typedef typename IteratorType::Entity EntityType;
182 typedef typename EntityType::Geometry GeometryType;
183 typedef typename RangeDiscreteFunctionSpaceType::DomainType DomainType;
184 typedef typename RangeDiscreteFunctionSpaceType::GridPartType GridPartType;
185 typedef typename GridPartType::IntersectionIteratorType IntersectionIteratorType;
186 typedef typename IntersectionIteratorType::Intersection IntersectionType;
187
189 // quadrature for faces - used for Neuman b.c.
191
196 : BaseType( space, space, model, parameter )
197 {}
200 const RangeDiscreteFunctionSpaceType &rSpace,
203 : BaseType( dSpace, rSpace, model, parameter )
204 {}
205
208 { assemble(u,jOp); }
209 template <class GridFunctionType>
210 void jacobian ( const GridFunctionType &u, JacobianOperatorType &jOp ) const
211 { assemble(u,jOp); }
212 template <class GridFunctionType>
213 void assemble ( const GridFunctionType &u, JacobianOperatorType &jOp ) const;
214 using BaseType::operator();
215
216 using BaseType::model;
219};
220
221
222
223// Implementation of EllipticOperator
224// ----------------------------------
225
226template< class DomainDiscreteFunction, class RangeDiscreteFunction, class Model >
227template<class GF>
229 ::apply( const GF &u, RangeFunctionType &w ) const
230{
231 w.clear();
232 // get discrete function space
233 const RangeDiscreteFunctionSpaceType &dfSpace = w.space();
234
235 // get local representation of the discrete functions
238
239 // iterate over grid
240 for( const EntityType &entity : dfSpace )
241 {
242 if( !model().init( entity ) )
243 continue;
244
245 // get elements geometry
246 const GeometryType &geometry = entity.geometry();
247
248 auto uGuard = Dune::Fem::bindGuard( uLocal, entity );
249 auto wGuard = Dune::Fem::bindGuard( wLocal, entity );
250
251 // obtain quadrature order
252 const int quadOrder = interiorOrder_==-1?
253 uLocal.order() + wLocal.order() : interiorOrder_;
254
255 { // element integral
256 QuadratureType quadrature( entity, quadOrder );
257 const size_t numQuadraturePoints = quadrature.nop();
258 for( size_t pt = 0; pt < numQuadraturePoints; ++pt )
259 {
261 const typename QuadratureType::CoordinateType &x = quadrature.point( pt );
262 const double weight = quadrature.weight( pt ) * geometry.integrationElement( x );
263
265 uLocal.evaluate( quadrature[ pt ], vu );
267 uLocal.jacobian( quadrature[ pt ], du );
268
269 // compute mass contribution (studying linear case so linearizing around zero)
270 RangeRangeType avu( 0 );
271 model().source( quadrature[ pt ], vu, du, avu );
272 avu *= weight;
273 // add to local functional wLocal.axpy( quadrature[ pt ], avu );
274
275 RangeJacobianRangeType adu( 0 );
276 // apply diffusive flux
277 model().flux( quadrature[ pt ], vu, du, adu );
278 adu *= weight;
279
280 // add to local function
281 wLocal.axpy( quadrature[ pt ], avu, adu );
283 }
284 }
285
286 if( model().hasNeumanBoundary() && entity.hasBoundaryIntersections() )
287 {
288 const IntersectionIteratorType iitend = dfSpace.gridPart().iend( entity );
289 for( IntersectionIteratorType iit = dfSpace.gridPart().ibegin( entity ); iit != iitend; ++iit )
290 {
291 const IntersectionType &intersection = *iit;
292 if( !intersection.boundary() )
293 continue;
294
295 std::array<int,RangeRangeType::dimension> components(0);
296
297 const bool hasDirichletComponent = model().isDirichletIntersection( intersection, components );
298
299 const auto &intersectionGeometry = intersection.geometry();
300 const int quadOrder = surfaceOrder_==-1?
301 uLocal.order() + wLocal.order() : surfaceOrder_;
302 FaceQuadratureType quadInside( dfSpace.gridPart(), intersection, quadOrder, FaceQuadratureType::INSIDE );
303 const std::size_t numQuadraturePoints = quadInside.nop();
304 for( std::size_t pt = 0; pt < numQuadraturePoints; ++pt )
305 {
306 const auto &x = quadInside.localPoint( pt );
307 double weight = quadInside.weight( pt ) * intersectionGeometry.integrationElement( x );
309 uLocal.evaluate( quadInside[ pt ], vu );
310 RangeRangeType alpha( 0 );
311 model().alpha( quadInside[ pt ], vu, alpha );
312 alpha *= weight;
313 for( int k = 0; k < RangeRangeType::dimension; ++k )
314 if( hasDirichletComponent && components[ k ] )
315 alpha[ k ] = 0;
316 wLocal.axpy( quadInside[ pt ], alpha );
317 }
318 }
319 }
320 }
321
322 w.communicate();
323}
324
325
326
327// Implementation of DifferentiableEllipticOperator
328// ------------------------------------------------
329
330template< class JacobianOperator, class Model >
331template<class GF>
333 ::assemble ( const GF &u, JacobianOperator &jOp ) const
334{
335 // std::cout << "starting assembly\n";
336 // Dune::Timer timer;
337 typedef typename JacobianOperator::LocalMatrixType LocalMatrixType;
338 typedef typename DomainDiscreteFunctionSpaceType::BasisFunctionSetType DomainBasisFunctionSetType;
339 typedef typename RangeDiscreteFunctionSpaceType::BasisFunctionSetType RangeBasisFunctionSetType;
340
341 const DomainDiscreteFunctionSpaceType &domainSpace = jOp.domainSpace();
342 const RangeDiscreteFunctionSpaceType &rangeSpace = jOp.rangeSpace();
343
345 jOp.reserve(stencil);
346 jOp.clear();
347
348 const int domainBlockSize = domainSpace.localBlockSize; // is equal to 1 for scalar functions
349 std::vector< typename DomainLocalFunctionType::RangeType > phi( domainSpace.blockMapper().maxNumDofs()*domainBlockSize );
350 std::vector< typename DomainLocalFunctionType::JacobianRangeType > dphi( domainSpace.blockMapper().maxNumDofs()*domainBlockSize );
351 const int rangeBlockSize = rangeSpace.localBlockSize; // is equal to 1 for scalar functions
352 std::vector< typename RangeLocalFunctionType::RangeType > rphi( rangeSpace.blockMapper().maxNumDofs()*rangeBlockSize );
353 std::vector< typename RangeLocalFunctionType::JacobianRangeType > rdphi( rangeSpace.blockMapper().maxNumDofs()*rangeBlockSize );
354
356 // std::cout << " in assembly: start element loop size=" << rangeSpace.gridPart().grid().size(0) << " time= " << timer.elapsed() << std::endl;;
357 for( const EntityType &entity : rangeSpace )
358 {
359 if( !model().init( entity ) )
360 continue;
361
362 const GeometryType &geometry = entity.geometry();
363
364 auto uGuard = Dune::Fem::bindGuard( uLocal, entity );
365 LocalMatrixType jLocal = jOp.localMatrix( entity, entity );
366
367 const DomainBasisFunctionSetType &domainBaseSet = jLocal.domainBasisFunctionSet();
368 const RangeBasisFunctionSetType &rangeBaseSet = jLocal.rangeBasisFunctionSet();
369 const std::size_t domainNumBasisFunctions = domainBaseSet.size();
370
371 const int quadOrder = interiorOrder_==-1?
372 domainSpace.order() + rangeSpace.order() : interiorOrder_;
373 QuadratureType quadrature( entity, quadOrder );
374 const size_t numQuadraturePoints = quadrature.nop();
375 for( size_t pt = 0; pt < numQuadraturePoints; ++pt )
376 {
378 const typename QuadratureType::CoordinateType &x = quadrature.point( pt );
379 const double weight = quadrature.weight( pt ) * geometry.integrationElement( x );
380
381 // evaluate all basis functions at given quadrature point
382 domainBaseSet.evaluateAll( quadrature[ pt ], phi );
383 rangeBaseSet.evaluateAll( quadrature[ pt ], rphi );
384
385 // evaluate jacobians of all basis functions at given quadrature point
386 domainBaseSet.jacobianAll( quadrature[ pt ], dphi );
387 rangeBaseSet.jacobianAll( quadrature[ pt ], rdphi );
388
389 // get value for linearization
392 uLocal.evaluate( quadrature[ pt ], u0 );
393 uLocal.jacobian( quadrature[ pt ], jacU0 );
394
395 RangeRangeType aphi( 0 );
396 RangeJacobianRangeType adphi( 0 );
397 for( std::size_t localCol = 0; localCol < domainNumBasisFunctions; ++localCol )
398 {
399 // if mass terms or right hand side is present
400 model().linSource( u0, jacU0, quadrature[ pt ], phi[ localCol ], dphi[localCol], aphi );
401
402 // if gradient term is present
403 model().linFlux( u0, jacU0, quadrature[ pt ], phi[ localCol ], dphi[ localCol ], adphi );
404
405 // get column object and call axpy method
406 jLocal.column( localCol ).axpy( rphi, rdphi, aphi, adphi, weight );
407 }
409 }
410
411 if( model().hasNeumanBoundary() && entity.hasBoundaryIntersections() )
412 {
413 const IntersectionIteratorType iitend = rangeSpace.gridPart().iend( entity );
414 for( IntersectionIteratorType iit = rangeSpace.gridPart().ibegin( entity ); iit != iitend; ++iit ) // looping over intersections
415 {
416 const IntersectionType &intersection = *iit;
417 if( !intersection.boundary() )
418 continue;
419
420 std::array<int,RangeRangeType::dimension> components(0);
421 bool hasDirichletComponent = model().isDirichletIntersection( intersection, components );
422
423 const auto &intersectionGeometry = intersection.geometry();
424 const int quadOrder = surfaceOrder_==-1?
425 domainSpace.order() + rangeSpace.order() : surfaceOrder_;
426 FaceQuadratureType quadInside( rangeSpace.gridPart(), intersection, quadOrder, FaceQuadratureType::INSIDE );
427 const std::size_t numQuadraturePoints = quadInside.nop();
428 for( std::size_t pt = 0; pt < numQuadraturePoints; ++pt )
429 {
430 const auto &x = quadInside.localPoint( pt );
431 const double weight = quadInside.weight( pt ) * intersectionGeometry.integrationElement( x );
433 uLocal.evaluate( quadInside[ pt ], u0 );
434 domainBaseSet.evaluateAll( quadInside[ pt ], phi );
435 for( std::size_t localCol = 0; localCol < domainNumBasisFunctions; ++localCol )
436 {
437 RangeRangeType alpha( 0 );
438 model().linAlpha( u0,quadInside[ pt ], phi[ localCol ], alpha );
439 for( int k = 0; k < RangeRangeType::dimension; ++k )
440 if( hasDirichletComponent && components[ k ] )
441 alpha[ k ] = 0;
442 jLocal.column( localCol ).axpy( phi, alpha, weight );
443 }
444 }
445 }
446 }
447 } // end grid traversal
448 // std::cout << " in assembly: final " << timer.elapsed() << std::endl;;
449 jOp.flushAssembly();
450}
451#endif // #ifndef DUNE_FEM_SCHEMES_ELLIPTIC_HH
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:604
static auto bindGuard(Object &object, Args &&... args) -> std::enable_if_t< isBindable< Object, Args... >::value, BindGuard< Object > >
Definition: bindguard.hh:67
Definition: common/localcontribution.hh:14
static ParameterContainer & container()
Definition: io/parameter.hh:193
abstract differentiable operator
Definition: differentiableoperator.hh:29
abstract operator
Definition: operator.hh:34
Stencil contaning the entries (en,en) for all entities in the space. Defailt for an operator over a L...
Definition: stencil.hh:348
quadrature class supporting base function caching
Definition: cachingquadrature.hh:41
[Class for elliptic operator]
Definition: elliptic.hh:74
DomainLocalFunctionType::JacobianRangeType DomainJacobianRangeType
Definition: elliptic.hh:83
RangeLocalFunctionType::JacobianRangeType RangeJacobianRangeType
Definition: elliptic.hh:87
Model ModelType
Definition: elliptic.hh:77
DomainFunctionType::LocalFunctionType DomainLocalFunctionType
Definition: elliptic.hh:81
int interiorOrder_
Definition: elliptic.hh:143
ModelType & model()
Definition: elliptic.hh:134
IteratorType::Entity EntityType
Definition: elliptic.hh:91
RangeDiscreteFunctionSpaceType::GridPartType GridPartType
Definition: elliptic.hh:94
RangeLocalFunctionType::RangeType RangeRangeType
Definition: elliptic.hh:86
RangeFunctionType::DiscreteFunctionSpaceType RangeDiscreteFunctionSpaceType
Definition: elliptic.hh:84
const DomainDiscreteFunctionSpaceType & domainSpace() const
Definition: elliptic.hh:125
RangeDiscreteFunctionSpaceType::IteratorType IteratorType
Definition: elliptic.hh:90
Dune::Fem::CachingQuadrature< GridPartType, 0, Dune::FemPy::FempyQuadratureTraits > QuadratureType
Definition: elliptic.hh:98
EllipticOperator(const RangeDiscreteFunctionSpaceType &rangeSpace, ModelType &model, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: elliptic.hh:101
void apply(const GF &u, RangeFunctionType &w) const
Definition: elliptic.hh:229
GridPartType::IntersectionIteratorType IntersectionIteratorType
Definition: elliptic.hh:95
void setQuadratureOrders(unsigned int interior, unsigned int surface)
Definition: elliptic.hh:136
Model DirichletModelType
Definition: elliptic.hh:78
void operator()(const GF &u, RangeFunctionType &w) const
Definition: elliptic.hh:120
DomainDiscreteFunction DomainFunctionType
Definition: elliptic.hh:75
RangeFunctionType::LocalFunctionType RangeLocalFunctionType
Definition: elliptic.hh:85
virtual void operator()(const DomainFunctionType &u, RangeFunctionType &w) const
application operator
Definition: elliptic.hh:117
RangeDiscreteFunctionSpaceType::DomainType DomainType
Definition: elliptic.hh:93
DomainFunctionType::DiscreteFunctionSpaceType DomainDiscreteFunctionSpaceType
Definition: elliptic.hh:80
Dune::Fem::CachingQuadrature< GridPartType, 1, Dune::FemPy::FempyQuadratureTraits > FaceQuadratureType
Definition: elliptic.hh:99
IntersectionIteratorType::Intersection IntersectionType
Definition: elliptic.hh:96
EntityType::Geometry GeometryType
Definition: elliptic.hh:92
const ModelType & model() const
Definition: elliptic.hh:133
int surfaceOrder_
Definition: elliptic.hh:143
const RangeDiscreteFunctionSpaceType & rangeSpace() const
Definition: elliptic.hh:129
RangeDiscreteFunction RangeFunctionType
Definition: elliptic.hh:76
DomainLocalFunctionType::RangeType DomainRangeType
Definition: elliptic.hh:82
EllipticOperator(const DomainDiscreteFunctionSpaceType &dSpace, const RangeDiscreteFunctionSpaceType &rSpace, ModelType &model, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: elliptic.hh:107
[Class for linearizable elliptic operator]
Definition: elliptic.hh:161
DomainFunctionType::DiscreteFunctionSpaceType DomainDiscreteFunctionSpaceType
Definition: elliptic.hh:170
DifferentiableEllipticOperator(const DomainDiscreteFunctionSpaceType &dSpace, const RangeDiscreteFunctionSpaceType &rSpace, ModelType &model, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
contructor
Definition: elliptic.hh:199
EllipticOperator< typename JacobianOperator::DomainFunctionType, typename JacobianOperator::RangeFunctionType, Model > BaseType
Definition: elliptic.hh:162
RangeDiscreteFunctionSpaceType::IteratorType IteratorType
Definition: elliptic.hh:180
DifferentiableEllipticOperator(const RangeDiscreteFunctionSpaceType &space, ModelType &model, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
contructor
Definition: elliptic.hh:193
BaseType::DomainFunctionType DomainFunctionType
Definition: elliptic.hh:166
void assemble(const GridFunctionType &u, JacobianOperatorType &jOp) const
IntersectionIteratorType::Intersection IntersectionType
Definition: elliptic.hh:186
DomainLocalFunctionType::JacobianRangeType DomainJacobianRangeType
Definition: elliptic.hh:173
BaseType::RangeFunctionType RangeFunctionType
Definition: elliptic.hh:167
RangeLocalFunctionType::JacobianRangeType RangeJacobianRangeType
Definition: elliptic.hh:177
void jacobian(const GridFunctionType &u, JacobianOperatorType &jOp) const
Definition: elliptic.hh:210
RangeDiscreteFunctionSpaceType::GridPartType GridPartType
Definition: elliptic.hh:184
RangeFunctionType::DiscreteFunctionSpaceType RangeDiscreteFunctionSpaceType
Definition: elliptic.hh:174
RangeDiscreteFunctionSpaceType::DomainType DomainType
Definition: elliptic.hh:183
BaseType::QuadratureType QuadratureType
Definition: elliptic.hh:188
BaseType::FaceQuadratureType FaceQuadratureType
Definition: elliptic.hh:190
DomainFunctionType::LocalFunctionType DomainLocalFunctionType
Definition: elliptic.hh:171
RangeFunctionType::LocalFunctionType RangeLocalFunctionType
Definition: elliptic.hh:175
JacobianOperator JacobianOperatorType
Definition: elliptic.hh:164
const ModelType & model() const
Definition: elliptic.hh:133
DomainLocalFunctionType::RangeType DomainRangeType
Definition: elliptic.hh:172
GridPartType::IntersectionIteratorType IntersectionIteratorType
Definition: elliptic.hh:185
EntityType::Geometry GeometryType
Definition: elliptic.hh:182
BaseType::ModelType ModelType
Definition: elliptic.hh:168
IteratorType::Entity EntityType
Definition: elliptic.hh:181
void jacobian(const DomainFunctionType &u, JacobianOperatorType &jOp) const
method to setup the jacobian of the operator for storage in a matrix
Definition: elliptic.hh:207
RangeLocalFunctionType::RangeType RangeRangeType
Definition: elliptic.hh:176