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>
35 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
36 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
39 typedef typename LFSU::Traits::FiniteElementType::
40 Traits::LocalBasisType::Traits::DomainFieldType DF;
41 typedef typename LFSU::Traits::FiniteElementType::
42 Traits::LocalBasisType::Traits::RangeFieldType RF;
45 Dune::FieldVector<DF,2> integrationpoint(1.0/3.0);
48 typedef typename LFSU::Traits::FiniteElementType::
49 Traits::LocalBasisType::Traits::JacobianType JT;
50 std::vector<JT> gradients(lfsu.size());
51 lfsu.finiteElement().localBasis().
52 evaluateJacobian(integrationpoint,gradients);
55 const typename EG::Geometry::JacobianInverseTransposed
56 jac = eg.geometry().jacobianInverseTransposed(integrationpoint);
57 Dune::FieldVector<RF,2> gradphi[3];
58 for (
int i=0; i<3; i++)
61 jac.umv(gradients[i][0],gradphi[i]);
65 Dune::FieldVector<RF,2> gradu(0.0);
66 for (
int i=0; i<3; i++)
67 gradu.axpy(x[i],gradphi[i]);
70 RF area = 0.5*eg.geometry().integrationElement(integrationpoint);
71 for (
int i=0; i<3; i++)
72 r.accumulate(lfsv, i, (gradu*gradphi[i])*area);
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
Definition: laplacedirichletp12d.hh:29
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:13
Definition: adaptivity.hh:27
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
Definition: laplacedirichletp12d.hh:22
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: laplacedirichletp12d.hh:36
Definition: laplacedirichletp12d.hh:32
const EG & eg
Definition: constraints.hh:280