4 #ifndef DUNE_PDELAB_POISSON_HH
5 #define DUNE_PDELAB_POISSON_HH
9 #include<dune/common/exceptions.hh>
10 #include<dune/common/fvector.hh>
11 #include<dune/common/static_assert.hh>
13 #include<dune/geometry/type.hh>
14 #include<dune/geometry/quadraturerules.hh>
16 #include <dune/localfunctions/common/interfaceswitch.hh>
43 template<
typename F,
typename B,
typename J>
61 Poisson (
const F& f_,
const B& bctype_,
const J& j_,
unsigned int quadOrder)
68 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
69 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
84 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
85 void jacobian_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, M & matrix)
const
90 template<
typename EG,
typename LFSV,
typename R>
94 typedef FiniteElementInterfaceSwitch<
95 typename LFSV::Traits::FiniteElementType
97 typedef BasisInterfaceSwitch<
98 typename FESwitch::Basis
100 typedef typename BasisSwitch::DomainField DF;
101 typedef typename BasisSwitch::RangeField RF;
102 typedef typename BasisSwitch::Range Range;
105 static const int dimLocal = EG::Geometry::mydimension;
108 Dune::GeometryType gt = eg.geometry().type();
109 const Dune::QuadratureRule<DF,dimLocal>& rule =
110 Dune::QuadratureRules<DF,dimLocal>::rule(gt,
quadOrder_);
113 for (
typename Dune::QuadratureRule<DF,dimLocal>::const_iterator it =
114 rule.begin(); it!=rule.end(); ++it)
117 std::vector<Range> phi(lfsv.size());
118 FESwitch::basis(lfsv.finiteElement()).
119 evaluateFunction(it->position(),phi);
122 typename F::Traits::RangeType y(0.0);
123 f.evaluate(eg.entity(),it->position(),y);
126 RF factor = - r.weight() * it->weight() * eg.geometry().integrationElement(it->position());
127 for (
size_t i=0; i<lfsv.size(); i++)
128 r.rawAccumulate(lfsv,i,y*phi[i]*factor);
133 template<
typename IG,
typename LFSV,
typename R>
137 typedef FiniteElementInterfaceSwitch<
138 typename LFSV::Traits::FiniteElementType
140 typedef BasisInterfaceSwitch<
141 typename FESwitch::Basis
143 typedef typename BasisSwitch::DomainField DF;
144 typedef typename BasisSwitch::DomainLocal DomainLocal;
145 typedef typename BasisSwitch::RangeField RF;
146 typedef typename BasisSwitch::Range Range;
149 static const int dimLocal = IG::Geometry::mydimension;
152 Dune::GeometryType gtface = ig.geometryInInside().type();
153 const Dune::QuadratureRule<DF,dimLocal>& rule =
154 Dune::QuadratureRules<DF,dimLocal>::rule(gtface,
quadOrder_);
157 for (
typename Dune::QuadratureRule<DF,dimLocal>::const_iterator it =
158 rule.begin(); it!=rule.end(); ++it)
162 if( !
bctype.isNeumann( ig,it->position() ) )
166 const DomainLocal& local =
167 ig.geometryInInside().global(it->position());
170 std::vector<Range> phi(lfsv.size());
171 FESwitch::basis(lfsv.finiteElement()).evaluateFunction(local,phi);
174 typename J::Traits::RangeType y(0.0);
175 j.evaluate(*(ig.inside()),local,y);
178 RF factor = r.weight() * it->weight()*ig.geometry().integrationElement(it->position());
179 for (
size_t i=0; i<lfsv.size(); i++)
180 r.rawAccumulate(lfsv,i,y*phi[i]*factor);
212 template<
typename Time,
typename F,
typename B,
typename J>
229 :
Base(f_, bctype_, j_, quadOrder)
Definition: poisson.hh:55
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &matrix) const
Compute the Laplace stiffness matrix for the element given in 'eg'.
Definition: laplace.hh:126
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: poisson.hh:69
Default flags for all local operators.
Definition: flags.hh:18
const J & j
Definition: poisson.hh:187
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Compute Laplace matrix times a given vector for one element.
Definition: laplace.hh:67
const B & bctype
Definition: poisson.hh:186
unsigned int quadOrder_
Definition: poisson.hh:193
a local operator for solving the Poisson equation in instationary problems
Definition: poisson.hh:213
J & j
Definition: poisson.hh:224
const IG & ig
Definition: common/constraints.hh:146
void lambda_boundary(const IG &ig, const LFSV &lfsv, R &r) const
Definition: poisson.hh:134
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &matrix) const
Compute the Laplace stiffness matrix for the element given in 'eg'.
Definition: poisson.hh:85
F & f
Definition: poisson.hh:222
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
Poisson(const F &f_, const B &bctype_, const J &j_, unsigned int quadOrder)
Constructor.
Definition: poisson.hh:61
const F & f
Definition: poisson.hh:185
Definition: poisson.hh:44
sparsity pattern generator
Definition: pattern.hh:14
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
Definition: poisson.hh:50
Definition: poisson.hh:54
const EG & eg
Definition: common/constraints.hh:277
InstationaryPoisson(F &f_, B &bctype_, J &j_, unsigned int quadOrder)
construct InstationaryPoisson
Definition: poisson.hh:228
Definition: poisson.hh:53
Laplace laplace_
Definition: poisson.hh:190
void setTime(Time t)
set the time for subsequent evaluation on the parameter functions
Definition: poisson.hh:236
B & bctype
Definition: poisson.hh:223
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: poisson.hh:91
Definition: laplace.hh:37
void setTime(Timet_)
set time for subsequent evaluation
Definition: idefault.hh:104