dune-pdelab  2.4.1
laplacedirichletp12d.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_LAPLACEDIRICHLETP12D_HH
3 #define DUNE_PDELAB_LAPLACEDIRICHLETP12D_HH
4 #warning This file is deprecated and will be removed after the Dune-PDELab 2.4 release! Use the ConvectionDiffusionFEM local operator from dune/pdelab/localoperator/convectiondiffusionfem.hh instead!
5 
6 #include <dune/common/exceptions.hh>
7 #include <dune/common/fvector.hh>
8 #include <dune/common/fmatrix.hh>
9 
11 
12 #include"pattern.hh"
13 #include"flags.hh"
14 
15 
16 namespace Dune {
17  namespace PDELab {
18 
19  // a local operator for solving the Laplace equation with Dirichlet boundary conditions
20  // - \Delta u = 0 in \Omega,
21  // u = g on \partial\Omega
22  // with P1 conforming finite elements on triangles
23  class DUNE_DEPRECATED_MSG("Deprecated in Dune-PDELab 2.4, use the local operator ConvectionDiffusionFEM instead!") LaplaceDirichletP12D
24  : public NumericalJacobianApplyVolume<LaplaceDirichletP12D>,
25  public NumericalJacobianVolume<LaplaceDirichletP12D>,
26  public FullVolumePattern,
28  {
29  public:
30  // pattern assembly flags
31  enum { doPatternVolume = true };
32 
33  // residual assembly flags
34  enum { doAlphaVolume = true };
35 
36  // volume integral depending on test and ansatz functions
37  template<typename EG, typename LFSU, typename X, typename LFSV, typename R>
38  void alpha_volume (const EG& eg, const LFSU& lfsu, const X& x, const LFSV& lfsv, R& r) const
39  {
40  // domain and range field type
41  typedef typename LFSU::Traits::FiniteElementType::
42  Traits::LocalBasisType::Traits::DomainFieldType DF;
43  typedef typename LFSU::Traits::FiniteElementType::
44  Traits::LocalBasisType::Traits::RangeFieldType RF;
45 
46  // define integration point (hard coded quadrature)
47  Dune::FieldVector<DF,2> integrationpoint(1.0/3.0);
48 
49  // gradient of shape functions at integration point
50  typedef typename LFSU::Traits::FiniteElementType::
51  Traits::LocalBasisType::Traits::JacobianType JT;
52  std::vector<JT> gradients(lfsu.size());
53  lfsu.finiteElement().localBasis().
54  evaluateJacobian(integrationpoint,gradients);
55 
56  // transformation of gradients to real element
57  const typename EG::Geometry::JacobianInverseTransposed
58  jac = eg.geometry().jacobianInverseTransposed(integrationpoint);
59  Dune::FieldVector<RF,2> gradphi[3];
60  for (int i=0; i<3; i++)
61  {
62  gradphi[i] = 0.0;
63  jac.umv(gradients[i][0],gradphi[i]);
64  }
65 
66  // compute gradient of solution at integration point
67  Dune::FieldVector<RF,2> gradu(0.0);
68  for (int i=0; i<3; i++)
69  gradu.axpy(x[i],gradphi[i]);
70 
71  // integrate grad u * grad phi_i (0.5 is the area of the reference element)
72  RF area = 0.5*eg.geometry().integrationElement(integrationpoint);
73  for (int i=0; i<3; i++)
74  r.accumulate(lfsv, i, (gradu*gradphi[i])*area);
75  }
76  };
77 
79  } // namespace PDELab
80 } // namespace Dune
81 
82 #endif
Definition: convectiondiffusionfem.hh:35
Definition: laplacedirichletp12d.hh:23
sparsity pattern generator
Definition: pattern.hh:13
Definition: adaptivity.hh:27
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: laplacedirichletp12d.hh:38
Default flags for all local operators.
Definition: flags.hh:18
Implement jacobian_volume() based on alpha_volume()
Definition: numericaljacobian.hh:31
Implements linear and nonlinear versions of jacobian_apply_volume() based on alpha_volume() ...
Definition: numericaljacobianapply.hh:32