4 #ifndef DUNE_PDELAB_LOCALOPERATOR_LAPLACE_HH 5 #define DUNE_PDELAB_LOCALOPERATOR_LAPLACE_HH 6 #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! 10 #include <dune/common/fvector.hh> 12 #include <dune/geometry/type.hh> 13 #include <dune/geometry/quadraturerules.hh> 15 #include <dune/localfunctions/common/interfaceswitch.hh> 52 DUNE_DEPRECATED_MSG(
"Deprecated in Dune-PDELab 2.4, use the local operator ConvectionDiffusionFEM instead!")
67 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
68 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 71 typedef FiniteElementInterfaceSwitch<
72 typename LFSU::Traits::FiniteElementType
74 typedef BasisInterfaceSwitch<
75 typename FESwitch::Basis
77 typedef typename BasisSwitch::DomainField DF;
78 typedef typename BasisSwitch::RangeField RF;
81 static const int dimLocal = EG::Geometry::mydimension;
82 static const int dimGlobal = EG::Geometry::coorddimension;
85 Dune::GeometryType gt = eg.geometry().type();
86 const Dune::QuadratureRule<DF,dimLocal>& rule =
87 Dune::QuadratureRules<DF,dimLocal>::rule(gt,
quadOrder_);
90 for(
typename Dune::QuadratureRule<DF,dimLocal>::const_iterator it =
91 rule.begin(); it!=rule.end(); ++it)
95 std::vector<Dune::FieldMatrix<RF,1,dimGlobal> >
96 gradphiu(lfsu.size());
97 BasisSwitch::gradient(FESwitch::basis(lfsu.finiteElement()),
98 eg.geometry(), it->position(), gradphiu);
99 std::vector<Dune::FieldMatrix<RF,1,dimGlobal> >
100 gradphiv(lfsv.size());
101 BasisSwitch::gradient(FESwitch::basis(lfsv.finiteElement()),
102 eg.geometry(), it->position(), gradphiv);
105 Dune::FieldVector<RF,dimGlobal> gradu(0.0);
106 for (
size_t i=0; i<lfsu.size(); i++)
107 gradu.axpy(x(lfsu,i),gradphiu[i][0]);
110 RF factor = r.weight() * it->weight() * eg.geometry().integrationElement(it->position());
111 for (
size_t i=0; i<lfsv.size(); i++)
112 r.rawAccumulate(lfsv,i,(gradu*gradphiv[i][0])*factor);
126 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
127 void jacobian_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, M & matrix)
const 130 typedef FiniteElementInterfaceSwitch<
131 typename LFSU::Traits::FiniteElementType
133 typedef BasisInterfaceSwitch<
134 typename FESwitch::Basis
138 typedef typename BasisSwitch::DomainField DF;
139 typedef typename BasisSwitch::RangeField RF;
140 typedef typename LFSU::Traits::SizeType size_type;
143 const int dim = EG::Entity::dimension;
146 Dune::GeometryType gt = eg.geometry().type();
147 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,
quadOrder_);
150 for (
typename QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
152 std::vector<Dune::FieldMatrix<RF,1,dim> > gradphi(lfsu.size());
153 BasisSwitch::gradient(FESwitch::basis(lfsu.finiteElement()),
154 eg.geometry(), it->position(), gradphi);
157 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
159 for (size_type i=0; i<lfsu.size(); i++)
161 for (size_type j=0; j<lfsv.size(); j++)
164 matrix.accumulate(lfsv,j,lfsu,i, gradphi[i][0] * gradphi[j][0] * factor);
Definition: laplace.hh:43
sparsity pattern generator
Definition: pattern.hh:13
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:127
static const int dim
Definition: adaptivity.hh:83
Definition: adaptivity.hh:27
Default flags for all local operators.
Definition: flags.hh:18
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:68
unsigned int quadOrder_
Definition: laplace.hh:172
Definition: laplace.hh:46
Definition: laplace.hh:37