2 #ifndef DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONFEM_HH 3 #define DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONFEM_HH 34 template<
typename T,
typename FiniteElementMap>
53 : param(param_), intorderadd(intorderadd_)
58 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
59 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 62 using RF =
typename LFSU::Traits::FiniteElementType::
63 Traits::LocalBasisType::Traits::RangeFieldType;
64 using size_type =
typename LFSU::Traits::SizeType;
67 const int dim = EG::Entity::dimension;
70 const auto& cell = eg.entity();
73 auto geo = eg.geometry();
77 auto localcenter = ref_el.position(0,0);
78 auto tensor = param.A(cell,localcenter);
81 std::vector<Dune::FieldVector<RF,dim> > gradphi(lfsu.size());
82 Dune::FieldVector<RF,dim> gradu(0.0);
83 Dune::FieldVector<RF,dim> Agradu(0.0);
86 typename EG::Geometry::JacobianInverseTransposed jac;
89 auto intorder = intorderadd+2*lfsu.finiteElement().localBasis().order();
93 auto& phi = cache.evaluateFunction(ip.position(),lfsu.finiteElement().localBasis());
97 for (size_type i=0; i<lfsu.size(); i++)
98 u += x(lfsu,i)*phi[i];
101 auto& js = cache.evaluateJacobian(ip.position(),lfsu.finiteElement().localBasis());
104 jac = geo.jacobianInverseTransposed(ip.position());
105 for (size_type i=0; i<lfsu.size(); i++)
106 jac.mv(js[i][0],gradphi[i]);
110 for (size_type i=0; i<lfsu.size(); i++)
111 gradu.axpy(x(lfsu,i),gradphi[i]);
114 tensor.mv(gradu,Agradu);
117 auto b = param.b(cell,ip.position());
118 auto c = param.c(cell,ip.position());
119 auto f = param.f(cell,ip.position());
122 RF factor = ip.weight() * geo.integrationElement(ip.position());
123 for (size_type i=0; i<lfsu.size(); i++)
124 r.accumulate(lfsu,i,( Agradu*gradphi[i] - u*(b*gradphi[i]) + (c*u-f)*phi[i] )*factor);
129 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
134 using RF =
typename LFSU::Traits::FiniteElementType::
135 Traits::LocalBasisType::Traits::RangeFieldType;
136 using size_type =
typename LFSU::Traits::SizeType;
139 const int dim = EG::Entity::dimension;
142 const auto& cell = eg.entity();
145 auto geo = eg.geometry();
149 auto localcenter = ref_el.position(0,0);
150 auto tensor = param.A(cell,localcenter);
153 std::vector<Dune::FieldVector<RF,dim> > gradphi(lfsu.size());
154 std::vector<Dune::FieldVector<RF,dim> > Agradphi(lfsu.size());
157 typename EG::Geometry::JacobianInverseTransposed jac;
160 auto intorder = intorderadd+2*lfsu.finiteElement().localBasis().order();
164 auto& js = cache.evaluateJacobian(ip.position(),lfsu.finiteElement().localBasis());
167 jac = geo.jacobianInverseTransposed(ip.position());
168 for (size_type i=0; i<lfsu.size(); i++)
170 jac.mv(js[i][0],gradphi[i]);
171 tensor.mv(gradphi[i],Agradphi[i]);
175 auto& phi = cache.evaluateFunction(ip.position(),lfsu.finiteElement().localBasis());
178 auto b = param.b(cell,ip.position());
179 auto c = param.c(cell,ip.position());
182 RF factor = ip.weight() * geo.integrationElement(ip.position());
183 for (size_type j=0; j<lfsu.size(); j++)
184 for (size_type i=0; i<lfsu.size(); i++)
185 mat.accumulate(lfsu,i,lfsu,j,( Agradphi[j]*gradphi[i]-phi[j]*(b*gradphi[i])+c*phi[j]*phi[i] )*factor);
190 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
192 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
196 using RF =
typename LFSV::Traits::FiniteElementType::
197 Traits::LocalBasisType::Traits::RangeFieldType;
198 using size_type =
typename LFSV::Traits::SizeType;
201 const auto& cell_inside = ig.inside();
204 auto geo = ig.geometry();
207 auto geo_in_inside = ig.geometryInInside();
211 auto local_face_center = ref_el.position(0,0);
212 auto intersection = ig.intersection();
213 auto bctype = param.bctype(intersection,local_face_center);
219 auto intorder = intorderadd+2*lfsu_s.finiteElement().localBasis().order();
223 auto local = geo_in_inside.global(ip.position());
226 auto& phi = cache.evaluateFunction(local,lfsu_s.finiteElement().localBasis());
231 auto j = param.j(intersection,ip.position());
234 auto factor = ip.weight()*geo.integrationElement(ip.position());
235 for (size_type i=0; i<lfsu_s.size(); i++)
236 r_s.accumulate(lfsu_s,i,j*phi[i]*factor);
243 for (size_type i=0; i<lfsu_s.size(); i++)
244 u += x_s(lfsu_s,i)*phi[i];
247 auto b = param.b(cell_inside,local);
248 auto n = ig.unitOuterNormal(ip.position());
251 auto o = param.o(intersection,ip.position());
254 auto factor = ip.weight()*geo.integrationElement(ip.position());
255 for (size_type i=0; i<lfsu_s.size(); i++)
256 r_s.accumulate(lfsu_s,i,( (b*n)*u + o)*phi[i]*factor);
262 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename M>
264 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
268 using size_type =
typename LFSV::Traits::SizeType;
271 const auto& cell_inside = ig.inside();
274 auto geo = ig.geometry();
277 auto geo_in_inside = ig.geometryInInside();
281 auto local_face_center = ref_el.position(0,0);
282 auto intersection = ig.intersection();
283 auto bctype = param.bctype(intersection,local_face_center);
290 auto intorder = intorderadd+2*lfsu_s.finiteElement().localBasis().order();
294 auto local = geo_in_inside.global(ip.position());
297 auto& phi = cache.evaluateFunction(local,lfsu_s.finiteElement().localBasis());
300 auto b = param.b(cell_inside,local);
301 auto n = ig.unitOuterNormal(ip.position());
304 auto factor = ip.weight()*geo.integrationElement(ip.position());
305 for (size_type j=0; j<lfsu_s.size(); j++)
306 for (size_type i=0; i<lfsu_s.size(); i++)
307 mat_s.accumulate(lfsu_s,i,lfsu_s,j,(b*n)*phi[j]*phi[i]*factor);
321 using LocalBasisType =
typename FiniteElementMap::Traits::FiniteElementType::Traits::LocalBasisType;
346 enum {
dim = T::Traits::GridViewType::dimension };
348 using Real =
typename T::Traits::RangeFieldType;
367 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
368 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 371 using RF =
typename LFSU::Traits::FiniteElementType::
372 Traits::LocalBasisType::Traits::RangeFieldType;
373 using RangeType =
typename LFSU::Traits::FiniteElementType::
374 Traits::LocalBasisType::Traits::RangeType;
375 using size_type =
typename LFSU::Traits::SizeType;
378 const auto& cell = eg.entity();
381 auto geo = eg.geometry();
384 std::vector<RangeType> phi(lfsu.size());
388 auto intorder = 2*lfsu.finiteElement().localBasis().order();
392 lfsu.finiteElement().localBasis().evaluateFunction(ip.position(),phi);
396 for (size_type i=0; i<lfsu.size(); i++)
397 u += x(lfsu,i)*phi[i];
400 auto c = param.c(cell,ip.position());
403 auto f = param.f(cell,ip.position());
406 auto factor = ip.weight() * geo.integrationElement(ip.position());
407 sum += (f*f-c*c*u*u)*factor;
411 auto h_T = diameter(geo);
412 r.accumulate(lfsv,0,h_T*h_T*sum);
418 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
420 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
421 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
422 R& r_s, R& r_n)
const 425 using RF =
typename LFSU::Traits::FiniteElementType::
426 Traits::LocalBasisType::Traits::RangeFieldType;
427 using JacobianType =
typename LFSU::Traits::FiniteElementType::
428 Traits::LocalBasisType::Traits::JacobianType;
429 using size_type =
typename LFSU::Traits::SizeType;
432 const int dim = IG::dimension;
435 const auto& cell_inside = ig.inside();
436 const auto& cell_outside = ig.outside();
439 auto geo = ig.geometry();
440 auto geo_inside = cell_inside.geometry();
441 auto geo_outside = cell_outside.geometry();
444 auto geo_in_inside = ig.geometryInInside();
445 auto geo_in_outside = ig.geometryInOutside();
450 auto local_inside = ref_el_inside.position(0,0);
451 auto local_outside = ref_el_outside.position(0,0);
452 auto A_s = param.A(cell_inside,local_inside);
453 auto A_n = param.A(cell_outside,local_outside);
456 auto n_F = ig.centerUnitOuterNormal();
457 Dune::FieldVector<RF,dim> An_F_s;
459 Dune::FieldVector<RF,dim> An_F_n;
463 std::vector<JacobianType> gradphi_s(lfsu_s.size());
464 std::vector<JacobianType> gradphi_n(lfsu_n.size());
465 std::vector<Dune::FieldVector<RF,dim> > tgradphi_s(lfsu_s.size());
466 std::vector<Dune::FieldVector<RF,dim> > tgradphi_n(lfsu_n.size());
467 Dune::FieldVector<RF,dim> gradu_s(0.0);
468 Dune::FieldVector<RF,dim> gradu_n(0.0);
471 typename IG::Entity::Geometry::JacobianInverseTransposed jac;
475 auto intorder = 2*lfsu_s.finiteElement().localBasis().order();
479 auto iplocal_s = geo_in_inside.global(ip.position());
480 auto iplocal_n = geo_in_outside.global(ip.position());
483 lfsu_s.finiteElement().localBasis().evaluateJacobian(iplocal_s,gradphi_s);
484 lfsu_n.finiteElement().localBasis().evaluateJacobian(iplocal_n,gradphi_n);
487 jac = geo_inside.jacobianInverseTransposed(iplocal_s);
488 for (size_type i=0; i<lfsu_s.size(); i++) jac.mv(gradphi_s[i][0],tgradphi_s[i]);
489 jac = geo_outside.jacobianInverseTransposed(iplocal_n);
490 for (size_type i=0; i<lfsu_n.size(); i++) jac.mv(gradphi_n[i][0],tgradphi_n[i]);
494 for (size_type i=0; i<lfsu_s.size(); i++)
495 gradu_s.axpy(x_s(lfsu_s,i),tgradphi_s[i]);
497 for (size_type i=0; i<lfsu_n.size(); i++)
498 gradu_n.axpy(x_n(lfsu_n,i),tgradphi_n[i]);
501 auto factor = ip.weight() * geo.integrationElement(ip.position());
502 auto jump = (An_F_s*gradu_s)-(An_F_n*gradu_n);
503 sum += 0.25*jump*jump*factor;
508 auto h_T = std::max(diameter(geo_inside),diameter(geo_outside));
509 r_s.accumulate(lfsv_s,0,h_T*sum);
510 r_n.accumulate(lfsv_n,0,h_T*sum);
516 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
518 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
522 using RF =
typename LFSU::Traits::FiniteElementType::
523 Traits::LocalBasisType::Traits::RangeFieldType;
524 using JacobianType =
typename LFSU::Traits::FiniteElementType::
525 Traits::LocalBasisType::Traits::JacobianType;
526 using size_type =
typename LFSU::Traits::SizeType;
529 const int dim = IG::dimension;
532 const auto& cell_inside = ig.inside();
535 auto geo = ig.geometry();
536 auto geo_inside = cell_inside.geometry();
540 auto local_inside = ref_el_inside.position(0,0);
541 auto A_s = param.A(cell_inside,local_inside);
544 auto n_F = ig.centerUnitOuterNormal();
545 Dune::FieldVector<RF,dim> An_F_s;
549 auto geo_in_inside = ig.geometryInInside();
553 auto face_local = ref_el_in_inside.position(0,0);
554 auto bctype = param.bctype(ig.intersection(),face_local);
559 std::vector<JacobianType> gradphi_s(lfsu_s.size());
560 std::vector<Dune::FieldVector<RF,dim> > tgradphi_s(lfsu_s.size());
561 Dune::FieldVector<RF,dim> gradu_s(0.0);
564 typename IG::Entity::Geometry::JacobianInverseTransposed jac;
568 auto intorder = 2*lfsu_s.finiteElement().localBasis().order();
572 auto iplocal_s = geo_in_inside.global(ip.position());
575 lfsu_s.finiteElement().localBasis().evaluateJacobian(iplocal_s,gradphi_s);
578 jac = geo_inside.jacobianInverseTransposed(iplocal_s);
579 for (size_type i=0; i<lfsu_s.size(); i++) jac.mv(gradphi_s[i][0],tgradphi_s[i]);
583 for (size_type i=0; i<lfsu_s.size(); i++)
584 gradu_s.axpy(x_s(lfsu_s,i),tgradphi_s[i]);
587 auto j = param.j(ig.intersection(),ip.position());
590 auto factor = ip.weight() * geo.integrationElement(ip.position());
591 auto jump = j+(An_F_s*gradu_s);
592 sum += jump*jump*factor;
597 auto h_T = diameter(geo_inside);
598 r_s.accumulate(lfsv_s,0,h_T*sum);
605 typename GEO::ctype diameter (
const GEO& geo)
const 607 using DF =
typename GEO::ctype;
609 const int dim = GEO::coorddimension;
610 for (
int i=0; i<geo.corners(); i++)
612 Dune::FieldVector<DF,dim> xi = geo.corner(i);
613 for (
int j=i+1; j<geo.corners(); j++)
615 Dune::FieldVector<DF,dim> xj = geo.corner(j);
617 hmax = std::max(hmax,xj.two_norm());
646 enum {
dim = T::Traits::GridViewType::dimension };
648 using Real =
typename T::Traits::RangeFieldType;
662 : param(param_), time(time_), dt(dt_), cmax(0)
666 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
667 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 670 using RF =
typename LFSU::Traits::FiniteElementType::
671 Traits::LocalBasisType::Traits::RangeFieldType;
672 using RangeType =
typename LFSU::Traits::FiniteElementType::
673 Traits::LocalBasisType::Traits::RangeType;
674 using size_type =
typename LFSU::Traits::SizeType;
677 const auto& cell = eg.entity();
680 auto geo = eg.geometry();
683 std::vector<RangeType> phi(lfsu.size());
690 auto intorder = 2*lfsu.finiteElement().localBasis().order();
694 lfsu.finiteElement().localBasis().evaluateFunction(ip.position(),phi);
698 for (size_type i=0; i<lfsu.size(); i++)
699 u += x(lfsu,i)*phi[i];
702 auto factor = ip.weight() * geo.integrationElement(ip.position());
707 auto f_down = param.f(cell,ip.position());
708 param.setTime(time+0.5*dt);
709 auto f_mid = param.f(cell,ip.position());
710 param.setTime(time+dt);
711 auto f_up = param.f(cell,ip.position());
712 auto f_average = (1.0/6.0)*f_down + (2.0/3.0)*f_mid + (1.0/6.0)*f_up;
715 fsum_down += (f_down-f_average)*(f_down-f_average)*factor;
716 fsum_mid += (f_mid-f_average)*(f_mid-f_average)*factor;
717 fsum_up += (f_up-f_average)*(f_up-f_average)*factor;
721 auto h_T = diameter(geo);
722 r.accumulate(lfsv,0,(h_T*h_T/dt)*sum);
723 r.accumulate(lfsv,0,h_T*h_T * dt*((1.0/6.0)*fsum_down+(2.0/3.0)*fsum_mid+(1.0/6.0)*fsum_up) );
743 typename GEO::ctype diameter (
const GEO& geo)
const 745 using DF =
typename GEO::ctype;
747 const int dim = GEO::coorddimension;
748 for (
int i=0; i<geo.corners(); i++)
750 Dune::FieldVector<DF,dim> xi = geo.corner(i);
751 for (
int j=i+1; j<geo.corners(); j++)
753 Dune::FieldVector<DF,dim> xj = geo.corner(j);
755 hmax = std::max(hmax,xj.two_norm());
764 template<
typename T,
typename EG>
771 template<
typename X,
typename Y>
774 y[0] = t.f(eg.entity(),x);
802 enum {
dim = T::Traits::GridViewType::dimension };
804 using Real =
typename T::Traits::RangeFieldType;
819 : param(param_), time(time_), dt(dt_)
823 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
824 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 827 using RF =
typename LFSU::Traits::FiniteElementType::
828 Traits::LocalBasisType::Traits::RangeFieldType;
829 using RangeType =
typename LFSU::Traits::FiniteElementType::
830 Traits::LocalBasisType::Traits::RangeType;
831 using size_type =
typename LFSU::Traits::SizeType;
832 using JacobianType =
typename LFSU::Traits::FiniteElementType::
833 Traits::LocalBasisType::Traits::JacobianType;
836 const int dim = EG::Geometry::mydimension;
839 auto geo = eg.geometry();
843 std::vector<RF> f_up, f_down, f_mid;
845 lfsu.finiteElement().localInterpolation().interpolate(f_adapter,f_down);
846 param.setTime(time+0.5*dt);
847 lfsu.finiteElement().localInterpolation().interpolate(f_adapter,f_mid);
848 param.setTime(time+dt);
849 lfsu.finiteElement().localInterpolation().interpolate(f_adapter,f_up);
852 std::vector<JacobianType> js(lfsu.size());
853 std::vector<Dune::FieldVector<RF,dim> > gradphi(lfsu.size());
854 Dune::FieldVector<RF,dim> gradu(0.0);
855 Dune::FieldVector<RF,dim> gradf_down(0.0);
856 Dune::FieldVector<RF,dim> gradf_mid(0.0);
857 Dune::FieldVector<RF,dim> gradf_up(0.0);
858 Dune::FieldVector<RF,dim> gradf_average(0.0);
861 typename EG::Geometry::JacobianInverseTransposed jac;
866 RF fsum_grad_up(0.0);
867 RF fsum_grad_mid(0.0);
868 RF fsum_grad_down(0.0);
869 auto intorder = 2*lfsu.finiteElement().localBasis().order();
873 std::vector<RangeType> phi(lfsu.size());
874 lfsu.finiteElement().localBasis().evaluateFunction(ip.position(),phi);
878 for (size_type i=0; i<lfsu.size(); i++)
879 u += x(lfsu,i)*phi[i];
882 auto factor = ip.weight() * geo.integrationElement(ip.position());
886 lfsu.finiteElement().localBasis().evaluateJacobian(ip.position(),js);
889 jac = geo.jacobianInverseTransposed(ip.position());
890 for (size_type i=0; i<lfsu.size(); i++)
891 jac.mv(js[i][0],gradphi[i]);
895 for (size_type i=0; i<lfsu.size(); i++)
896 gradu.axpy(x(lfsu,i),gradphi[i]);
899 sum_grad += (gradu*gradu)*factor;
903 for (size_type i=0; i<lfsu.size(); i++) gradf_down.axpy(f_down[i],gradphi[i]);
905 for (size_type i=0; i<lfsu.size(); i++) gradf_mid.axpy(f_mid[i],gradphi[i]);
907 for (size_type i=0; i<lfsu.size(); i++) gradf_up.axpy(f_up[i],gradphi[i]);
909 for (size_type i=0; i<lfsu.size(); i++)
910 gradf_average.axpy((1.0/6.0)*f_down[i]+(2.0/3.0)*f_mid[i]+(1.0/6.0)*f_up[i],gradphi[i]);
913 gradf_down -= gradf_average;
914 fsum_grad_down += (gradf_down*gradf_down)*factor;
915 gradf_mid -= gradf_average;
916 fsum_grad_mid += (gradf_mid*gradf_mid)*factor;
917 gradf_up -= gradf_average;
918 fsum_grad_up += (gradf_up*gradf_up)*factor;
922 auto h_T = diameter(geo);
923 r.accumulate(lfsv,0,dt * sum_grad);
924 r.accumulate(lfsv,0,dt*dt * dt*((1.0/6.0)*fsum_grad_down+(2.0/3.0)*fsum_grad_mid+(1.0/6.0)*fsum_grad_up));
929 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
931 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
935 using RF =
typename LFSU::Traits::FiniteElementType::
936 Traits::LocalBasisType::Traits::RangeFieldType;
939 const auto& cell_inside = ig.inside();
942 auto geo = ig.geometry();
943 auto geo_inside = cell_inside.geometry();
946 auto geo_in_inside = ig.geometryInInside();
950 auto face_local = ref_el_in_inside.position(0,0);
951 auto bctype = param.bctype(ig.intersection(),face_local);
958 const int intorder = 2*lfsu_s.finiteElement().localBasis().order();
963 auto j_down = param.j(ig.intersection(),ip.position());
964 param.setTime(time+0.5*dt);
965 auto j_mid = param.j(ig.intersection(),ip.position());
966 param.setTime(time+dt);
967 auto j_up = param.j(ig.intersection(),ip.position());
970 auto factor = ip.weight() * geo.integrationElement(ip.position());
971 sum_down += (j_down-j_mid)*(j_down-j_mid)*factor;
972 sum_up += (j_up-j_mid)*(j_up-j_mid)*factor;
977 auto h_T = diameter(geo_inside);
978 r_s.accumulate(lfsv_s,0,(h_T+dt*dt/h_T)*dt*0.5*(sum_down+sum_up));
987 typename GEO::ctype diameter (
const GEO& geo)
const 989 using DF =
typename GEO::ctype;
991 const int dim = GEO::coorddimension;
992 for (
int i=0; i<geo.corners(); i++)
994 Dune::FieldVector<DF,dim> xi = geo.corner(i);
995 for (
int j=i+1; j<geo.corners(); j++)
997 Dune::FieldVector<DF,dim> xj = geo.corner(j);
999 hmax = std::max(hmax,xj.two_norm());
1009 #endif // DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONFEM_HH void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:824
const IG & ig
Definition: constraints.hh:148
static const int dim
Definition: adaptivity.hh:83
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:667
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:368
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: convectiondiffusionfem.hh:930
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: convectiondiffusionfem.hh:191
Implements linear and nonlinear versions of jacobian_apply_boundary() based on alpha_boundary() ...
Definition: numericaljacobianapply.hh:285
CD_RHS_LocalAdapter(const T &t_, const EG &eg_)
Definition: convectiondiffusionfem.hh:768
void evaluate(const X &x, Y &y) const
Definition: convectiondiffusionfem.hh:772
void clearCmax()
Definition: convectiondiffusionfem.hh:726
Definition: convectiondiffusionfem.hh:50
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: convectiondiffusionfem.hh:130
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:59
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionfem.hh:765
Definition: convectiondiffusionfem.hh:46
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: convectiondiffusionfem.hh:517
ConvectionDiffusionTemporalResidualEstimator1(T ¶m_, double time_, double dt_)
constructor: pass parameter object
Definition: convectiondiffusionfem.hh:661
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
Type
Definition: convectiondiffusionparameter.hh:65
void jacobian_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, M &mat_s) const
Definition: convectiondiffusionfem.hh:263
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionfem.hh:35
void setTime(double t)
set time in parameter class
Definition: convectiondiffusionfem.hh:313
ConvectionDiffusionFEMResidualEstimator(T ¶m_)
constructor: pass parameter object
Definition: convectiondiffusionfem.hh:362
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: convectiondiffusionfem.hh:419
Definition: convectiondiffusionfem.hh:642
Implements linear and nonlinear versions of jacobian_apply_volume() based on alpha_volume() ...
Definition: numericaljacobianapply.hh:32
double getCmax() const
Definition: convectiondiffusionfem.hh:731
Definition: convectiondiffusionfem.hh:49
ConvectionDiffusionTemporalResidualEstimator(T ¶m_, double time_, double dt_)
constructor: pass parameter object
Definition: convectiondiffusionfem.hh:818
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionfem.hh:798
sparsity pattern generator
Definition: pattern.hh:13
Definition: convectiondiffusionfem.hh:343
QuadratureRuleWrapper< QuadratureRule< typename Geometry::ctype, Geometry::mydimension > > quadratureRule(const Geometry &geo, std::size_t order, QuadratureType::Enum quadrature_type=QuadratureType::GaussLegendre)
Returns a quadrature rule for the given geometry.
Definition: quadraturerules.hh:117
ReferenceElementWrapper< ReferenceElement< typename Geometry::ctype, Geometry::mydimension > > referenceElement(const Geometry &geo)
Returns the reference element for the given geometry.
Definition: referenceelements.hh:144
store values of basis functions and gradients in a cache
Definition: localbasiscache.hh:17
Default flags for all local operators.
Definition: flags.hh:18
Implement jacobian_boundary() based on alpha_boundary()
Definition: numericaljacobian.hh:250
Implement jacobian_volume() based on alpha_volume()
Definition: numericaljacobian.hh:31
ConvectionDiffusionFEM(T ¶m_, int intorderadd_=0)
Definition: convectiondiffusionfem.hh:52