2 #ifndef DUNE_PDELAB_MAXWELLDG_HH
3 #define DUNE_PDELAB_MAXWELLDG_HH
7 #include<dune/common/exceptions.hh>
8 #include<dune/common/fvector.hh>
9 #include<dune/common/static_assert.hh>
11 #include <dune/geometry/referenceelements.hh>
49 template<
typename T1,
typename T2,
typename T3>
50 static void eigenvalues (T1 eps, T1 mu,
const Dune::FieldVector<T2,2*dim>&
e)
52 T1
s = 1.0/sqrt(mu*eps);
72 template<
typename T1,
typename T2,
typename T3>
73 static void eigenvectors (T1 eps, T1 mu,
const Dune::FieldVector<T2,dim>& n, Dune::FieldMatrix<T3,2*dim,2*dim>& R)
75 T1 a=n[0], b=n[1], c=n[2];
77 Dune::FieldVector<T2,dim> re, im;
80 re[0]=a*c; re[1]=b*c; re[2]=c*c-1;
81 im[0]=-b; im[1]=a; im[2]=0;
85 re[0]=a*b; re[1]=b*b-1; re[2]=b*c;
86 im[0]=c; im[1]=0.0; im[2]=-a;
90 R[0][0] = re[0]; R[0][1] = -im[0];
91 R[1][0] = re[1]; R[1][1] = -im[1];
92 R[2][0] = re[2]; R[2][1] = -im[2];
93 R[3][0] = im[0]; R[3][1] = re[0];
94 R[4][0] = im[1]; R[4][1] = re[1];
95 R[5][0] = im[2]; R[5][1] = re[2];
98 R[0][2] = im[0]; R[0][3] = re[0];
99 R[1][2] = im[1]; R[1][3] = re[1];
100 R[2][2] = im[2]; R[2][3] = re[2];
101 R[3][2] = re[0]; R[3][3] = -im[0];
102 R[4][2] = re[1]; R[4][3] = -im[1];
103 R[5][2] = re[2]; R[5][3] = -im[2];
106 R[0][4] = a; R[0][5] = 0;
107 R[1][4] = b; R[1][5] = 0;
108 R[2][4] = c; R[2][5] = 0;
109 R[3][4] = 0; R[3][5] = a;
110 R[4][4] = 0; R[4][5] = b;
111 R[5][4] = 0; R[5][5] = c;
116 for (std::size_t i=0; i<3; i++)
117 for (std::size_t j=0; j<6; j++)
119 for (std::size_t i=3; i<6; i++)
120 for (std::size_t j=0; j<6; j++)
300 template<
typename T,
typename FEM>
313 enum { dim = T::Traits::GridViewType::dimension };
328 : param(param_), overintegration(overintegration_), cache(20)
333 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
334 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
337 typedef typename LFSV::template Child<0>::Type DGSpace;
338 typedef typename DGSpace::Traits::FiniteElementType::
339 Traits::LocalBasisType::Traits::DomainFieldType DF;
340 typedef typename DGSpace::Traits::FiniteElementType::
341 Traits::LocalBasisType::Traits::RangeFieldType RF;
342 typedef typename DGSpace::Traits::FiniteElementType::
343 Traits::LocalBasisType::Traits::RangeType RangeType;
344 typedef typename DGSpace::Traits::FiniteElementType::
345 Traits::LocalBasisType::Traits::JacobianType JacobianType;
346 typedef typename DGSpace::Traits::SizeType size_type;
349 if (LFSV::CHILDREN!=dim*2)
350 DUNE_THROW(Dune::Exception,
"need exactly dim*2 components!");
353 const DGSpace& dgspace = lfsv.template child<0>();
356 const int order = dgspace.finiteElement().localBasis().order();
357 const int intorder = overintegration+2*order;
358 Dune::GeometryType gt = eg.geometry().type();
359 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
362 typename EG::Geometry::JacobianInverseTransposed jac;
365 Dune::FieldVector<DF,dim> localcenter = Dune::ReferenceElements<DF,dim>::general(gt).position(0,0);
366 RF mu = param.mu(eg.entity(),localcenter);
367 RF eps = param.eps(eg.entity(),localcenter);
368 RF sigma = param.sigma(eg.entity(),localcenter);
375 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
378 const std::vector<RangeType>& phi = cache[order].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
381 Dune::FieldVector<RF,dim*2> u(0.0);
382 for (size_type k=0; k<dim*2; k++)
383 for (size_type j=0; j<dgspace.size(); j++)
384 u[k] += x(lfsv.child(k),j)*phi[j];
388 const std::vector<JacobianType>& js = cache[order].evaluateJacobian(it->position(),dgspace.finiteElement().localBasis());
391 jac = eg.geometry().jacobianInverseTransposed(it->position());
392 std::vector<Dune::FieldVector<RF,dim> > gradphi(dgspace.size());
393 for (size_type i=0; i<dgspace.size(); i++)
394 jac.mv(js[i][0],gradphi[i]);
397 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
399 Dune::FieldMatrix<RF,dim*2,dim> F;
400 F[0][0] = 0; F[0][1] = -muinv*u[5]; F[0][2] = muinv*u[4];
401 F[1][0] = muinv*u[5]; F[1][1] = 0; F[1][2] = -muinv*u[3];
402 F[2][0] =-muinv*u[4]; F[2][1] = muinv*u[3]; F[2][2] = 0;
403 F[3][0] = 0; F[3][1] = epsinv*u[2]; F[3][2] = -epsinv*u[1];
404 F[4][0] = -epsinv*u[2]; F[4][1] = 0; F[4][2] = epsinv*u[0];
405 F[5][0] = epsinv*u[1]; F[5][1] = -epsinv*u[0]; F[5][2] = 0;
408 for (size_type i=0; i<dim*2; i++)
410 for (size_type k=0; k<dgspace.size(); k++)
412 for (size_type j=0; j<dim; j++)
413 r.accumulate(lfsv.child(i),k,-F[i][j]*gradphi[k][j]*factor);
416 for (size_type i=0; i<dim; i++)
418 for (size_type k=0; k<dgspace.size(); k++)
419 r.accumulate(lfsv.child(i),k,(sigma/eps)*u[i]*phi[k]*factor);
429 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
431 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
432 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
433 R& r_s, R& r_n)
const
436 typedef typename LFSV::template Child<0>::Type DGSpace;
437 typedef typename DGSpace::Traits::FiniteElementType::
438 Traits::LocalBasisType::Traits::DomainFieldType DF;
439 typedef typename DGSpace::Traits::FiniteElementType::
440 Traits::LocalBasisType::Traits::RangeFieldType RF;
441 typedef typename DGSpace::Traits::FiniteElementType::
442 Traits::LocalBasisType::Traits::RangeType RangeType;
443 typedef typename DGSpace::Traits::SizeType size_type;
446 const DGSpace& dgspace_s = lfsv_s.template child<0>();
447 const DGSpace& dgspace_n = lfsv_n.template child<0>();
450 const Dune::FieldVector<DF,dim> n_F = ig.centerUnitOuterNormal();
453 const Dune::FieldVector<DF,dim>&
454 inside_local = Dune::ReferenceElements<DF,dim>::general(ig.inside()->type()).position(0,0);
455 const Dune::FieldVector<DF,dim>&
456 outside_local = Dune::ReferenceElements<DF,dim>::general(ig.outside()->type()).position(0,0);
457 RF mu_s = param.mu(*(ig.inside()),inside_local);
458 RF mu_n = param.mu(*(ig.outside()),outside_local);
459 RF eps_s = param.eps(*(ig.inside()),inside_local);
460 RF eps_n = param.eps(*(ig.outside()),outside_local);
465 Dune::FieldMatrix<DF,dim*2,dim*2> R_s;
467 Dune::FieldMatrix<DF,dim*2,dim*2> Dplus_s(0.0);
468 Dplus_s[0][0] = 1.0/sqrt(eps_s*mu_s);
469 Dplus_s[1][1] = 1.0/sqrt(eps_s*mu_s);
470 Dune::FieldMatrix<DF,dim*2,dim*2> Aplus_s(R_s);
471 Aplus_s.rightmultiply(Dplus_s);
473 Aplus_s.rightmultiply(R_s);
476 Dune::FieldMatrix<DF,dim*2,dim*2> R_n;
478 Dune::FieldMatrix<DF,dim*2,dim*2> Dminus_n(0.0);
479 Dminus_n[2][2] = -1.0/sqrt(eps_n*mu_n);
480 Dminus_n[3][3] = -1.0/sqrt(eps_n*mu_n);
481 Dune::FieldMatrix<DF,dim*2,dim*2> Aminus_n(R_n);
482 Aminus_n.rightmultiply(Dminus_n);
484 Aminus_n.rightmultiply(R_n);
487 const int order_s = dgspace_s.finiteElement().localBasis().order();
488 const int order_n = dgspace_n.finiteElement().localBasis().order();
489 const int intorder = overintegration+1+2*std::max(order_s,order_n);
490 Dune::GeometryType gtface = ig.geometry().type();
491 const Dune::QuadratureRule<DF,dim-1>& rule = Dune::QuadratureRules<DF,dim-1>::rule(gtface,intorder);
496 for (
typename Dune::QuadratureRule<DF,dim-1>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
499 Dune::FieldVector<DF,dim> iplocal_s = ig.geometryInInside().global(it->position());
500 Dune::FieldVector<DF,dim> iplocal_n = ig.geometryInOutside().global(it->position());
503 const std::vector<RangeType>& phi_s = cache[order_s].evaluateFunction(iplocal_s,dgspace_s.finiteElement().localBasis());
504 const std::vector<RangeType>& phi_n = cache[order_n].evaluateFunction(iplocal_n,dgspace_n.finiteElement().localBasis());
507 Dune::FieldVector<RF,dim*2> u_s(0.0);
508 for (size_type i=0; i<dim*2; i++)
509 for (size_type k=0; k<dgspace_s.size(); k++)
510 u_s[i] += x_s(lfsv_s.child(i),k)*phi_s[k];
511 Dune::FieldVector<RF,dim*2> u_n(0.0);
512 for (size_type i=0; i<dim*2; i++)
513 for (size_type k=0; k<dgspace_n.size(); k++)
514 u_n[i] += x_n(lfsv_n.child(i),k)*phi_n[k];
517 Dune::FieldVector<RF,dim*2>
f(0.0);
526 RF factor = it->weight() * ig.geometry().integrationElement(it->position());
527 for (size_type k=0; k<dgspace_s.size(); k++)
528 for (size_type i=0; i<dim*2; i++)
529 r_s.accumulate(lfsv_s.child(i),k,f[i]*phi_s[k]*factor);
530 for (size_type k=0; k<dgspace_n.size(); k++)
531 for (size_type i=0; i<dim*2; i++)
532 r_n.accumulate(lfsv_n.child(i),k,-f[i]*phi_n[k]*factor);
544 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
546 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
550 typedef typename LFSV::template Child<0>::Type DGSpace;
551 typedef typename DGSpace::Traits::FiniteElementType::
552 Traits::LocalBasisType::Traits::DomainFieldType DF;
553 typedef typename DGSpace::Traits::FiniteElementType::
554 Traits::LocalBasisType::Traits::RangeFieldType RF;
555 typedef typename DGSpace::Traits::FiniteElementType::
556 Traits::LocalBasisType::Traits::RangeType RangeType;
557 typedef typename DGSpace::Traits::SizeType size_type;
560 const DGSpace& dgspace_s = lfsv_s.template child<0>();
563 const Dune::FieldVector<DF,dim> n_F = ig.centerUnitOuterNormal();
566 const Dune::FieldVector<DF,dim>&
567 inside_local = Dune::ReferenceElements<DF,dim>::general(ig.inside()->type()).position(0,0);
568 RF mu_s = param.mu(*(ig.inside()),inside_local);
569 RF eps_s = param.eps(*(ig.inside()),inside_local);
573 Dune::FieldMatrix<DF,dim*2,dim*2> R_s;
575 Dune::FieldMatrix<DF,dim*2,dim*2> Dplus_s(0.0);
576 Dplus_s[0][0] = 1.0/sqrt(eps_s*mu_s);
577 Dplus_s[1][1] = 1.0/sqrt(eps_s*mu_s);
578 Dune::FieldMatrix<DF,dim*2,dim*2> Aplus_s(R_s);
579 Aplus_s.rightmultiply(Dplus_s);
581 Aplus_s.rightmultiply(R_s);
584 Dune::FieldMatrix<DF,dim*2,dim*2> R_n;
586 Dune::FieldMatrix<DF,dim*2,dim*2> Dminus_n(0.0);
587 Dminus_n[2][2] = -1.0/sqrt(eps_s*mu_s);
588 Dminus_n[3][3] = -1.0/sqrt(eps_s*mu_s);
589 Dune::FieldMatrix<DF,dim*2,dim*2> Aminus_n(R_n);
590 Aminus_n.rightmultiply(Dminus_n);
592 Aminus_n.rightmultiply(R_n);
595 const int order_s = dgspace_s.finiteElement().localBasis().order();
596 const int intorder = overintegration+1+2*order_s;
597 Dune::GeometryType gtface = ig.geometry().type();
598 const Dune::QuadratureRule<DF,dim-1>& rule = Dune::QuadratureRules<DF,dim-1>::rule(gtface,intorder);
603 for (
typename Dune::QuadratureRule<DF,dim-1>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
606 Dune::FieldVector<DF,dim> iplocal_s = ig.geometryInInside().global(it->position());
609 const std::vector<RangeType>& phi_s = cache[order_s].evaluateFunction(iplocal_s,dgspace_s.finiteElement().localBasis());
612 Dune::FieldVector<RF,dim*2> u_s(0.0);
613 for (size_type i=0; i<dim*2; i++)
614 for (size_type k=0; k<dgspace_s.size(); k++)
615 u_s[i] += x_s(lfsv_s.child(i),k)*phi_s[k];
619 Dune::FieldVector<RF,dim*2> u_n(param.g(ig.intersection(),it->position(),u_s));
623 Dune::FieldVector<RF,dim*2>
f(0.0);
630 RF factor = it->weight() * ig.geometry().integrationElement(it->position());
631 for (size_type k=0; k<dgspace_s.size(); k++)
632 for (size_type i=0; i<dim*2; i++)
633 r_s.accumulate(lfsv_s.child(i),k,f[i]*phi_s[k]*factor);
642 template<
typename EG,
typename LFSV,
typename R>
646 typedef typename LFSV::template Child<0>::Type DGSpace;
647 typedef typename DGSpace::Traits::FiniteElementType::
648 Traits::LocalBasisType::Traits::DomainFieldType DF;
649 typedef typename DGSpace::Traits::FiniteElementType::
650 Traits::LocalBasisType::Traits::RangeFieldType RF;
651 typedef typename DGSpace::Traits::FiniteElementType::
652 Traits::LocalBasisType::Traits::RangeType RangeType;
653 typedef typename DGSpace::Traits::SizeType size_type;
656 const DGSpace& dgspace = lfsv.template child<0>();
659 const int order_s = dgspace.finiteElement().localBasis().order();
660 const int intorder = overintegration+2*order_s;
661 Dune::GeometryType gt = eg.geometry().type();
662 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
665 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
668 Dune::FieldVector<RF,dim*2> j(param.j(eg.entity(),it->position()));
671 const std::vector<RangeType>& phi = cache[order_s].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
674 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
675 for (size_type k=0; k<dim*2; k++)
676 for (size_type i=0; i<dgspace.size(); i++)
677 r.accumulate(lfsv.child(k),i,-j[k]*phi[i]*factor);
682 void setTime (
typename T::Traits::RangeFieldType t)
687 void preStep (
typename T::Traits::RangeFieldType time,
typename T::Traits::RangeFieldType dt,
693 void preStage (
typename T::Traits::RangeFieldType time,
int r)
703 typename T::Traits::RangeFieldType
suggestTimestep (
typename T::Traits::RangeFieldType dt)
const
711 typedef typename FEM::Traits::FiniteElementType::Traits::LocalBasisType LocalBasisType;
713 std::vector<Cache> cache;
724 template<
typename T,
typename FEM>
730 enum { dim = T::Traits::GridViewType::dimension };
739 : param(param_), overintegration(overintegration_), cache(20)
743 template<
typename LFSU,
typename LFSV,
typename LocalPattern>
745 LocalPattern& pattern)
const
748 if (LFSU::CHILDREN!=LFSV::CHILDREN)
749 DUNE_THROW(Dune::Exception,
"need U=V!");
750 if (LFSV::CHILDREN!=dim*2)
751 DUNE_THROW(Dune::Exception,
"need exactly dim*2 components!");
753 for (
size_t k=0; k<LFSV::CHILDREN; k++)
754 for (
size_t i=0; i<lfsv.child(k).size(); ++i)
755 for (
size_t j=0; j<lfsu.child(k).size(); ++j)
756 pattern.addLink(lfsv.child(k),i,lfsu.child(k),j);
760 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
761 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
764 typedef typename LFSV::template Child<0>::Type DGSpace;
765 typedef typename DGSpace::Traits::FiniteElementType::
766 Traits::LocalBasisType::Traits::DomainFieldType DF;
767 typedef typename DGSpace::Traits::FiniteElementType::
768 Traits::LocalBasisType::Traits::RangeFieldType RF;
769 typedef typename DGSpace::Traits::FiniteElementType::
770 Traits::LocalBasisType::Traits::RangeType RangeType;
771 typedef typename DGSpace::Traits::SizeType size_type;
774 const DGSpace& dgspace = lfsv.template child<0>();
777 const int order = dgspace.finiteElement().localBasis().order();
778 const int intorder = overintegration+2*order;
779 Dune::GeometryType gt = eg.geometry().type();
780 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
783 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
786 const std::vector<RangeType>& phi = cache[order].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
789 Dune::FieldVector<RF,dim*2> u(0.0);
790 for (size_type k=0; k<dim*2; k++)
791 for (size_type j=0; j<dgspace.size(); j++)
792 u[k] += x(lfsv.child(k),j)*phi[j];
795 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
796 for (size_type k=0; k<dim*2; k++)
797 for (size_type i=0; i<dgspace.size(); i++)
798 r.accumulate(lfsv.child(k),i,u[k]*phi[i]*factor);
803 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
808 typedef typename LFSV::template Child<0>::Type DGSpace;
809 typedef typename DGSpace::Traits::FiniteElementType::
810 Traits::LocalBasisType::Traits::DomainFieldType DF;
811 typedef typename DGSpace::Traits::FiniteElementType::
812 Traits::LocalBasisType::Traits::RangeFieldType RF;
813 typedef typename DGSpace::Traits::FiniteElementType::
814 Traits::LocalBasisType::Traits::RangeType RangeType;
815 typedef typename DGSpace::Traits::SizeType size_type;
818 const DGSpace& dgspace = lfsv.template child<0>();
821 const int order = dgspace.finiteElement().localBasis().order();
822 const int intorder = overintegration+2*order;
823 Dune::GeometryType gt = eg.geometry().type();
824 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
827 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
830 const std::vector<RangeType>& phi = cache[order].evaluateFunction(it->position(),dgspace.finiteElement().localBasis());
833 RF factor = it->weight() * eg.geometry().integrationElement(it->position());
834 for (size_type k=0; k<dim*2; k++)
835 for (size_type i=0; i<dgspace.size(); i++)
836 for (size_type j=0; j<dgspace.size(); j++)
837 mat.accumulate(lfsv.child(k),i,lfsu.child(k),j,phi[j]*phi[i]*factor);
844 typedef typename FEM::Traits::FiniteElementType::Traits::LocalBasisType LocalBasisType;
846 std::vector<Cache> cache;
Definition: maxwelldg.hh:324
Definition: maxwelldg.hh:323
sparsity pattern generator
Definition: pattern.hh:30
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
Definition: maxwelldg.hh:725
void preStage(typename T::Traits::RangeFieldType time, int r)
to be called once before each stage
Definition: maxwelldg.hh:693
Definition: maxwelldg.hh:733
Definition: maxwelldg.hh:317
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538
T::Traits::RangeFieldType suggestTimestep(typename T::Traits::RangeFieldType dt) const
to be called once before each stage
Definition: maxwelldg.hh:703
static const int dim
Definition: adaptivity.hh:82
Definition: maxwelldg.hh:318
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: maxwelldg.hh:761
Definition: maxwelldg.hh:321
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: maxwelldg.hh:334
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: maxwelldg.hh:545
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
DGMaxwellSpatialOperator(T ¶m_, int overintegration_=0)
Definition: maxwelldg.hh:327
const IG & ig
Definition: common/constraints.hh:146
Definition: maxwelldg.hh:322
static void eigenvectors(T1 eps, T1 mu, const Dune::FieldVector< T2, dim > &n, Dune::FieldMatrix< T3, 2 *dim, 2 *dim > &R)
Definition: maxwelldg.hh:73
store values of basis functions and gradients in a cache
Definition: localbasiscache.hh:16
Definition: maxwelldg.hh:301
static void eigenvalues(T1 eps, T1 mu, const Dune::FieldVector< T2, 2 *dim > &e)
Definition: maxwelldg.hh:50
Definition: maxwelldg.hh:736
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
const F & f
Definition: common/constraints.hh:145
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: maxwelldg.hh:430
Definition: maxwelldg.hh:28
sparsity pattern generator
Definition: pattern.hh:14
DGMaxwellTemporalOperator(T ¶m_, int overintegration_=0)
Definition: maxwelldg.hh:738
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 preStep(typename T::Traits::RangeFieldType time, typename T::Traits::RangeFieldType dt, int stages)
to be called once before each time step
Definition: maxwelldg.hh:687
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: maxwelldg.hh:804
void pattern_volume(const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
Definition: maxwelldg.hh:744
void postStage()
to be called once at the end of each stage
Definition: maxwelldg.hh:698
const E & e
Definition: interpolate.hh:172
void setTime(typename T::Traits::RangeFieldType t)
set time in parameter class
Definition: maxwelldg.hh:682
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: maxwelldg.hh:643
const std::string s
Definition: function.hh:1103