2 #ifndef DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONFEM_HH 3 #define DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONFEM_HH 7 #include<dune/geometry/referenceelements.hh> 35 template<
typename T,
typename FiniteElementMap>
54 : param(param_), intorderadd(intorderadd_)
59 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
60 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 63 using RF =
typename LFSU::Traits::FiniteElementType::
64 Traits::LocalBasisType::Traits::RangeFieldType;
65 using size_type =
typename LFSU::Traits::SizeType;
68 const int dim = EG::Entity::dimension;
71 const auto& cell = eg.entity();
74 auto geo = eg.geometry();
77 auto ref_el = referenceElement(geo);
78 auto localcenter = ref_el.position(0,0);
79 auto tensor = param.A(cell,localcenter);
82 std::vector<Dune::FieldVector<RF,dim> > gradphi(lfsu.size());
83 Dune::FieldVector<RF,dim> gradu(0.0);
84 Dune::FieldVector<RF,dim> Agradu(0.0);
87 typename EG::Geometry::JacobianInverseTransposed jac;
90 auto intorder = intorderadd+2*lfsu.finiteElement().localBasis().order();
94 auto& phi = cache.evaluateFunction(ip.position(),lfsu.finiteElement().localBasis());
98 for (size_type i=0; i<lfsu.size(); i++)
99 u += x(lfsu,i)*phi[i];
102 auto& js = cache.evaluateJacobian(ip.position(),lfsu.finiteElement().localBasis());
105 jac = geo.jacobianInverseTransposed(ip.position());
106 for (size_type i=0; i<lfsu.size(); i++)
107 jac.mv(js[i][0],gradphi[i]);
111 for (size_type i=0; i<lfsu.size(); i++)
112 gradu.axpy(x(lfsu,i),gradphi[i]);
115 tensor.mv(gradu,Agradu);
118 auto b = param.b(cell,ip.position());
119 auto c = param.c(cell,ip.position());
120 auto f = param.f(cell,ip.position());
123 RF factor = ip.weight() * geo.integrationElement(ip.position());
124 for (size_type i=0; i<lfsu.size(); i++)
125 r.accumulate(lfsu,i,( Agradu*gradphi[i] - u*(b*gradphi[i]) + (c*u-f)*phi[i] )*factor);
130 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
135 using RF =
typename LFSU::Traits::FiniteElementType::
136 Traits::LocalBasisType::Traits::RangeFieldType;
137 using size_type =
typename LFSU::Traits::SizeType;
140 const int dim = EG::Entity::dimension;
143 const auto& cell = eg.entity();
146 auto geo = eg.geometry();
149 auto ref_el = referenceElement(geo);
150 auto localcenter = ref_el.position(0,0);
151 auto tensor = param.A(cell,localcenter);
154 std::vector<Dune::FieldVector<RF,dim> > gradphi(lfsu.size());
155 std::vector<Dune::FieldVector<RF,dim> > Agradphi(lfsu.size());
158 typename EG::Geometry::JacobianInverseTransposed jac;
161 auto intorder = intorderadd+2*lfsu.finiteElement().localBasis().order();
165 auto& js = cache.evaluateJacobian(ip.position(),lfsu.finiteElement().localBasis());
168 jac = geo.jacobianInverseTransposed(ip.position());
169 for (size_type i=0; i<lfsu.size(); i++)
171 jac.mv(js[i][0],gradphi[i]);
172 tensor.mv(gradphi[i],Agradphi[i]);
176 auto& phi = cache.evaluateFunction(ip.position(),lfsu.finiteElement().localBasis());
179 auto b = param.b(cell,ip.position());
180 auto c = param.c(cell,ip.position());
183 RF factor = ip.weight() * geo.integrationElement(ip.position());
184 for (size_type j=0; j<lfsu.size(); j++)
185 for (size_type i=0; i<lfsu.size(); i++)
186 mat.accumulate(lfsu,i,lfsu,j,( Agradphi[j]*gradphi[i]-phi[j]*(b*gradphi[i])+c*phi[j]*phi[i] )*factor);
191 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
193 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
197 using RF =
typename LFSV::Traits::FiniteElementType::
198 Traits::LocalBasisType::Traits::RangeFieldType;
199 using size_type =
typename LFSV::Traits::SizeType;
202 const auto& cell_inside = ig.inside();
205 auto geo = ig.geometry();
208 auto geo_in_inside = ig.geometryInInside();
211 auto ref_el = referenceElement(geo_in_inside);
212 auto local_face_center = ref_el.position(0,0);
213 auto intersection = ig.intersection();
214 auto bctype = param.bctype(intersection,local_face_center);
220 auto intorder = intorderadd+2*lfsu_s.finiteElement().localBasis().order();
224 auto local = geo_in_inside.global(ip.position());
227 auto& phi = cache.evaluateFunction(local,lfsu_s.finiteElement().localBasis());
232 auto j = param.j(intersection,ip.position());
235 auto factor = ip.weight()*geo.integrationElement(ip.position());
236 for (size_type i=0; i<lfsu_s.size(); i++)
237 r_s.accumulate(lfsu_s,i,j*phi[i]*factor);
244 for (size_type i=0; i<lfsu_s.size(); i++)
245 u += x_s(lfsu_s,i)*phi[i];
248 auto b = param.b(cell_inside,local);
249 auto n = ig.unitOuterNormal(ip.position());
252 auto o = param.o(intersection,ip.position());
255 auto factor = ip.weight()*geo.integrationElement(ip.position());
256 for (size_type i=0; i<lfsu_s.size(); i++)
257 r_s.accumulate(lfsu_s,i,( (b*n)*u + o)*phi[i]*factor);
263 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename M>
265 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
269 using size_type =
typename LFSV::Traits::SizeType;
272 const auto& cell_inside = ig.inside();
275 auto geo = ig.geometry();
278 auto geo_in_inside = ig.geometryInInside();
281 auto ref_el = referenceElement(geo_in_inside);
282 auto local_face_center = ref_el.position(0,0);
283 auto intersection = ig.intersection();
284 auto bctype = param.bctype(intersection,local_face_center);
291 auto intorder = intorderadd+2*lfsu_s.finiteElement().localBasis().order();
295 auto local = geo_in_inside.global(ip.position());
298 auto& phi = cache.evaluateFunction(local,lfsu_s.finiteElement().localBasis());
301 auto b = param.b(cell_inside,local);
302 auto n = ig.unitOuterNormal(ip.position());
305 auto factor = ip.weight()*geo.integrationElement(ip.position());
306 for (size_type j=0; j<lfsu_s.size(); j++)
307 for (size_type i=0; i<lfsu_s.size(); i++)
308 mat_s.accumulate(lfsu_s,i,lfsu_s,j,(b*n)*phi[j]*phi[i]*factor);
322 using LocalBasisType =
typename FiniteElementMap::Traits::FiniteElementType::Traits::LocalBasisType;
347 enum {
dim = T::Traits::GridViewType::dimension };
349 using Real =
typename T::Traits::RangeFieldType;
368 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
369 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 372 using RF =
typename LFSU::Traits::FiniteElementType::
373 Traits::LocalBasisType::Traits::RangeFieldType;
374 using RangeType =
typename LFSU::Traits::FiniteElementType::
375 Traits::LocalBasisType::Traits::RangeType;
376 using size_type =
typename LFSU::Traits::SizeType;
379 const auto& cell = eg.entity();
382 auto geo = eg.geometry();
385 std::vector<RangeType> phi(lfsu.size());
389 auto intorder = 2*lfsu.finiteElement().localBasis().order();
393 lfsu.finiteElement().localBasis().evaluateFunction(ip.position(),phi);
397 for (size_type i=0; i<lfsu.size(); i++)
398 u += x(lfsu,i)*phi[i];
401 auto c = param.c(cell,ip.position());
404 auto f = param.f(cell,ip.position());
407 auto factor = ip.weight() * geo.integrationElement(ip.position());
408 sum += (f*f-c*c*u*u)*factor;
412 auto h_T = diameter(geo);
413 r.accumulate(lfsv,0,h_T*h_T*sum);
419 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
421 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
422 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
423 R& r_s, R& r_n)
const 426 using RF =
typename LFSU::Traits::FiniteElementType::
427 Traits::LocalBasisType::Traits::RangeFieldType;
428 using JacobianType =
typename LFSU::Traits::FiniteElementType::
429 Traits::LocalBasisType::Traits::JacobianType;
430 using size_type =
typename LFSU::Traits::SizeType;
433 const int dim = IG::Entity::dimension;
436 const auto& cell_inside = ig.inside();
437 const auto& cell_outside = ig.outside();
440 auto geo = ig.geometry();
441 auto geo_inside = cell_inside.geometry();
442 auto geo_outside = cell_outside.geometry();
445 auto geo_in_inside = ig.geometryInInside();
446 auto geo_in_outside = ig.geometryInOutside();
449 auto ref_el_inside = referenceElement(geo_inside);
450 auto ref_el_outside = referenceElement(geo_outside);
451 auto local_inside = ref_el_inside.position(0,0);
452 auto local_outside = ref_el_outside.position(0,0);
453 auto A_s = param.A(cell_inside,local_inside);
454 auto A_n = param.A(cell_outside,local_outside);
457 auto n_F = ig.centerUnitOuterNormal();
458 Dune::FieldVector<RF,dim> An_F_s;
460 Dune::FieldVector<RF,dim> An_F_n;
464 std::vector<JacobianType> gradphi_s(lfsu_s.size());
465 std::vector<JacobianType> gradphi_n(lfsu_n.size());
466 std::vector<Dune::FieldVector<RF,dim> > tgradphi_s(lfsu_s.size());
467 std::vector<Dune::FieldVector<RF,dim> > tgradphi_n(lfsu_n.size());
468 Dune::FieldVector<RF,dim> gradu_s(0.0);
469 Dune::FieldVector<RF,dim> gradu_n(0.0);
472 typename IG::Entity::Geometry::JacobianInverseTransposed jac;
476 auto intorder = 2*lfsu_s.finiteElement().localBasis().order();
480 auto iplocal_s = geo_in_inside.global(ip.position());
481 auto iplocal_n = geo_in_outside.global(ip.position());
484 lfsu_s.finiteElement().localBasis().evaluateJacobian(iplocal_s,gradphi_s);
485 lfsu_n.finiteElement().localBasis().evaluateJacobian(iplocal_n,gradphi_n);
488 jac = geo_inside.jacobianInverseTransposed(iplocal_s);
489 for (size_type i=0; i<lfsu_s.size(); i++) jac.mv(gradphi_s[i][0],tgradphi_s[i]);
490 jac = geo_outside.jacobianInverseTransposed(iplocal_n);
491 for (size_type i=0; i<lfsu_n.size(); i++) jac.mv(gradphi_n[i][0],tgradphi_n[i]);
495 for (size_type i=0; i<lfsu_s.size(); i++)
496 gradu_s.axpy(x_s(lfsu_s,i),tgradphi_s[i]);
498 for (size_type i=0; i<lfsu_n.size(); i++)
499 gradu_n.axpy(x_n(lfsu_n,i),tgradphi_n[i]);
502 auto factor = ip.weight() * geo.integrationElement(ip.position());
503 auto jump = (An_F_s*gradu_s)-(An_F_n*gradu_n);
504 sum += 0.25*jump*jump*factor;
509 auto h_T = std::max(diameter(geo_inside),diameter(geo_outside));
510 r_s.accumulate(lfsv_s,0,h_T*sum);
511 r_n.accumulate(lfsv_n,0,h_T*sum);
517 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
519 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
523 using RF =
typename LFSU::Traits::FiniteElementType::
524 Traits::LocalBasisType::Traits::RangeFieldType;
525 using JacobianType =
typename LFSU::Traits::FiniteElementType::
526 Traits::LocalBasisType::Traits::JacobianType;
527 using size_type =
typename LFSU::Traits::SizeType;
530 const int dim = IG::Entity::dimension;
533 const auto& cell_inside = ig.inside();
536 auto geo = ig.geometry();
537 auto geo_inside = cell_inside.geometry();
540 auto ref_el_inside = referenceElement(geo_inside);
541 auto local_inside = ref_el_inside.position(0,0);
542 auto A_s = param.A(cell_inside,local_inside);
545 auto n_F = ig.centerUnitOuterNormal();
546 Dune::FieldVector<RF,dim> An_F_s;
550 auto geo_in_inside = ig.geometryInInside();
553 auto ref_el_in_inside = referenceElement(geo_in_inside);
554 auto face_local = ref_el_in_inside.position(0,0);
555 auto bctype = param.bctype(ig.intersection(),face_local);
560 std::vector<JacobianType> gradphi_s(lfsu_s.size());
561 std::vector<Dune::FieldVector<RF,dim> > tgradphi_s(lfsu_s.size());
562 Dune::FieldVector<RF,dim> gradu_s(0.0);
565 typename IG::Entity::Geometry::JacobianInverseTransposed jac;
569 auto intorder = 2*lfsu_s.finiteElement().localBasis().order();
573 auto iplocal_s = geo_in_inside.global(ip.position());
576 lfsu_s.finiteElement().localBasis().evaluateJacobian(iplocal_s,gradphi_s);
579 jac = geo_inside.jacobianInverseTransposed(iplocal_s);
580 for (size_type i=0; i<lfsu_s.size(); i++) jac.mv(gradphi_s[i][0],tgradphi_s[i]);
584 for (size_type i=0; i<lfsu_s.size(); i++)
585 gradu_s.axpy(x_s(lfsu_s,i),tgradphi_s[i]);
588 auto j = param.j(ig.intersection(),ip.position());
591 auto factor = ip.weight() * geo.integrationElement(ip.position());
592 auto jump = j+(An_F_s*gradu_s);
593 sum += jump*jump*factor;
598 auto h_T = diameter(geo_inside);
599 r_s.accumulate(lfsv_s,0,h_T*sum);
606 typename GEO::ctype diameter (
const GEO& geo)
const 608 using DF =
typename GEO::ctype;
610 const int dim = GEO::coorddimension;
611 for (
int i=0; i<geo.corners(); i++)
613 Dune::FieldVector<DF,dim> xi = geo.corner(i);
614 for (
int j=i+1; j<geo.corners(); j++)
616 Dune::FieldVector<DF,dim> xj = geo.corner(j);
618 hmax = std::max(hmax,xj.two_norm());
647 enum {
dim = T::Traits::GridViewType::dimension };
649 using Real =
typename T::Traits::RangeFieldType;
663 : param(param_), time(time_), dt(dt_), cmax(0)
667 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
668 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 671 using RF =
typename LFSU::Traits::FiniteElementType::
672 Traits::LocalBasisType::Traits::RangeFieldType;
673 using RangeType =
typename LFSU::Traits::FiniteElementType::
674 Traits::LocalBasisType::Traits::RangeType;
675 using size_type =
typename LFSU::Traits::SizeType;
678 const auto& cell = eg.entity();
681 auto geo = eg.geometry();
684 std::vector<RangeType> phi(lfsu.size());
691 auto intorder = 2*lfsu.finiteElement().localBasis().order();
695 lfsu.finiteElement().localBasis().evaluateFunction(ip.position(),phi);
699 for (size_type i=0; i<lfsu.size(); i++)
700 u += x(lfsu,i)*phi[i];
703 auto factor = ip.weight() * geo.integrationElement(ip.position());
708 auto f_down = param.f(cell,ip.position());
709 param.setTime(time+0.5*dt);
710 auto f_mid = param.f(cell,ip.position());
711 param.setTime(time+dt);
712 auto f_up = param.f(cell,ip.position());
713 auto f_average = (1.0/6.0)*f_down + (2.0/3.0)*f_mid + (1.0/6.0)*f_up;
716 fsum_down += (f_down-f_average)*(f_down-f_average)*factor;
717 fsum_mid += (f_mid-f_average)*(f_mid-f_average)*factor;
718 fsum_up += (f_up-f_average)*(f_up-f_average)*factor;
722 auto h_T = diameter(geo);
723 r.accumulate(lfsv,0,(h_T*h_T/dt)*sum);
724 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) );
744 typename GEO::ctype diameter (
const GEO& geo)
const 746 using DF =
typename GEO::ctype;
748 const int dim = GEO::coorddimension;
749 for (
int i=0; i<geo.corners(); i++)
751 Dune::FieldVector<DF,dim> xi = geo.corner(i);
752 for (
int j=i+1; j<geo.corners(); j++)
754 Dune::FieldVector<DF,dim> xj = geo.corner(j);
756 hmax = std::max(hmax,xj.two_norm());
765 template<
typename T,
typename EG>
772 template<
typename X,
typename Y>
775 y[0] = t.f(eg.entity(),x);
803 enum {
dim = T::Traits::GridViewType::dimension };
805 using Real =
typename T::Traits::RangeFieldType;
820 : param(param_), time(time_), dt(dt_)
824 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
825 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 828 using RF =
typename LFSU::Traits::FiniteElementType::
829 Traits::LocalBasisType::Traits::RangeFieldType;
830 using RangeType =
typename LFSU::Traits::FiniteElementType::
831 Traits::LocalBasisType::Traits::RangeType;
832 using size_type =
typename LFSU::Traits::SizeType;
833 using JacobianType =
typename LFSU::Traits::FiniteElementType::
834 Traits::LocalBasisType::Traits::JacobianType;
837 const int dim = EG::Geometry::mydimension;
840 auto geo = eg.geometry();
844 std::vector<RF> f_up, f_down, f_mid;
846 lfsu.finiteElement().localInterpolation().interpolate(f_adapter,f_down);
847 param.setTime(time+0.5*dt);
848 lfsu.finiteElement().localInterpolation().interpolate(f_adapter,f_mid);
849 param.setTime(time+dt);
850 lfsu.finiteElement().localInterpolation().interpolate(f_adapter,f_up);
853 std::vector<JacobianType> js(lfsu.size());
854 std::vector<Dune::FieldVector<RF,dim> > gradphi(lfsu.size());
855 Dune::FieldVector<RF,dim> gradu(0.0);
856 Dune::FieldVector<RF,dim> gradf_down(0.0);
857 Dune::FieldVector<RF,dim> gradf_mid(0.0);
858 Dune::FieldVector<RF,dim> gradf_up(0.0);
859 Dune::FieldVector<RF,dim> gradf_average(0.0);
862 typename EG::Geometry::JacobianInverseTransposed jac;
867 RF fsum_grad_up(0.0);
868 RF fsum_grad_mid(0.0);
869 RF fsum_grad_down(0.0);
870 auto intorder = 2*lfsu.finiteElement().localBasis().order();
874 std::vector<RangeType> phi(lfsu.size());
875 lfsu.finiteElement().localBasis().evaluateFunction(ip.position(),phi);
879 for (size_type i=0; i<lfsu.size(); i++)
880 u += x(lfsu,i)*phi[i];
883 auto factor = ip.weight() * geo.integrationElement(ip.position());
887 lfsu.finiteElement().localBasis().evaluateJacobian(ip.position(),js);
890 jac = geo.jacobianInverseTransposed(ip.position());
891 for (size_type i=0; i<lfsu.size(); i++)
892 jac.mv(js[i][0],gradphi[i]);
896 for (size_type i=0; i<lfsu.size(); i++)
897 gradu.axpy(x(lfsu,i),gradphi[i]);
900 sum_grad += (gradu*gradu)*factor;
904 for (size_type i=0; i<lfsu.size(); i++) gradf_down.axpy(f_down[i],gradphi[i]);
906 for (size_type i=0; i<lfsu.size(); i++) gradf_mid.axpy(f_mid[i],gradphi[i]);
908 for (size_type i=0; i<lfsu.size(); i++) gradf_up.axpy(f_up[i],gradphi[i]);
910 for (size_type i=0; i<lfsu.size(); i++)
911 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]);
914 gradf_down -= gradf_average;
915 fsum_grad_down += (gradf_down*gradf_down)*factor;
916 gradf_mid -= gradf_average;
917 fsum_grad_mid += (gradf_mid*gradf_mid)*factor;
918 gradf_up -= gradf_average;
919 fsum_grad_up += (gradf_up*gradf_up)*factor;
923 auto h_T = diameter(geo);
924 r.accumulate(lfsv,0,dt * sum_grad);
925 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));
930 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
932 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
936 using RF =
typename LFSU::Traits::FiniteElementType::
937 Traits::LocalBasisType::Traits::RangeFieldType;
940 const auto& cell_inside = ig.inside();
943 auto geo = ig.geometry();
944 auto geo_inside = cell_inside.geometry();
947 auto geo_in_inside = ig.geometryInInside();
950 auto ref_el_in_inside = referenceElement(geo_in_inside);
951 auto face_local = ref_el_in_inside.position(0,0);
952 auto bctype = param.bctype(ig.intersection(),face_local);
959 const int intorder = 2*lfsu_s.finiteElement().localBasis().order();
964 auto j_down = param.j(ig.intersection(),ip.position());
965 param.setTime(time+0.5*dt);
966 auto j_mid = param.j(ig.intersection(),ip.position());
967 param.setTime(time+dt);
968 auto j_up = param.j(ig.intersection(),ip.position());
971 auto factor = ip.weight() * geo.integrationElement(ip.position());
972 sum_down += (j_down-j_mid)*(j_down-j_mid)*factor;
973 sum_up += (j_up-j_mid)*(j_up-j_mid)*factor;
978 auto h_T = diameter(geo_inside);
979 r_s.accumulate(lfsv_s,0,(h_T+dt*dt/h_T)*dt*0.5*(sum_down+sum_up));
988 typename GEO::ctype diameter (
const GEO& geo)
const 990 using DF =
typename GEO::ctype;
992 const int dim = GEO::coorddimension;
993 for (
int i=0; i<geo.corners(); i++)
995 Dune::FieldVector<DF,dim> xi = geo.corner(i);
996 for (
int j=i+1; j<geo.corners(); j++)
998 Dune::FieldVector<DF,dim> xj = geo.corner(j);
1000 hmax = std::max(hmax,xj.two_norm());
1010 #endif // DUNE_PDELAB_LOCALOPERATOR_CONVECTIONDIFFUSIONFEM_HH
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: convectiondiffusionfem.hh:131
static const int dim
Definition: adaptivity.hh:84
Default flags for all local operators.
Definition: flags.hh:18
Definition: convectiondiffusionfem.hh:766
Definition: convectiondiffusionfem.hh:36
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:931
ConvectionDiffusionTemporalResidualEstimator1(T ¶m_, double time_, double dt_)
constructor: pass parameter object
Definition: convectiondiffusionfem.hh:662
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:825
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
Definition: convectiondiffusionfem.hh:643
ConvectionDiffusionTemporalResidualEstimator(T ¶m_, double time_, double dt_)
constructor: pass parameter object
Definition: convectiondiffusionfem.hh:819
sparsity pattern generator
Definition: pattern.hh:13
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:60
CD_RHS_LocalAdapter(const T &t_, const EG &eg_)
Definition: convectiondiffusionfem.hh:769
void setTime(double t)
set time in parameter class
Definition: convectiondiffusionfem.hh:314
void clearCmax()
Definition: convectiondiffusionfem.hh:727
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionparameter.hh:65
Implement jacobian_boundary() based on alpha_boundary()
Definition: numericaljacobian.hh:250
Implement jacobian_volume() based on alpha_volume()
Definition: numericaljacobian.hh:31
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:264
store values of basis functions and gradients in a cache
Definition: localbasiscache.hh:17
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
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
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:420
Implements linear and nonlinear versions of jacobian_apply_boundary() based on alpha_boundary() ...
Definition: numericaljacobianapply.hh:285
double getCmax() const
Definition: convectiondiffusionfem.hh:732
const IG & ig
Definition: constraints.hh:148
Implements linear and nonlinear versions of jacobian_apply_volume() based on alpha_volume() ...
Definition: numericaljacobianapply.hh:32
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:518
Definition: convectiondiffusionfem.hh:47
Definition: convectiondiffusionfem.hh:51
Definition: convectiondiffusionfem.hh:799
Type
Definition: convectiondiffusionparameter.hh:65
Definition: convectiondiffusionfem.hh:344
void evaluate(const X &x, Y &y) const
Definition: convectiondiffusionfem.hh:773
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:668
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: convectiondiffusionfem.hh:369
Definition: convectiondiffusionparameter.hh:65
ConvectionDiffusionFEMResidualEstimator(T ¶m_)
constructor: pass parameter object
Definition: convectiondiffusionfem.hh:363
Definition: convectiondiffusionfem.hh:50
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:192
ConvectionDiffusionFEM(T ¶m_, int intorderadd_=0)
Definition: convectiondiffusionfem.hh:53