dune-pdelab  2.4-dev
poisson.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_POISSON_HH
5 #define DUNE_PDELAB_POISSON_HH
6 #warning "The file dune/pdelab/localoperator/poisson.hh is deprecated. Please use the ConvectionDiffusionFEM local operator from dune/pdelab/localoperator/convectiondiffusionfem.hh instead."
7 
8 #include<vector>
9 
10 #include<dune/common/deprecated.hh>
11 #include<dune/common/exceptions.hh>
12 #include<dune/common/fvector.hh>
13 
14 #include<dune/geometry/type.hh>
15 #include<dune/geometry/quadraturerules.hh>
16 
17 #include <dune/localfunctions/common/interfaceswitch.hh>
18 
20 
21 #include"defaultimp.hh"
22 #include"idefault.hh"
23 #include"pattern.hh"
24 #include"flags.hh"
25 
26 namespace Dune {
27  namespace PDELab {
31 
44  template<typename F, typename B, typename J>
45  class DUNE_DEPRECATED_MSG("Deprecated in DUNE-PDELab 2.4. Please use ConvectionDiffusionFEM instead!") Poisson
46  : public NumericalJacobianApplyVolume<Poisson<F,B,J> >,
47  public FullVolumePattern,
49  {
50  public:
51  // pattern assembly flags
52  enum { doPatternVolume = true };
53 
54  // residual assembly flags
55  enum { doAlphaVolume = true };
56  enum { doLambdaVolume = true };
57  enum { doLambdaBoundary = true };
58 
63  DUNE_DEPRECATED_MSG("Deprecated in DUNE-PDELab 2.4, use the local operator ConvectionDiffusionFEM instead!")
64  Poisson (const F& f_, const B& bctype_, const J& j_, unsigned int quadOrder)
65  : f(f_), bctype(bctype_), j(j_),
66  laplace_(quadOrder),
67  quadOrder_(quadOrder)
68  {}
69 
70  // volume integral depending on test and ansatz functions
71  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
72  void alpha_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, R& r) const
73  {
74  laplace_.alpha_volume(eg, lfsu, x, lfsv, r);
75  }
76 
87  template<typename EG, typename LFSU, typename X, typename LFSV, typename M>
88  void jacobian_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, M & matrix) const
89  {
90  laplace_.jacobian_volume(eg, lfsu, x, lfsv, matrix);
91  }
92  // volume integral depending only on test functions
93  template<typename EG, typename LFSV, typename R>
94  void lambda_volume (const EG& eg, const LFSV& lfsv, R& r) const
95  {
96  // domain and range field type
97  typedef FiniteElementInterfaceSwitch<
98  typename LFSV::Traits::FiniteElementType
99  > FESwitch;
100  typedef BasisInterfaceSwitch<
101  typename FESwitch::Basis
102  > BasisSwitch;
103  typedef typename BasisSwitch::DomainField DF;
104  typedef typename BasisSwitch::RangeField RF;
105  typedef typename BasisSwitch::Range Range;
106 
107  // dimensions
108  static const int dimLocal = EG::Geometry::mydimension;
109 
110  // select quadrature rule
111  Dune::GeometryType gt = eg.geometry().type();
112  const Dune::QuadratureRule<DF,dimLocal>& rule =
113  Dune::QuadratureRules<DF,dimLocal>::rule(gt,quadOrder_);
114 
115  // loop over quadrature points
116  for (typename Dune::QuadratureRule<DF,dimLocal>::const_iterator it =
117  rule.begin(); it!=rule.end(); ++it)
118  {
119  // evaluate shape functions
120  std::vector<Range> phi(lfsv.size());
121  FESwitch::basis(lfsv.finiteElement()).
122  evaluateFunction(it->position(),phi);
123 
124  // evaluate right hand side parameter function
125  typename F::Traits::RangeType y(0.0);
126  f.evaluate(eg.entity(),it->position(),y);
127 
128  // integrate f
129  RF factor = - r.weight() * it->weight() * eg.geometry().integrationElement(it->position());
130  for (size_t i=0; i<lfsv.size(); i++)
131  r.rawAccumulate(lfsv,i,y*phi[i]*factor);
132  }
133  }
134 
135  // boundary integral independen of ansatz functions
136  template<typename IG, typename LFSV, typename R>
137  void lambda_boundary (const IG& ig, const LFSV& lfsv, R& r) const
138  {
139  // domain and range field type
140  typedef FiniteElementInterfaceSwitch<
141  typename LFSV::Traits::FiniteElementType
142  > FESwitch;
143  typedef BasisInterfaceSwitch<
144  typename FESwitch::Basis
145  > BasisSwitch;
146  typedef typename BasisSwitch::DomainField DF;
147  typedef typename BasisSwitch::DomainLocal DomainLocal;
148  typedef typename BasisSwitch::RangeField RF;
149  typedef typename BasisSwitch::Range Range;
150 
151  // dimensions
152  static const int dimLocal = IG::Geometry::mydimension;
153 
154  // select quadrature rule
155  Dune::GeometryType gtface = ig.geometryInInside().type();
156  const Dune::QuadratureRule<DF,dimLocal>& rule =
157  Dune::QuadratureRules<DF,dimLocal>::rule(gtface,quadOrder_);
158 
159  // loop over quadrature points and integrate normal flux
160  for (typename Dune::QuadratureRule<DF,dimLocal>::const_iterator it =
161  rule.begin(); it!=rule.end(); ++it)
162  {
163  // evaluate boundary condition type
164  // skip rest if we are on Dirichlet boundary
165  if( !bctype.isNeumann( ig,it->position() ) )
166  continue;
167 
168  // position of quadrature point in local coordinates of element
169  const DomainLocal& local =
170  ig.geometryInInside().global(it->position());
171 
172  // evaluate test shape functions
173  std::vector<Range> phi(lfsv.size());
174  FESwitch::basis(lfsv.finiteElement()).evaluateFunction(local,phi);
175 
176  // evaluate flux boundary condition
177  typename J::Traits::RangeType y(0.0);
178  j.evaluate(*(ig.inside()),local,y);
179 
180  // integrate J
181  RF factor = r.weight() * it->weight()*ig.geometry().integrationElement(it->position());
182  for (size_t i=0; i<lfsv.size(); i++)
183  r.rawAccumulate(lfsv,i,y*phi[i]*factor);
184  }
185  }
186 
187  protected:
188  const F& f;
189  const B& bctype;
190  const J& j;
191 
192  // Laplace assembler to handle the matrix assembly
194 
195  // Quadrature rule order
196  unsigned int quadOrder_;
197  };
198 
201 
215  template<typename Time, typename F, typename B, typename J>
216  class DUNE_DEPRECATED_MSG("Deprecated in DUNE-PDELab 2.4. Please use ConvectionDiffusionFEM instead!") InstationaryPoisson
217  : public Poisson<F,B,J>,
219  {
220  typedef Poisson<F,B,J> Base;
222 
223  protected:
224  // need non-const references for setTime()
225  F& f;
226  B& bctype;
227  J& j;
228 
229  public:
231  InstationaryPoisson(F& f_, B& bctype_, J& j_, unsigned int quadOrder)
232  : Base(f_, bctype_, j_, quadOrder)
233  , f(f_)
234  , bctype(bctype_)
235  , j(j_)
236  {}
237 
239  void setTime (Time t) {
240  f.setTime(t);
241  bctype.setTime(t);
242  j.setTime(t);
243  IDefault::setTime(t);
244  }
245  };
246 
248  } // namespace PDELab
249 } // namespace Dune
250 
251 #endif
F & f
Definition: poisson.hh:225
const F & f
Definition: poisson.hh:188
J & j
Definition: poisson.hh:227
InstationaryPoisson(F &f_, B &bctype_, J &j_, unsigned int quadOrder)
construct InstationaryPoisson
Definition: poisson.hh:231
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
Laplace laplace_
Definition: poisson.hh:193
void setTime(Time t)
set the time for subsequent evaluation on the parameter functions
Definition: poisson.hh:239
B & bctype
Definition: poisson.hh:226
const IG & ig
Definition: constraints.hh:147
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:88
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:13
void lambda_boundary(const IG &ig, const LFSV &lfsv, R &r) const
Definition: poisson.hh:137
Definition: poisson.hh:45
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: poisson.hh:94
Definition: convectiondiffusionfem.hh:39
a local operator for solving the Poisson equation in instationary problems
Definition: poisson.hh:216
unsigned int quadOrder_
Definition: poisson.hh:196
Definition: adaptivity.hh:27
const J & j
Definition: poisson.hh:190
Definition: laplace.hh:36
const EG & eg
Definition: constraints.hh:280
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: poisson.hh:72
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
const B & bctype
Definition: poisson.hh:189