2 #ifndef DUNE_PDELAB_LAPLACEDIRICHLETP12D_HH
3 #define DUNE_PDELAB_LAPLACEDIRICHLETP12D_HH
5 #include <dune/common/exceptions.hh>
6 #include <dune/common/fvector.hh>
7 #include <dune/common/fmatrix.hh>
8 #include <dune/common/static_assert.hh>
36 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
37 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
40 typedef typename LFSU::Traits::FiniteElementType::
41 Traits::LocalBasisType::Traits::DomainFieldType DF;
42 typedef typename LFSU::Traits::FiniteElementType::
43 Traits::LocalBasisType::Traits::RangeFieldType RF;
46 Dune::FieldVector<DF,2> integrationpoint(1.0/3.0);
49 typedef typename LFSU::Traits::FiniteElementType::
50 Traits::LocalBasisType::Traits::JacobianType JT;
51 std::vector<JT> gradients(lfsu.size());
52 lfsu.finiteElement().localBasis().
53 evaluateJacobian(integrationpoint,gradients);
56 const typename EG::Geometry::JacobianInverseTransposed
57 jac = eg.geometry().jacobianInverseTransposed(integrationpoint);
58 Dune::FieldVector<RF,2> gradphi[3];
59 for (
int i=0; i<3; i++)
62 jac.umv(gradients[i][0],gradphi[i]);
66 Dune::FieldVector<RF,2> gradu(0.0);
67 for (
int i=0; i<3; i++)
68 gradu.axpy(x[i],gradphi[i]);
71 RF area = 0.5*eg.geometry().integrationElement(integrationpoint);
72 for (
int i=0; i<3; i++)
73 r.accumulate(lfsv, i, (gradu*gradphi[i])*area);
Default flags for all local operators.
Definition: flags.hh:18
Definition: laplacedirichletp12d.hh:23
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: laplacedirichletp12d.hh:37
Definition: laplacedirichletp12d.hh:33
sparsity pattern generator
Definition: pattern.hh:14
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
const EG & eg
Definition: common/constraints.hh:277
Definition: laplacedirichletp12d.hh:30