2 #ifndef DUNE_PDELAB_STOKESDG_VECFEM_HH
3 #define DUNE_PDELAB_STOKESDG_VECFEM_HH
5 #include <dune/common/exceptions.hh>
6 #include <dune/common/fvector.hh>
7 #include <dune/common/static_assert.hh>
8 #include <dune/common/parametertree.hh>
10 #include <dune/geometry/quadraturerules.hh>
12 #include <dune/localfunctions/common/interfaceswitch.hh>
23 #define PBLOCK (- VBLOCK + 1)
28 template<
class Basis,
class Dummy =
void>
35 static const std::size_t
dimRange = Basis::Traits::dimRange;
38 template<
typename Geometry>
39 static void jacobian(
const Basis& basis,
const Geometry& geometry,
42 Geometry::coorddimension> >& jac)
44 jac.resize(basis.size());
45 basis.evaluateJacobian(xl, jac);
52 Basis, typename enable_if<Basis::Traits::dimDomain>::type
58 typedef typename Basis::Traits::RangeFieldType
RangeField;
60 static const std::size_t
dimRange = Basis::Traits::dimRange;
63 template<
typename Geometry>
64 static void jacobian(
const Basis& basis,
const Geometry& geometry,
67 Geometry::coorddimension> >& jac)
69 jac.resize(basis.size());
71 std::vector<FieldMatrix<
73 basis.evaluateJacobian(xl, ljac);
75 const typename Geometry::Jacobian& geojac =
76 geometry.jacobianInverseTransposed(xl);
78 for(std::size_t i = 0; i < basis.size(); ++i)
79 for(std::size_t row=0; row <
dimRange; ++row)
80 geojac.mv(ljac[i][row], jac[i][row]);
94 template<
typename F,
typename B,
typename V,
typename P,
95 typename IP = DefaultInteriorPenalty<typename V::Traits::RangeFieldType> >
108 typedef typename V::Traits::RangeFieldType RF;
122 StokesDGVectorFEM (
const Dune::ParameterTree & configuration,
const IP & ip_factor_,
const RF mu_,
123 const F & _f,
const B & _b,
const V & _v,
const P & _p,
int _qorder=4) :
124 f(_f), b(_b), v(_v), p(_p), qorder(_qorder), mu(mu_), ip_factor(ip_factor_)
126 epsilon = configuration.get<
int>(
"epsilon");
130 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
131 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
134 static const unsigned int dim = EG::Geometry::dimension;
135 static const unsigned int dimw = EG::Geometry::dimensionworld;
139 ((LFSV::CHILDREN == 2),
"You seem to use the wrong function space for StokesDGVectorFEM");
141 typedef typename LFSV::template Child<VBLOCK>::Type LFSV_V;
142 const LFSV_V& lfsv_v = lfsv.template child<VBLOCK>();
143 const LFSV_V& lfsu_v = lfsu.template child<VBLOCK>();
145 typedef typename LFSV::template Child<PBLOCK>::Type LFSV_P;
146 const LFSV_P& lfsv_p = lfsv.template child<PBLOCK>();
147 const LFSV_P& lfsu_p = lfsu.template child<PBLOCK>();
150 typedef FiniteElementInterfaceSwitch<typename LFSV_V::Traits::FiniteElementType > FESwitch_V;
151 typedef BasisInterfaceSwitch<typename FESwitch_V::Basis > BasisSwitch_V;
153 typedef FiniteElementInterfaceSwitch<typename LFSV_P::Traits::FiniteElementType > FESwitch_P;
154 typedef BasisInterfaceSwitch<typename FESwitch_P::Basis > BasisSwitch_P;
155 typedef typename BasisSwitch_V::DomainField DF;
156 typedef typename BasisSwitch_V::RangeField RF;
157 typedef typename BasisSwitch_V::Range Range_V;
158 typedef typename BasisSwitch_P::Range Range_P;
159 typedef typename LFSV::Traits::SizeType size_type;
162 Dune::GeometryType gt = eg.geometry().type();
163 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,qorder);
166 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
169 const Dune::FieldVector<DF,dim> local = it->position();
170 const Dune::FieldVector<DF,dimw> global = eg.geometry().global(local);
171 const RF weight = it->weight() * eg.geometry().integrationElement(it->position());
174 std::vector<Range_V> phi_v(lfsv_v.size());
175 FESwitch_V::basis(lfsv_v.finiteElement()).evaluateFunction(local,phi_v);
178 std::vector<Range_P> phi_p(lfsv_p.size());
179 FESwitch_P::basis(lfsv_p.finiteElement()).evaluateFunction(local,phi_p);
182 std::vector<Dune::FieldMatrix<RF,dim,dim> > jac_phi_v(lfsu_v.size());
183 VectorBasisSwitch_V::jacobian
184 (FESwitch_V::basis(lfsv_v.finiteElement()), eg.geometry(), it->position(), jac_phi_v);
187 std::vector<RF> div_phi_v(lfsv_v.size(),0.0);
188 for (size_type i=0; i<lfsv_v.size(); i++)
189 for (size_type d=0; d<
dim; d++)
190 div_phi_v[i] += jac_phi_v[i][d][d];
194 Dune::FieldMatrix<RF,dim,dim> jac_u(0.0);
196 for (size_type i=0; i<lfsu_v.size(); i++){
197 val_u.axpy(x(lfsu_v,i),phi_v[i]);
198 jac_u.axpy(x(lfsu_v,i),jac_phi_v[i]);
199 div_u += x(lfsu_v,i) * div_phi_v[i];
204 for (size_type i=0; i<lfsu_p.size(); i++)
205 val_p.axpy(x(lfsu_p,i),phi_p[i]);
208 typename F::Traits::RangeType fval;
209 f.evaluateGlobal(global,fval);
212 const RF factor = weight;
213 for (size_type i=0; i<lfsv_v.size(); i++)
214 r.accumulate(lfsv_v, i, (fval * phi_v[i]) * factor );
218 const RF factor = mu * weight;
219 for (size_type i=0; i<lfsv_v.size(); i++){
220 RF dvdu(0); contraction(jac_u,jac_phi_v[i],dvdu);
221 r.accumulate(lfsv_v, i, dvdu * factor);
226 for (size_type i=0; i<lfsv_v.size(); i++)
227 r.accumulate(lfsv_v, i, - div_phi_v[i] * val_p * weight);
231 for (size_type i=0; i<lfsv_p.size(); i++)
232 r.accumulate(lfsv_p, i, - div_u * phi_p[i] * weight);
238 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
240 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
241 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
242 R& r_s, R& r_n)
const
245 static const unsigned int dim = IG::Geometry::dimension;
246 static const unsigned int dimw = IG::Geometry::dimensionworld;
250 ((LFSV::CHILDREN == 2),
"You seem to use the wrong function space for StokesDGVectorFEM");
252 typedef typename LFSV::template Child<VBLOCK>::Type LFSV_V;
253 const LFSV_V& lfsv_v_s = lfsv_s.template child<VBLOCK>();
254 const LFSV_V& lfsu_v_s = lfsu_s.template child<VBLOCK>();
255 const LFSV_V& lfsv_v_n = lfsv_n.template child<VBLOCK>();
256 const LFSV_V& lfsu_v_n = lfsu_n.template child<VBLOCK>();
258 typedef typename LFSV::template Child<PBLOCK>::Type LFSV_P;
259 const LFSV_P& lfsv_p_s = lfsv_s.template child<PBLOCK>();
260 const LFSV_P& lfsu_p_s = lfsu_s.template child<PBLOCK>();
261 const LFSV_P& lfsv_p_n = lfsv_n.template child<PBLOCK>();
262 const LFSV_P& lfsu_p_n = lfsu_n.template child<PBLOCK>();
265 typedef FiniteElementInterfaceSwitch<typename LFSV_V::Traits::FiniteElementType > FESwitch_V;
266 typedef BasisInterfaceSwitch<typename FESwitch_V::Basis > BasisSwitch_V;
268 typedef FiniteElementInterfaceSwitch<typename LFSV_P::Traits::FiniteElementType > FESwitch_P;
269 typedef BasisInterfaceSwitch<typename FESwitch_P::Basis > BasisSwitch_P;
270 typedef typename BasisSwitch_V::DomainField DF;
271 typedef typename BasisSwitch_V::RangeField RF;
272 typedef typename BasisSwitch_V::Range Range_V;
273 typedef typename BasisSwitch_P::Range Range_P;
274 typedef typename LFSV::Traits::SizeType size_type;
277 Dune::GeometryType gt = ig.geometry().type();
278 const Dune::QuadratureRule<DF,dim-1>& rule
279 = Dune::QuadratureRules<DF,dim-1>::rule(gt,qorder);
281 const typename IG::EntityPointer
self = ig.inside();
282 const typename IG::EntityPointer neighbor = ig.outside();
285 for (
typename Dune::QuadratureRule<DF,dim-1>::const_iterator it=rule.begin(); it!=rule.end(); ++it){
288 const Dune::FieldVector<DF,dim> local_s = ig.geometryInInside().global(it->position());
289 const Dune::FieldVector<DF,dim> local_n = ig.geometryInOutside().global(it->position());
290 const Dune::FieldVector<DF,dim> global = ig.geometry().global(it->position());
291 const RF weight = it->weight() * ig.geometry().integrationElement(it->position());
294 std::vector<Range_V> phi_v_s(lfsv_v_s.size());
295 std::vector<Range_V> phi_v_n(lfsv_v_n.size());
296 FESwitch_V::basis(lfsv_v_s.finiteElement()).evaluateFunction(local_s,phi_v_s);
297 FESwitch_V::basis(lfsv_v_n.finiteElement()).evaluateFunction(local_n,phi_v_n);
300 std::vector<Range_P> phi_p_s(lfsv_p_s.size());
301 std::vector<Range_P> phi_p_n(lfsv_p_n.size());
302 FESwitch_P::basis(lfsv_p_s.finiteElement()).evaluateFunction(local_s,phi_p_s);
303 FESwitch_P::basis(lfsv_p_n.finiteElement()).evaluateFunction(local_n,phi_p_n);
306 std::vector<Dune::FieldMatrix<RF,dim,dim> > jac_phi_v_s(lfsu_v_s.size());
307 std::vector<Dune::FieldMatrix<RF,dim,dim> > jac_phi_v_n(lfsu_v_n.size());
308 VectorBasisSwitch_V::jacobian
309 (FESwitch_V::basis(lfsv_v_s.finiteElement()), ig.inside()->geometry(), local_s, jac_phi_v_s);
310 VectorBasisSwitch_V::jacobian
311 (FESwitch_V::basis(lfsv_v_n.finiteElement()), ig.outside()->geometry(), local_n, jac_phi_v_n);
314 std::vector<RF> div_phi_v_s(lfsv_v_s.size(),0.0);
315 std::vector<RF> div_phi_v_n(lfsv_v_s.size(),0.0);
316 for (size_type d=0; d<
dim; d++){
317 for (size_type i=0; i<lfsv_v_s.size(); i++)
318 div_phi_v_s[i] += jac_phi_v_s[i][d][d];
319 for (size_type i=0; i<lfsv_v_n.size(); i++)
320 div_phi_v_n[i] += jac_phi_v_n[i][d][d];
324 Range_V val_u_s(0.0);
325 Range_V val_u_n(0.0);
326 Dune::FieldMatrix<RF,dim,dim> jac_u_s(0.0);
327 Dune::FieldMatrix<RF,dim,dim> jac_u_n(0.0);
330 for (size_type i=0; i<lfsu_v_s.size(); i++){
331 val_u_s.axpy(x_s(lfsu_v_s,i),phi_v_s[i]);
332 jac_u_s.axpy(x_s(lfsu_v_s,i),jac_phi_v_s[i]);
333 div_u_s += x_s(lfsu_v_s,i) * div_phi_v_s[i];
335 for (size_type i=0; i<lfsu_v_n.size(); i++){
336 val_u_n.axpy(x_n(lfsu_v_n,i),phi_v_n[i]);
337 jac_u_n.axpy(x_n(lfsu_v_n,i),jac_phi_v_n[i]);
338 div_u_n += x_n(lfsu_v_n,i) * div_phi_v_n[i];
342 Range_P val_p_s(0.0);
343 Range_P val_p_n(0.0);
344 for (size_type i=0; i<lfsu_p_s.size(); i++)
345 val_p_s.axpy(x_s(lfsu_p_s,i),phi_p_s[i]);
346 for (size_type i=0; i<lfsu_p_n.size(); i++)
347 val_p_n.axpy(x_n(lfsu_p_n,i),phi_p_n[i]);
350 const Dune::FieldVector<DF,dimw> normal = ig.unitOuterNormal(it->position());
353 const Dune::FieldVector<DF,dimw> jump = val_u_s - val_u_n;
358 Dune::FieldVector<DF,dimw> flux(0.0);
359 add_compute_flux(jac_u_s,normal,flux);
360 add_compute_flux(jac_u_n,normal,flux);
363 const RF factor = weight * mu;
364 for (
size_t i=0; i<lfsv_v_s.size(); i++)
365 r_s.accumulate(lfsv_v_s, i, -(flux * phi_v_s[i]) * factor);
366 for (
size_t i=0; i<lfsv_v_n.size(); i++)
367 r_n.accumulate(lfsv_v_n, i, (flux * phi_v_n[i]) * factor);
371 const RF factor = weight * mu;
373 for (
size_t i=0; i<lfsv_v_s.size(); i++){
374 Dune::FieldVector<DF,dimw> flux(0.0);
375 add_compute_flux(jac_phi_v_s[i],normal,flux);
376 r_s.accumulate(lfsv_v_s, i, epsilon * 0.5 * (flux * jump) * factor);
378 for (
size_t i=0; i<lfsv_v_n.size(); i++){
379 Dune::FieldVector<DF,dimw> flux(0.0);
380 add_compute_flux(jac_phi_v_n[i],normal,flux);
381 r_n.accumulate(lfsv_v_n, i, epsilon * 0.5 * (flux * jump) * factor);
386 const RF factor = weight;
387 const RF gamma = ip_factor.getFaceIP(ig);
388 for (
size_t i=0; i<lfsv_v_s.size(); i++)
389 r_s.accumulate(lfsv_v_s,i, gamma * (jump * phi_v_s[i]) * factor);
390 for (
size_t i=0; i<lfsv_v_n.size(); i++)
391 r_n.accumulate(lfsv_v_n,i, - gamma * (jump * phi_v_n[i]) * factor);
395 const RF factor = weight;
396 const RF val_p = 0.5 * (val_p_s + val_p_n);
397 for (
size_t i=0; i<lfsv_v_s.size(); i++)
398 r_s.accumulate(lfsv_v_s, i, val_p * (phi_v_s[i] * normal) * factor);
399 for (
size_t i=0; i<lfsv_v_n.size(); i++)
400 r_n.accumulate(lfsv_v_n, i, - val_p * (phi_v_n[i] * normal) * factor);
402 for (
size_t i=0; i<lfsv_p_s.size(); i++)
403 r_s.accumulate(lfsv_p_s, i, 0.5 * phi_p_s[i] * (jump*normal) * factor);
404 for (
size_t i=0; i<lfsv_p_n.size(); i++)
405 r_n.accumulate(lfsv_p_n, i, 0.5 * phi_p_n[i] * (jump*normal) * factor);
412 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
414 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
418 static const unsigned int dim = IG::Geometry::dimension;
419 static const unsigned int dimw = IG::Geometry::dimensionworld;
423 ((LFSV::CHILDREN == 2),
"You seem to use the wrong function space for StokesDGVectorFEM");
425 typedef typename LFSV::template Child<VBLOCK>::Type LFSV_V;
426 const LFSV_V& lfsv_v_s = lfsv_s.template child<VBLOCK>();
427 const LFSV_V& lfsu_v_s = lfsu_s.template child<VBLOCK>();
429 typedef typename LFSV::template Child<PBLOCK>::Type LFSV_P;
430 const LFSV_P& lfsv_p_s = lfsv_s.template child<PBLOCK>();
431 const LFSV_P& lfsu_p_s = lfsu_s.template child<PBLOCK>();
434 typedef FiniteElementInterfaceSwitch<typename LFSV_V::Traits::FiniteElementType > FESwitch_V;
435 typedef BasisInterfaceSwitch<typename FESwitch_V::Basis > BasisSwitch_V;
437 typedef FiniteElementInterfaceSwitch<typename LFSV_P::Traits::FiniteElementType > FESwitch_P;
438 typedef BasisInterfaceSwitch<typename FESwitch_P::Basis > BasisSwitch_P;
439 typedef typename BasisSwitch_V::DomainField DF;
440 typedef typename BasisSwitch_V::RangeField RF;
441 typedef typename BasisSwitch_V::Range Range_V;
442 typedef typename BasisSwitch_P::Range Range_P;
443 typedef typename LFSV::Traits::SizeType size_type;
446 Dune::GeometryType gt = ig.geometry().type();
447 const Dune::QuadratureRule<DF,dim-1>& rule
448 = Dune::QuadratureRules<DF,dim-1>::rule(gt,qorder);
450 const typename IG::EntityPointer
self = ig.inside();
453 for (
typename Dune::QuadratureRule<DF,dim-1>::const_iterator it=rule.begin(); it!=rule.end(); ++it){
456 const Dune::FieldVector<DF,dim> local_s = ig.geometryInInside().global(it->position());
457 const Dune::FieldVector<DF,dim> global = ig.geometry().global(it->position());
458 const RF weight = it->weight() * ig.geometry().integrationElement(it->position());
461 typename B::Traits::RangeType bctype;
462 b.evaluate(ig,it->position(),bctype);
465 std::vector<Range_V> phi_v_s(lfsv_v_s.size());
466 FESwitch_V::basis(lfsv_v_s.finiteElement()).evaluateFunction(local_s,phi_v_s);
469 std::vector<Range_P> phi_p_s(lfsv_p_s.size());
470 FESwitch_P::basis(lfsv_p_s.finiteElement()).evaluateFunction(local_s,phi_p_s);
473 std::vector<Dune::FieldMatrix<RF,dim,dim> > jac_phi_v_s(lfsu_v_s.size());
474 VectorBasisSwitch_V::jacobian
475 (FESwitch_V::basis(lfsv_v_s.finiteElement()), ig.inside()->geometry(), local_s, jac_phi_v_s);
478 std::vector<RF> div_phi_v_s(lfsv_v_s.size(),0.0);
480 for (size_type d=0; d<
dim; d++){
481 for (size_type i=0; i<lfsv_v_s.size(); i++)
482 div_phi_v_s[i] += jac_phi_v_s[i][d][d];
486 Range_V val_u_s(0.0);
487 Dune::FieldMatrix<RF,dim,dim> jac_u_s(0.0);
489 for (size_type i=0; i<lfsu_v_s.size(); i++){
490 val_u_s.axpy(x_s(lfsu_v_s,i),phi_v_s[i]);
491 jac_u_s.axpy(x_s(lfsu_v_s,i),jac_phi_v_s[i]);
492 div_u_s += x_s(lfsu_v_s,i) * div_phi_v_s[i];
496 Range_P val_p_s(0.0);
497 for (size_type i=0; i<lfsu_p_s.size(); i++)
498 val_p_s.axpy(x_s(lfsu_p_s,i),phi_p_s[i]);
501 const Dune::FieldVector<DF,dimw> normal = ig.unitOuterNormal(it->position());
506 typename V::Traits::RangeType u0;
507 v.evaluateGlobal(global,u0);
508 const Dune::FieldVector<DF,dimw> jump = val_u_s - u0;
513 Dune::FieldVector<DF,dimw> flux(0.0);
514 add_compute_flux(jac_u_s,normal,flux);
516 const RF factor = weight * mu;
517 for (
size_t i=0; i<lfsv_v_s.size(); i++)
518 r_s.accumulate(lfsv_v_s, i, -(flux * phi_v_s[i]) * factor);
522 const RF factor = weight * mu;
524 for (
size_t i=0; i<lfsv_v_s.size(); i++){
525 Dune::FieldVector<DF,dimw> flux(0.0);
526 add_compute_flux(jac_phi_v_s[i],normal,flux);
527 r_s.accumulate(lfsv_v_s, i, epsilon * 0.5 * (flux * jump) * factor);
532 const RF factor = weight;
533 const RF gamma = ip_factor.getFaceIP(ig);
534 for (
size_t i=0; i<lfsv_v_s.size(); i++)
535 r_s.accumulate(lfsv_v_s,i, gamma * (jump * phi_v_s[i]) * factor);
539 const RF factor = weight;
540 const RF val_p = val_p_s;
541 for (
size_t i=0; i<lfsv_v_s.size(); i++)
542 r_s.accumulate(lfsv_v_s, i, val_p * (phi_v_s[i] * normal) * factor);
543 for (
size_t i=0; i<lfsv_p_s.size(); i++)
544 r_s.accumulate(lfsv_p_s, i, 0.5 * phi_p_s[i] * (jump*normal) * factor);
550 typename P::Traits::RangeType p0;
551 p.evaluateGlobal(global,p0);
553 for (
size_t i=0; i<lfsv_v_s.size(); i++){
554 const RF val = p0 * (normal*phi_v_s[i]) * weight;
555 r_s.accumulate(lfsv_v_s, i, val);
565 template<
class M,
class RF>
566 static void contraction(
const M & a,
const M & b, RF & v)
569 const int an = a.N();
570 const int am = a.M();
571 for(
int r=0; r<an; ++r)
572 for(
int c=0; c<am; ++c)
573 v += a[r][c] * b[r][c];
576 template<
class DU,
class R>
577 static void add_compute_flux(
const DU & du,
const R & n, R & result)
579 const int N = du.N();
580 const int M = du.M();
581 for(
int r=0; r<N; ++r)
582 for(
int c=0; c<M; ++c)
583 result[r] += du[r][c] * n[c];
595 const IP & ip_factor;
Definition: stokesdg_vecfem.hh:120
sparsity pattern generator
Definition: pattern.hh:30
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: stokesdg_vecfem.hh:413
Implement jacobian_apply_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:445
Default flags for all local operators.
Definition: flags.hh:18
Implement jacobian_boundary() based on alpha_boundary()
Definition: defaultimp.hh:251
static void jacobian(const Basis &basis, const Geometry &geometry, const DomainLocal &xl, std::vector< FieldMatrix< RangeField, dimRange, Geometry::coorddimension > > &jac)
Compute global jacobian matrix for vector valued bases.
Definition: stokesdg_vecfem.hh:39
void alpha_skeleton(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const X &x_n, const LFSV &lfsv_n, R &r_s, R &r_n) const
Definition: stokesdg_vecfem.hh:239
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538
StokesDGVectorFEM(const Dune::ParameterTree &configuration, const IP &ip_factor_, const RF mu_, const F &_f, const B &_b, const V &_v, const P &_p, int _qorder=4)
Definition: stokesdg_vecfem.hh:122
Basis::Traits::DomainLocal DomainLocal
export vector type of the local coordinates
Definition: stokesdg_vecfem.hh:31
static const int dim
Definition: adaptivity.hh:82
Definition: stokesparameter.hh:32
Definition: stokesdg_vecfem.hh:118
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
static void jacobian(const Basis &basis, const Geometry &geometry, const DomainLocal &xl, std::vector< FieldMatrix< RangeField, dimRange, Geometry::coorddimension > > &jac)
Compute global jacobian matrix for vector valued bases.
Definition: stokesdg_vecfem.hh:64
const IG & ig
Definition: common/constraints.hh:146
static const std::size_t dimRange
export dimension of the values
Definition: stokesdg_vecfem.hh:35
Basis::Traits::RangeField RangeField
export field type of the values
Definition: stokesdg_vecfem.hh:33
Definition: stokesdg_vecfem.hh:115
A local operator for solving the stokes equation using a DG discretization and a vector valued finite...
Definition: stokesdg_vecfem.hh:96
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: stokesdg_vecfem.hh:131
sparsity pattern generator
Definition: pattern.hh:14
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
Definition: stokesdg_vecfem.hh:114
Definition: stokesparameter.hh:36
const EG & eg
Definition: common/constraints.hh:277
IP InteriorPenaltyFactor
Definition: stokesdg_vecfem.hh:111
Definition: stokesdg_vecfem.hh:29
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:157
Basis::Traits::DomainType DomainLocal
export vector type of the local coordinates
Definition: stokesdg_vecfem.hh:56
Definition: stokesparameter.hh:35
Basis::Traits::RangeFieldType RangeField
export field type of the values
Definition: stokesdg_vecfem.hh:58
Definition: stokesdg_vecfem.hh:119