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/geometry/referenceelements.hh>
41 template<
typename T1,
typename T2,
typename T3>
44 RT[0][0] = 1; RT[0][1] = c*n[0];
45 RT[1][0] = -1; RT[1][1] = c*n[0];
48 template<
typename T1,
typename T2,
typename T3>
49 static void eigenvectors (T1 c,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,dim+1,dim+1>& RT)
51 RT[0][0] = 1; RT[1][0] = c*n[0];
52 RT[0][1] = -1; RT[1][1] = c*n[0];
70 template<
typename T1,
typename T2,
typename T3>
73 RT[0][0] = 0; RT[0][1] = -n[1]; RT[0][2] = n[0];
74 RT[1][0] = 1; RT[1][1] = c*n[0]; RT[1][2] = c*n[1];
75 RT[2][0] = -1; RT[2][1] = c*n[0]; RT[2][2] = c*n[1];
78 template<
typename T1,
typename T2,
typename T3>
79 static void eigenvectors (T1 c,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,dim+1,dim+1>& RT)
81 RT[0][0] = 0; RT[1][0] = -n[1]; RT[2][0] = n[0];
82 RT[0][1] = 1; RT[1][1] = c*n[0]; RT[2][1] = c*n[1];
83 RT[0][2] = -1; RT[1][2] = c*n[0]; RT[2][2] = c*n[1];
100 template<
typename T1,
typename T2,
typename T3>
103 Dune::FieldVector<T2,dim>
s;
104 s[0] = n[1]-n[2]; s[1] = n[2]-n[0]; s[2] = n[0]-n[1];
105 if (s.two_norm()<1
e-5)
107 s[0] = n[1]+n[2]; s[1] = n[2]-n[0]; s[2] = -(n[0]+n[1]);
110 Dune::FieldVector<T2,dim> t;
111 t[0] = n[1]*s[2] - n[2]*s[1];
112 t[1] = n[2]*s[0] - n[0]*s[2];
113 t[2] = n[0]*s[1] - n[1]*s[0];
115 RT[0][0] = 0; RT[0][1] = s[0]; RT[0][2] = s[1]; RT[0][3] = s[2];
116 RT[1][0] = 0; RT[1][1] = t[0]; RT[1][2] = t[1]; RT[1][3] = t[2];
117 RT[2][0] = 1; RT[2][1] = c*n[0]; RT[2][2] = c*n[1]; RT[2][3] = c*n[2];
118 RT[3][0] = -1; RT[3][1] = c*n[0]; RT[3][2] = c*n[1]; RT[3][3] = c*n[2];
121 template<
typename T1,
typename T2,
typename T3>
122 static void eigenvectors (T1 c,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,dim+1,dim+1>& RT)
124 Dune::FieldVector<T2,dim>
s;
125 s[0] = n[1]-n[2]; s[1] = n[2]-n[0]; s[2] = n[0]-n[1];
126 if (s.two_norm()<1
e-5)
128 s[0] = n[1]+n[2]; s[1] = n[2]-n[0]; s[2] = -(n[0]+n[1]);
131 Dune::FieldVector<T2,dim> t;
132 t[0] = n[1]*s[2] - n[2]*s[1];
133 t[1] = n[2]*s[0] - n[0]*s[2];
134 t[2] = n[0]*s[1] - n[1]*s[0];
136 RT[0][0] = 0; RT[1][0] = s[0]; RT[2][0] = s[1]; RT[3][0] = s[2];
137 RT[0][1] = 0; RT[1][1] = t[0]; RT[2][1] = t[1]; RT[3][1] = t[2];
138 RT[0][2] = 1; RT[1][2] = c*n[0]; RT[2][2] = c*n[1]; RT[3][2] = c*n[2];
139 RT[0][3] = -1; RT[1][3] = c*n[0]; RT[2][3] = c*n[1]; RT[3][3] = c*n[2];
159 template<
typename T,
typename FEM>
172 enum { dim = T::Traits::GridViewType::dimension };
187 : param(param_), overintegration(overintegration_), cache(20)
192 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
193 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
196 typedef typename LFSV::template Child<0>::Type DGSpace;
197 typedef typename DGSpace::Traits::FiniteElementType::
198 Traits::LocalBasisType::Traits::DomainFieldType DF;
199 typedef typename DGSpace::Traits::FiniteElementType::
200 Traits::LocalBasisType::Traits::RangeFieldType RF;
201 typedef typename DGSpace::Traits::FiniteElementType::
202 Traits::LocalBasisType::Traits::RangeType RangeType;
203 typedef typename DGSpace::Traits::FiniteElementType::
204 Traits::LocalBasisType::Traits::JacobianType JacobianType;
205 typedef typename DGSpace::Traits::SizeType size_type;
208 if (LFSV::CHILDREN!=dim+1)
209 DUNE_THROW(Dune::Exception,
"need exactly dim+1 components!");
212 const DGSpace& dgspace = lfsv.template child<0>();
214 auto geometry = eg.geometry();
217 const int order = dgspace.finiteElement().localBasis().order();
218 const int intorder = overintegration+2*order;
219 Dune::GeometryType gt = geometry.type();
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 (
const auto& ip : Dune::QuadratureRules<DF,dim>::rule(gt,intorder))
235 const std::vector<RangeType>& phi = cache[order].evaluateFunction(ip.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(ip.position(),dgspace.finiteElement().localBasis());
248 jac = geometry.jacobianInverseTransposed(ip.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 = ip.weight() * geometry.integrationElement(ip.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 auto inside_cell = ig.inside();
297 auto outside_cell = ig.outside();
300 const Dune::FieldVector<DF,dim>&
301 inside_local = Dune::ReferenceElements<DF,dim>::general(inside_cell.type()).position(0,0);
302 const Dune::FieldVector<DF,dim>&
303 outside_local = Dune::ReferenceElements<DF,dim>::general(outside_cell.type()).position(0,0);
304 RF c_s = param.c(inside_cell,inside_local);
305 RF c_n = param.c(outside_cell,outside_local);
308 Dune::FieldMatrix<DF,dim+1,dim+1> RT;
310 Dune::FieldVector<DF,dim+1> alpha;
311 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim-1][i];
312 Dune::FieldVector<DF,dim+1> unit(0.0);
314 Dune::FieldVector<DF,dim+1> beta;
316 Dune::FieldMatrix<DF,dim+1,dim+1> A_plus_s;
317 for (
int i=0; i<=dim; i++)
318 for (
int j=0; j<=dim; j++)
319 A_plus_s[i][j] = c_s*alpha[i]*beta[j];
323 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim][i];
327 Dune::FieldMatrix<DF,dim+1,dim+1> A_minus_n;
328 for (
int i=0; i<=dim; i++)
329 for (
int j=0; j<=dim; j++)
330 A_minus_n[i][j] = -c_n*alpha[i]*beta[j];
332 auto geometry = ig.geometry();
335 const int order_s = dgspace_s.finiteElement().localBasis().order();
336 const int order_n = dgspace_n.finiteElement().localBasis().order();
337 const int intorder = overintegration+1+2*std::max(order_s,order_n);
338 Dune::GeometryType gtface = geometry.type();
343 for (
const auto& ip : Dune::QuadratureRules<DF,dim-1>::rule(gtface,intorder))
346 Dune::FieldVector<DF,dim> iplocal_s = ig.geometryInInside().global(ip.position());
347 Dune::FieldVector<DF,dim> iplocal_n = ig.geometryInOutside().global(ip.position());
350 const std::vector<RangeType>& phi_s = cache[order_s].evaluateFunction(iplocal_s,dgspace_s.finiteElement().localBasis());
351 const std::vector<RangeType>& phi_n = cache[order_n].evaluateFunction(iplocal_n,dgspace_n.finiteElement().localBasis());
354 Dune::FieldVector<RF,dim+1> u_s(0.0);
355 for (size_type i=0; i<=dim; i++)
356 for (size_type k=0; k<dgspace_s.size(); k++)
357 u_s[i] += x_s(lfsv_s.child(i),k)*phi_s[k];
358 Dune::FieldVector<RF,dim+1> u_n(0.0);
359 for (size_type i=0; i<=dim; i++)
360 for (size_type k=0; k<dgspace_n.size(); k++)
361 u_n[i] += x_n(lfsv_n.child(i),k)*phi_n[k];
364 Dune::FieldVector<RF,dim+1> f(0.0);
367 A_minus_n.umv(u_n,f);
371 RF factor = ip.weight() * geometry.integrationElement(ip.position());
372 for (size_type k=0; k<dgspace_s.size(); k++)
373 for (size_type i=0; i<=dim; i++)
374 r_s.accumulate(lfsv_s.child(i),k, f[i]*phi_s[k]*factor);
375 for (size_type k=0; k<dgspace_n.size(); k++)
376 for (size_type i=0; i<=dim; i++)
377 r_n.accumulate(lfsv_n.child(i),k, - f[i]*phi_n[k]*factor);
390 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
392 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
396 typedef typename LFSV::template Child<0>::Type DGSpace;
397 typedef typename DGSpace::Traits::FiniteElementType::
398 Traits::LocalBasisType::Traits::DomainFieldType DF;
399 typedef typename DGSpace::Traits::FiniteElementType::
400 Traits::LocalBasisType::Traits::RangeFieldType RF;
401 typedef typename DGSpace::Traits::FiniteElementType::
402 Traits::LocalBasisType::Traits::RangeType RangeType;
403 typedef typename DGSpace::Traits::SizeType size_type;
406 const DGSpace& dgspace_s = lfsv_s.template child<0>();
409 const Dune::FieldVector<DF,dim> n_F = ig.centerUnitOuterNormal();
411 auto cell = ig.inside();
414 const Dune::FieldVector<DF,dim>&
415 inside_local = Dune::ReferenceElements<DF,dim>::general(cell.type()).position(0,0);
416 RF c_s = param.c(cell,inside_local);
419 Dune::FieldMatrix<DF,dim+1,dim+1> RT;
421 Dune::FieldVector<DF,dim+1> alpha;
422 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim-1][i];
423 Dune::FieldVector<DF,dim+1> unit(0.0);
425 Dune::FieldVector<DF,dim+1> beta;
427 Dune::FieldMatrix<DF,dim+1,dim+1> A_plus_s;
428 for (
int i=0; i<=dim; i++)
429 for (
int j=0; j<=dim; j++)
430 A_plus_s[i][j] = c_s*alpha[i]*beta[j];
434 for (
int i=0; i<=dim; i++) alpha[i] = RT[dim][i];
438 Dune::FieldMatrix<DF,dim+1,dim+1> A_minus_n;
439 for (
int i=0; i<=dim; i++)
440 for (
int j=0; j<=dim; j++)
441 A_minus_n[i][j] = -c_s*alpha[i]*beta[j];
443 auto geometry = ig.geometry();
446 const int order_s = dgspace_s.finiteElement().localBasis().order();
447 const int intorder = overintegration+1+2*order_s;
448 Dune::GeometryType gtface = geometry.type();
453 for (
const auto& ip : Dune::QuadratureRules<DF,dim-1>::rule(gtface,intorder))
456 Dune::FieldVector<DF,dim> iplocal_s = ig.geometryInInside().global(ip.position());
459 const std::vector<RangeType>& phi_s = cache[order_s].evaluateFunction(iplocal_s,dgspace_s.finiteElement().localBasis());
462 Dune::FieldVector<RF,dim+1> u_s(0.0);
463 for (size_type i=0; i<=dim; i++)
464 for (size_type k=0; k<dgspace_s.size(); k++)
465 u_s[i] += x_s(lfsv_s.child(i),k)*phi_s[k];
469 Dune::FieldVector<RF,dim+1> u_n(param.g(ig.intersection(),ip.position(),u_s));
473 Dune::FieldVector<RF,dim+1> f(0.0);
476 A_minus_n.umv(u_n,f);
480 RF factor = ip.weight() * geometry.integrationElement(ip.position());
481 for (size_type k=0; k<dgspace_s.size(); k++)
482 for (size_type i=0; i<=dim; i++)
483 r_s.accumulate(lfsv_s.child(i),k, f[i]*phi_s[k]*factor);
492 template<
typename EG,
typename LFSV,
typename R>
496 typedef typename LFSV::template Child<0>::Type DGSpace;
497 typedef typename DGSpace::Traits::FiniteElementType::
498 Traits::LocalBasisType::Traits::DomainFieldType DF;
499 typedef typename DGSpace::Traits::FiniteElementType::
500 Traits::LocalBasisType::Traits::RangeFieldType RF;
501 typedef typename DGSpace::Traits::FiniteElementType::
502 Traits::LocalBasisType::Traits::RangeType RangeType;
503 typedef typename DGSpace::Traits::SizeType size_type;
506 const DGSpace& dgspace = lfsv.template child<0>();
508 auto geometry = eg.geometry();
511 const int order_s = dgspace.finiteElement().localBasis().order();
512 const int intorder = overintegration+2*order_s;
513 Dune::GeometryType gt = geometry.type();
516 for (
const auto& ip : Dune::QuadratureRules<DF,dim>::rule(gt,intorder))
519 Dune::FieldVector<RF,dim+1> q(param.q(eg.entity(),ip.position()));
522 const std::vector<RangeType>& phi = cache[order_s].evaluateFunction(ip.position(),dgspace.finiteElement().localBasis());
525 RF factor = ip.weight() * geometry.integrationElement(ip.position());
526 for (size_type k=0; k<=dim; k++)
527 for (size_type i=0; i<dgspace.size(); i++)
528 r.accumulate(lfsv.child(k),i, - q[k]*phi[i]*factor);
533 void setTime (
typename T::Traits::RangeFieldType t)
538 void preStep (
typename T::Traits::RangeFieldType time,
typename T::Traits::RangeFieldType dt,
544 void preStage (
typename T::Traits::RangeFieldType time,
int r)
554 typename T::Traits::RangeFieldType
suggestTimestep (
typename T::Traits::RangeFieldType dt)
const
562 typedef typename FEM::Traits::FiniteElementType::Traits::LocalBasisType LocalBasisType;
564 std::vector<Cache> cache;
575 template<
typename T,
typename FEM>
581 enum { dim = T::Traits::GridViewType::dimension };
590 : param(param_), overintegration(overintegration_), cache(20)
594 template<
typename LFSU,
typename LFSV,
typename LocalPattern>
596 LocalPattern& pattern)
const
599 if (LFSU::CHILDREN!=LFSV::CHILDREN)
600 DUNE_THROW(Dune::Exception,
"need U=V!");
601 if (LFSV::CHILDREN!=dim+1)
602 DUNE_THROW(Dune::Exception,
"need exactly dim+1 components!");
604 for (
size_t k=0; k<LFSV::CHILDREN; k++)
605 for (
size_t i=0; i<lfsv.child(k).size(); ++i)
606 for (
size_t j=0; j<lfsu.child(k).size(); ++j)
607 pattern.addLink(lfsv.child(k),i,lfsu.child(k),j);
611 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
612 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
615 typedef typename LFSV::template Child<0>::Type DGSpace;
616 typedef typename DGSpace::Traits::FiniteElementType::
617 Traits::LocalBasisType::Traits::DomainFieldType DF;
618 typedef typename DGSpace::Traits::FiniteElementType::
619 Traits::LocalBasisType::Traits::RangeFieldType RF;
620 typedef typename DGSpace::Traits::FiniteElementType::
621 Traits::LocalBasisType::Traits::RangeType RangeType;
622 typedef typename DGSpace::Traits::SizeType size_type;
625 const DGSpace& dgspace = lfsv.template child<0>();
627 auto geometry = eg.geometry();
630 const int order = dgspace.finiteElement().localBasis().order();
631 const int intorder = overintegration+2*order;
632 Dune::GeometryType gt = geometry.type();
635 for (
const auto& ip : Dune::QuadratureRules<DF,dim>::rule(gt,intorder))
638 const std::vector<RangeType>& phi = cache[order].evaluateFunction(ip.position(),dgspace.finiteElement().localBasis());
641 Dune::FieldVector<RF,dim+1> u(0.0);
642 for (size_type k=0; k<=dim; k++)
643 for (size_type j=0; j<dgspace.size(); j++)
644 u[k] += x(lfsv.child(k),j)*phi[j];
647 RF factor = ip.weight() * geometry.integrationElement(ip.position());
648 for (size_type k=0; k<=dim; k++)
649 for (size_type i=0; i<dgspace.size(); i++)
650 r.accumulate(lfsv.child(k),i, u[k]*phi[i]*factor);
655 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
660 typedef typename LFSV::template Child<0>::Type DGSpace;
661 typedef typename DGSpace::Traits::FiniteElementType::
662 Traits::LocalBasisType::Traits::DomainFieldType DF;
663 typedef typename DGSpace::Traits::FiniteElementType::
664 Traits::LocalBasisType::Traits::RangeFieldType RF;
665 typedef typename DGSpace::Traits::FiniteElementType::
666 Traits::LocalBasisType::Traits::RangeType RangeType;
667 typedef typename DGSpace::Traits::SizeType size_type;
670 const DGSpace& dgspace = lfsv.template child<0>();
672 auto geometry = eg.geometry();
675 const int order = dgspace.finiteElement().localBasis().order();
676 const int intorder = overintegration+2*order;
677 Dune::GeometryType gt = geometry.type();
680 for (
const auto& ip : Dune::QuadratureRules<DF,dim>::rule(gt,intorder))
683 const std::vector<RangeType>& phi = cache[order].evaluateFunction(ip.position(),dgspace.finiteElement().localBasis());
686 RF factor = ip.weight() * geometry.integrationElement(ip.position());
687 for (size_type k=0; k<=dim; k++)
688 for (size_type i=0; i<dgspace.size(); i++)
689 for (size_type j=0; j<dgspace.size(); j++)
690 mat.accumulate(lfsv.child(k),i,lfsu.child(k),j, phi[j]*phi[i]*factor);
697 typedef typename FEM::Traits::FiniteElementType::Traits::LocalBasisType LocalBasisType;
699 std::vector<Cache> cache;
DGLinearAcousticsSpatialOperator(T ¶m_, int overintegration_=0)
Definition: linearacousticsdg.hh:186
static void eigenvectors_transposed(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:42
void setTime(typename T::Traits::RangeFieldType t)
set time in parameter class
Definition: linearacousticsdg.hh:533
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: linearacousticsdg.hh:493
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: linearacousticsdg.hh:193
void postStage()
to be called once at the end of each stage
Definition: linearacousticsdg.hh:549
Definition: linearacousticsdg.hh:160
static void eigenvectors_transposed(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:101
const E & e
Definition: interpolate.hh:172
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
Definition: linearacousticsdg.hh:177
Definition: linearacousticsdg.hh:587
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
const IG & ig
Definition: constraints.hh:147
Definition: linearacousticsdg.hh:183
Definition: linearacousticsdg.hh:180
static void eigenvectors(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:79
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:13
static void eigenvectors(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:122
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: linearacousticsdg.hh:656
static void eigenvectors_transposed(T1 c, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, dim+1, dim+1 > &RT)
Definition: linearacousticsdg.hh:71
Definition: linearacousticsdg.hh:182
Definition: linearacousticsdg.hh:584
Implement jacobian_apply_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:445
Definition: adaptivity.hh:27
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
Implement jacobian_boundary() based on alpha_boundary()
Definition: defaultimp.hh:251
static const int dim
Definition: adaptivity.hh:83
T::Traits::RangeFieldType suggestTimestep(typename T::Traits::RangeFieldType dt) const
to be called once before each stage
Definition: linearacousticsdg.hh:554
Definition: linearacousticsdg.hh:176
const std::string s
Definition: function.hh:1102
Definition: linearacousticsdg.hh:25
sparsity pattern generator
Definition: pattern.hh:29
void pattern_volume(const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
Definition: linearacousticsdg.hh:595
Definition: linearacousticsdg.hh:576
void preStage(typename T::Traits::RangeFieldType time, int r)
to be called once before each stage
Definition: linearacousticsdg.hh:544
Definition: linearacousticsdg.hh:181
DGLinearAcousticsTemporalOperator(T ¶m_, int overintegration_=0)
Definition: linearacousticsdg.hh:589
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:538
store values of basis functions and gradients in a cache
Definition: localbasiscache.hh:15
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:157
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:391
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: linearacousticsdg.hh:612
const EG & eg
Definition: constraints.hh:280
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
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:49