2 #ifndef DUNE_PDELAB_LAPLACEDIRICHLETCCFV_HH
3 #define DUNE_PDELAB_LAPLACEDIRICHLETCCFV_HH
5 #include<dune/common/exceptions.hh>
6 #include<dune/common/fvector.hh>
7 #include<dune/geometry/referenceelements.hh>
45 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
47 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
48 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
52 typedef typename LFSU::Traits::FiniteElementType::
53 Traits::LocalBasisType::Traits::DomainFieldType DF;
54 typedef typename LFSU::Traits::FiniteElementType::
55 Traits::LocalBasisType::Traits::RangeFieldType RF;
58 const Dune::FieldVector<DF,IG::dimension-1>&
59 face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
62 RF face_volume = ig.geometry().integrationElement(face_local)
63 *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
66 const Dune::FieldVector<DF,IG::dimension>&
67 inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
68 const Dune::FieldVector<DF,IG::dimension>&
69 outside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.outside()->type()).position(0,0);
72 Dune::FieldVector<DF,IG::dimension>
73 inside_global = ig.inside()->geometry().global(inside_local);
74 Dune::FieldVector<DF,IG::dimension>
75 outside_global = ig.outside()->geometry().global(outside_local);
76 inside_global -= outside_global;
77 RF distance = inside_global.two_norm();
80 r_s.accumulate(lfsu_s,0,(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
81 r_n.accumulate(lfsu_n,0,-(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
86 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
88 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
92 typedef typename LFSU::Traits::FiniteElementType::
93 Traits::LocalBasisType::Traits::DomainFieldType DF;
94 typedef typename LFSU::Traits::FiniteElementType::
95 Traits::LocalBasisType::Traits::RangeFieldType RF;
98 const Dune::FieldVector<DF,IG::dimension-1>&
99 face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
102 RF face_volume = ig.geometry().integrationElement(face_local)
103 *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
106 const Dune::FieldVector<DF,IG::dimension>&
107 inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
110 Dune::FieldVector<DF,IG::dimension>
111 inside_global = ig.inside()->geometry().global(inside_local);
112 Dune::FieldVector<DF,IG::dimension>
113 outside_global = ig.geometry().global(face_local);
114 inside_global -= outside_global;
115 RF distance = inside_global.two_norm();
118 typename G::Traits::DomainType x = ig.geometryInInside().global(face_local);
119 typename G::Traits::RangeType y;
120 g.evaluate(*(ig.inside()),x,y);
123 r_s.accumulate(lfsu_s,0,(x_s(lfsu_s,0)-y[0])*face_volume/distance);
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: laplacedirichletccfv.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: laplacedirichletccfv.hh:87
Definition: laplacedirichletccfv.hh:35
Definition: laplacedirichletccfv.hh:38
const IG & ig
Definition: constraints.hh:147
LaplaceDirichletCCFV(const G &g_)
Definition: laplacedirichletccfv.hh:41
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:13
Definition: laplacedirichletccfv.hh:24
Definition: laplacedirichletccfv.hh:34
Implement jacobian_apply_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:445
Definition: adaptivity.hh:27
Implement jacobian_boundary() based on alpha_boundary()
Definition: defaultimp.hh:251
sparsity pattern generator
Definition: pattern.hh:29
Definition: laplacedirichletccfv.hh:39
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:157
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538