2 #ifndef DUNE_PDELAB_LINEARACOUSTICSDG_HH
3 #define DUNE_PDELAB_LINEARACOUSTICSDG_HH
7 #include<dune/common/exceptions.hh>
8 #include<dune/common/fvector.hh>
9 #include<dune/common/static_assert.hh>
10 #include<dune/geometry/referenceelements.hh>
42 template<
typename T1,
typename T2,
typename T3>
45 RT[0][0] = 1; RT[0][1] = c*n[0];
46 RT[1][0] = -1; RT[1][1] = c*n[0];
49 template<
typename T1,
typename T2,
typename T3>
50 static void eigenvectors (T1 c,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,dim+1,dim+1>& RT)
52 RT[0][0] = 1; RT[1][0] = c*n[0];
53 RT[0][1] = -1; RT[1][1] = c*n[0];
71 template<
typename T1,
typename T2,
typename T3>
74 RT[0][0] = 0; RT[0][1] = -n[1]; RT[0][2] = n[0];
75 RT[1][0] = 1; RT[1][1] = c*n[0]; RT[1][2] = c*n[1];
76 RT[2][0] = -1; RT[2][1] = c*n[0]; RT[2][2] = c*n[1];
79 template<
typename T1,
typename T2,
typename T3>
80 static void eigenvectors (T1 c,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,dim+1,dim+1>& RT)
82 RT[0][0] = 0; RT[1][0] = -n[1]; RT[2][0] = n[0];
83 RT[0][1] = 1; RT[1][1] = c*n[0]; RT[2][1] = c*n[1];
84 RT[0][2] = -1; RT[1][2] = c*n[0]; RT[2][2] = c*n[1];
101 template<
typename T1,
typename T2,
typename T3>
104 Dune::FieldVector<T2,dim>
s;
105 s[0] = n[1]-n[2]; s[1] = n[2]-n[0]; s[2] = n[0]-n[1];
106 if (s.two_norm()<1
e-5)
108 s[0] = n[1]+n[2]; s[1] = n[2]-n[0]; s[2] = -(n[0]+n[1]);
111 Dune::FieldVector<T2,dim> t;
112 t[0] = n[1]*s[2] - n[2]*s[1];
113 t[1] = n[2]*s[0] - n[0]*s[2];
114 t[2] = n[0]*s[1] - n[1]*s[0];
116 RT[0][0] = 0; RT[0][1] = s[0]; RT[0][2] = s[1]; RT[0][3] = s[2];
117 RT[1][0] = 0; RT[1][1] = t[0]; RT[1][2] = t[1]; RT[1][3] = t[2];
118 RT[2][0] = 1; RT[2][1] = c*n[0]; RT[2][2] = c*n[1]; RT[2][3] = c*n[2];
119 RT[3][0] = -1; RT[3][1] = c*n[0]; RT[3][2] = c*n[1]; RT[3][3] = c*n[2];
122 template<
typename T1,
typename T2,
typename T3>
123 static void eigenvectors (T1 c,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,dim+1,dim+1>& RT)
125 Dune::FieldVector<T2,dim>
s;
126 s[0] = n[1]-n[2]; s[1] = n[2]-n[0]; s[2] = n[0]-n[1];
127 if (s.two_norm()<1
e-5)
129 s[0] = n[1]+n[2]; s[1] = n[2]-n[0]; s[2] = -(n[0]+n[1]);
132 Dune::FieldVector<T2,dim> t;
133 t[0] = n[1]*s[2] - n[2]*s[1];
134 t[1] = n[2]*s[0] - n[0]*s[2];
135 t[2] = n[0]*s[1] - n[1]*s[0];
137 RT[0][0] = 0; RT[1][0] = s[0]; RT[2][0] = s[1]; RT[3][0] = s[2];
138 RT[0][1] = 0; RT[1][1] = t[0]; RT[2][1] = t[1]; RT[3][1] = t[2];
139 RT[0][2] = 1; RT[1][2] = c*n[0]; RT[2][2] = c*n[1]; RT[3][2] = c*n[2];
140 RT[0][3] = -1; RT[1][3] = c*n[0]; RT[2][3] = c*n[1]; RT[3][3] = c*n[2];
160 template<
typename T,
typename FEM>
173 enum { dim = T::Traits::GridViewType::dimension };
188 : param(param_), overintegration(overintegration_), cache(20)
193 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
194 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
197 typedef typename LFSV::template Child<0>::Type DGSpace;
198 typedef typename DGSpace::Traits::FiniteElementType::
199 Traits::LocalBasisType::Traits::DomainFieldType DF;
200 typedef typename DGSpace::Traits::FiniteElementType::
201 Traits::LocalBasisType::Traits::RangeFieldType RF;
202 typedef typename DGSpace::Traits::FiniteElementType::
203 Traits::LocalBasisType::Traits::RangeType RangeType;
204 typedef typename DGSpace::Traits::FiniteElementType::
205 Traits::LocalBasisType::Traits::JacobianType JacobianType;
206 typedef typename DGSpace::Traits::SizeType size_type;
209 if (LFSV::CHILDREN!=dim+1)
210 DUNE_THROW(Dune::Exception,
"need exactly dim+1 components!");
213 const DGSpace& dgspace = lfsv.template child<0>();
216 const int order = dgspace.finiteElement().localBasis().order();
217 const int intorder = overintegration+2*order;
218 Dune::GeometryType gt = eg.geometry().type();
219 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
222 typename EG::Geometry::JacobianInverseTransposed jac;
225 Dune::FieldVector<DF,dim> localcenter = Dune::ReferenceElements<DF,dim>::general(gt).position(0,0);
226 RF c2 = param.c(eg.entity(),localcenter);
232 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
235 const std::vector<RangeType>& phi = cache[order].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
238 Dune::FieldVector<RF,dim+1> u(0.0);
239 for (size_type k=0; k<=dim; k++)
240 for (size_type j=0; j<dgspace.size(); j++)
241 u[k] += x(lfsv.child(k),j)*phi[j];
245 const std::vector<JacobianType>& js = cache[order].evaluateJacobian(it->position(),dgspace.finiteElement().localBasis());
248 jac = eg.geometry().jacobianInverseTransposed(it->position());
249 std::vector<Dune::FieldVector<RF,dim> > gradphi(dgspace.size());
250 for (size_type i=0; i<dgspace.size(); i++)
251 jac.mv(js[i][0],gradphi[i]);
254 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
255 for (size_type k=0; k<dgspace.size(); k++)
258 for (size_type j=0; j<dim; j++)
259 r.accumulate(lfsv.child(0),k, - u[j+1]*gradphi[k][j]*factor);
261 for (size_type i=1; i<=dim; i++)
262 r.accumulate(lfsv.child(i),k, - c2*u[0]*gradphi[k][i-1]*factor);
272 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
274 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
275 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
276 R& r_s, R& r_n)
const
279 typedef typename LFSV::template Child<0>::Type DGSpace;
280 typedef typename DGSpace::Traits::FiniteElementType::
281 Traits::LocalBasisType::Traits::DomainFieldType DF;
282 typedef typename DGSpace::Traits::FiniteElementType::
283 Traits::LocalBasisType::Traits::RangeFieldType RF;
284 typedef typename DGSpace::Traits::FiniteElementType::
285 Traits::LocalBasisType::Traits::RangeType RangeType;
286 typedef typename DGSpace::Traits::SizeType size_type;
289 const DGSpace& dgspace_s = lfsv_s.template child<0>();
290 const DGSpace& dgspace_n = lfsv_n.template child<0>();
293 const Dune::FieldVector<DF,dim> n_F = ig.centerUnitOuterNormal();
296 const Dune::FieldVector<DF,dim>&
297 inside_local = Dune::ReferenceElements<DF,dim>::general(ig.inside()->type()).position(0,0);
298 const Dune::FieldVector<DF,dim>&
299 outside_local = Dune::ReferenceElements<DF,dim>::general(ig.outside()->type()).position(0,0);
300 RF c_s = param.c(*(ig.inside()),inside_local);
301 RF c_n = param.c(*(ig.outside()),outside_local);
304 Dune::FieldMatrix<DF,dim+1,dim+1> RT;
306 Dune::FieldVector<DF,dim+1> alpha;
307 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim-1][i];
308 Dune::FieldVector<DF,dim+1> unit(0.0);
310 Dune::FieldVector<DF,dim+1> beta;
312 Dune::FieldMatrix<DF,dim+1,dim+1> A_plus_s;
313 for (
int i=0; i<=dim; i++)
314 for (
int j=0; j<=dim; j++)
315 A_plus_s[i][j] = c_s*alpha[i]*beta[j];
319 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim][i];
323 Dune::FieldMatrix<DF,dim+1,dim+1> A_minus_n;
324 for (
int i=0; i<=dim; i++)
325 for (
int j=0; j<=dim; j++)
326 A_minus_n[i][j] = -c_n*alpha[i]*beta[j];
329 const int order_s = dgspace_s.finiteElement().localBasis().order();
330 const int order_n = dgspace_n.finiteElement().localBasis().order();
331 const int intorder = overintegration+1+2*std::max(order_s,order_n);
332 Dune::GeometryType gtface = ig.geometry().type();
333 const Dune::QuadratureRule<DF,dim-1>& rule = Dune::QuadratureRules<DF,dim-1>::rule(gtface,intorder);
338 for (
typename Dune::QuadratureRule<DF,dim-1>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
341 Dune::FieldVector<DF,dim> iplocal_s = ig.geometryInInside().global(it->position());
342 Dune::FieldVector<DF,dim> iplocal_n = ig.geometryInOutside().global(it->position());
345 const std::vector<RangeType>& phi_s = cache[order_s].evaluateFunction(iplocal_s,dgspace_s.finiteElement().localBasis());
346 const std::vector<RangeType>& phi_n = cache[order_n].evaluateFunction(iplocal_n,dgspace_n.finiteElement().localBasis());
349 Dune::FieldVector<RF,dim+1> u_s(0.0);
350 for (size_type i=0; i<=dim; i++)
351 for (size_type k=0; k<dgspace_s.size(); k++)
352 u_s[i] += x_s(lfsv_s.child(i),k)*phi_s[k];
353 Dune::FieldVector<RF,dim+1> u_n(0.0);
354 for (size_type i=0; i<=dim; i++)
355 for (size_type k=0; k<dgspace_n.size(); k++)
356 u_n[i] += x_n(lfsv_n.child(i),k)*phi_n[k];
359 Dune::FieldVector<RF,dim+1>
f(0.0);
362 A_minus_n.umv(u_n,f);
366 RF factor = it->weight() * ig.geometry().integrationElement(it->position());
367 for (size_type k=0; k<dgspace_s.size(); k++)
368 for (size_type i=0; i<=dim; i++)
369 r_s.accumulate(lfsv_s.child(i),k, f[i]*phi_s[k]*factor);
370 for (size_type k=0; k<dgspace_n.size(); k++)
371 for (size_type i=0; i<=dim; i++)
372 r_n.accumulate(lfsv_n.child(i),k, - f[i]*phi_n[k]*factor);
385 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
387 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
391 typedef typename LFSV::template Child<0>::Type DGSpace;
392 typedef typename DGSpace::Traits::FiniteElementType::
393 Traits::LocalBasisType::Traits::DomainFieldType DF;
394 typedef typename DGSpace::Traits::FiniteElementType::
395 Traits::LocalBasisType::Traits::RangeFieldType RF;
396 typedef typename DGSpace::Traits::FiniteElementType::
397 Traits::LocalBasisType::Traits::RangeType RangeType;
398 typedef typename DGSpace::Traits::SizeType size_type;
401 const DGSpace& dgspace_s = lfsv_s.template child<0>();
404 const Dune::FieldVector<DF,dim> n_F = ig.centerUnitOuterNormal();
407 const Dune::FieldVector<DF,dim>&
408 inside_local = Dune::ReferenceElements<DF,dim>::general(ig.inside()->type()).position(0,0);
409 RF c_s = param.c(*(ig.inside()),inside_local);
412 Dune::FieldMatrix<DF,dim+1,dim+1> RT;
414 Dune::FieldVector<DF,dim+1> alpha;
415 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim-1][i];
416 Dune::FieldVector<DF,dim+1> unit(0.0);
418 Dune::FieldVector<DF,dim+1> beta;
420 Dune::FieldMatrix<DF,dim+1,dim+1> A_plus_s;
421 for (
int i=0; i<=dim; i++)
422 for (
int j=0; j<=dim; j++)
423 A_plus_s[i][j] = c_s*alpha[i]*beta[j];
427 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim][i];
431 Dune::FieldMatrix<DF,dim+1,dim+1> A_minus_n;
432 for (
int i=0; i<=dim; i++)
433 for (
int j=0; j<=dim; j++)
434 A_minus_n[i][j] = -c_s*alpha[i]*beta[j];
437 const int order_s = dgspace_s.finiteElement().localBasis().order();
438 const int intorder = overintegration+1+2*order_s;
439 Dune::GeometryType gtface = ig.geometry().type();
440 const Dune::QuadratureRule<DF,dim-1>& rule = Dune::QuadratureRules<DF,dim-1>::rule(gtface,intorder);
445 for (
typename Dune::QuadratureRule<DF,dim-1>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
448 Dune::FieldVector<DF,dim> iplocal_s = ig.geometryInInside().global(it->position());
451 const std::vector<RangeType>& phi_s = cache[order_s].evaluateFunction(iplocal_s,dgspace_s.finiteElement().localBasis());
454 Dune::FieldVector<RF,dim+1> u_s(0.0);
455 for (size_type i=0; i<=dim; i++)
456 for (size_type k=0; k<dgspace_s.size(); k++)
457 u_s[i] += x_s(lfsv_s.child(i),k)*phi_s[k];
461 Dune::FieldVector<RF,dim+1> u_n(param.g(ig.intersection(),it->position(),u_s));
465 Dune::FieldVector<RF,dim+1>
f(0.0);
468 A_minus_n.umv(u_n,f);
472 RF factor = it->weight() * ig.geometry().integrationElement(it->position());
473 for (size_type k=0; k<dgspace_s.size(); k++)
474 for (size_type i=0; i<=dim; i++)
475 r_s.accumulate(lfsv_s.child(i),k, f[i]*phi_s[k]*factor);
484 template<
typename EG,
typename LFSV,
typename R>
488 typedef typename LFSV::template Child<0>::Type DGSpace;
489 typedef typename DGSpace::Traits::FiniteElementType::
490 Traits::LocalBasisType::Traits::DomainFieldType DF;
491 typedef typename DGSpace::Traits::FiniteElementType::
492 Traits::LocalBasisType::Traits::RangeFieldType RF;
493 typedef typename DGSpace::Traits::FiniteElementType::
494 Traits::LocalBasisType::Traits::RangeType RangeType;
495 typedef typename DGSpace::Traits::SizeType size_type;
498 const DGSpace& dgspace = lfsv.template child<0>();
501 const int order_s = dgspace.finiteElement().localBasis().order();
502 const int intorder = overintegration+2*order_s;
503 Dune::GeometryType gt = eg.geometry().type();
504 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
507 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
510 Dune::FieldVector<RF,dim+1> q(param.q(eg.entity(),it->position()));
513 const std::vector<RangeType>& phi = cache[order_s].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
516 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
517 for (size_type k=0; k<=dim; k++)
518 for (size_type i=0; i<dgspace.size(); i++)
519 r.accumulate(lfsv.child(k),i, - q[k]*phi[i]*factor);
524 void setTime (
typename T::Traits::RangeFieldType t)
529 void preStep (
typename T::Traits::RangeFieldType time,
typename T::Traits::RangeFieldType dt,
535 void preStage (
typename T::Traits::RangeFieldType time,
int r)
545 typename T::Traits::RangeFieldType
suggestTimestep (
typename T::Traits::RangeFieldType dt)
const
553 typedef typename FEM::Traits::FiniteElementType::Traits::LocalBasisType LocalBasisType;
555 std::vector<Cache> cache;
566 template<
typename T,
typename FEM>
572 enum { dim = T::Traits::GridViewType::dimension };
581 : param(param_), overintegration(overintegration_), cache(20)
585 template<
typename LFSU,
typename LFSV,
typename LocalPattern>
587 LocalPattern& pattern)
const
590 if (LFSU::CHILDREN!=LFSV::CHILDREN)
591 DUNE_THROW(Dune::Exception,
"need U=V!");
592 if (LFSV::CHILDREN!=dim+1)
593 DUNE_THROW(Dune::Exception,
"need exactly dim+1 components!");
595 for (
size_t k=0; k<LFSV::CHILDREN; k++)
596 for (
size_t i=0; i<lfsv.child(k).size(); ++i)
597 for (
size_t j=0; j<lfsu.child(k).size(); ++j)
598 pattern.addLink(lfsv.child(k),i,lfsu.child(k),j);
602 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
603 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
606 typedef typename LFSV::template Child<0>::Type DGSpace;
607 typedef typename DGSpace::Traits::FiniteElementType::
608 Traits::LocalBasisType::Traits::DomainFieldType DF;
609 typedef typename DGSpace::Traits::FiniteElementType::
610 Traits::LocalBasisType::Traits::RangeFieldType RF;
611 typedef typename DGSpace::Traits::FiniteElementType::
612 Traits::LocalBasisType::Traits::RangeType RangeType;
613 typedef typename DGSpace::Traits::SizeType size_type;
616 const DGSpace& dgspace = lfsv.template child<0>();
619 const int order = dgspace.finiteElement().localBasis().order();
620 const int intorder = overintegration+2*order;
621 Dune::GeometryType gt = eg.geometry().type();
622 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
625 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
628 const std::vector<RangeType>& phi = cache[order].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
631 Dune::FieldVector<RF,dim+1> u(0.0);
632 for (size_type k=0; k<=dim; k++)
633 for (size_type j=0; j<dgspace.size(); j++)
634 u[k] += x(lfsv.child(k),j)*phi[j];
637 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
638 for (size_type k=0; k<=dim; k++)
639 for (size_type i=0; i<dgspace.size(); i++)
640 r.accumulate(lfsv.child(k),i, u[k]*phi[i]*factor);
645 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
650 typedef typename LFSV::template Child<0>::Type DGSpace;
651 typedef typename DGSpace::Traits::FiniteElementType::
652 Traits::LocalBasisType::Traits::DomainFieldType DF;
653 typedef typename DGSpace::Traits::FiniteElementType::
654 Traits::LocalBasisType::Traits::RangeFieldType RF;
655 typedef typename DGSpace::Traits::FiniteElementType::
656 Traits::LocalBasisType::Traits::RangeType RangeType;
657 typedef typename DGSpace::Traits::SizeType size_type;
660 const DGSpace& dgspace = lfsv.template child<0>();
663 const int order = dgspace.finiteElement().localBasis().order();
664 const int intorder = overintegration+2*order;
665 Dune::GeometryType gt = eg.geometry().type();
666 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
669 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
672 const std::vector<RangeType>& phi = cache[order].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
675 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
676 for (size_type k=0; k<=dim; k++)
677 for (size_type i=0; i<dgspace.size(); i++)
678 for (size_type j=0; j<dgspace.size(); j++)
679 mat.accumulate(lfsv.child(k),i,lfsu.child(k),j, phi[j]*phi[i]*factor);
686 typedef typename FEM::Traits::FiniteElementType::Traits::LocalBasisType LocalBasisType;
688 std::vector<Cache> cache;
DGLinearAcousticsSpatialOperator(T ¶m_, int overintegration_=0)
Definition: linearacousticsdg.hh:187
sparsity pattern generator
Definition: pattern.hh:30
Implement jacobian_apply_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:445
static void eigenvectors_transposed(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:43
void setTime(typename T::Traits::RangeFieldType t)
set time in parameter class
Definition: linearacousticsdg.hh:524
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: linearacousticsdg.hh:485
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: linearacousticsdg.hh:194
void postStage()
to be called once at the end of each stage
Definition: linearacousticsdg.hh:540
Default flags for all local operators.
Definition: flags.hh:18
Implement jacobian_boundary() based on alpha_boundary()
Definition: defaultimp.hh:251
Definition: linearacousticsdg.hh:161
Definition: linearacousticsdg.hh:177
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538
static void eigenvectors(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:80
static const int dim
Definition: adaptivity.hh:82
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: linearacousticsdg.hh:273
Definition: linearacousticsdg.hh:183
Definition: linearacousticsdg.hh:182
Definition: linearacousticsdg.hh:181
Definition: linearacousticsdg.hh:178
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: linearacousticsdg.hh:646
Definition: linearacousticsdg.hh:578
Definition: linearacousticsdg.hh:26
const IG & ig
Definition: common/constraints.hh:146
static void eigenvectors(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:123
static void eigenvectors_transposed(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:102
store values of basis functions and gradients in a cache
Definition: localbasiscache.hh:16
Definition: linearacousticsdg.hh:575
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
const F & f
Definition: common/constraints.hh:145
Definition: linearacousticsdg.hh:184
T::Traits::RangeFieldType suggestTimestep(typename T::Traits::RangeFieldType dt) const
to be called once before each stage
Definition: linearacousticsdg.hh:545
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
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:157
void pattern_volume(const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
Definition: linearacousticsdg.hh:586
Definition: linearacousticsdg.hh:567
void preStage(typename T::Traits::RangeFieldType time, int r)
to be called once before each stage
Definition: linearacousticsdg.hh:535
DGLinearAcousticsTemporalOperator(T ¶m_, int overintegration_=0)
Definition: linearacousticsdg.hh:580
void preStep(typename T::Traits::RangeFieldType time, typename T::Traits::RangeFieldType dt, int stages)
to be called once before each time step
Definition: linearacousticsdg.hh:529
const E & e
Definition: interpolate.hh:172
static void eigenvectors(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:50
static void eigenvectors_transposed(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:72
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: linearacousticsdg.hh:386
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: linearacousticsdg.hh:603
const std::string s
Definition: function.hh:1103