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/common/static_assert.hh>
8 #include<dune/geometry/referenceelements.hh>
46 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
48 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
49 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
53 typedef typename LFSU::Traits::FiniteElementType::
54 Traits::LocalBasisType::Traits::DomainFieldType DF;
55 typedef typename LFSU::Traits::FiniteElementType::
56 Traits::LocalBasisType::Traits::RangeFieldType RF;
59 const Dune::FieldVector<DF,IG::dimension-1>&
60 face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
63 RF face_volume = ig.geometry().integrationElement(face_local)
64 *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
67 const Dune::FieldVector<DF,IG::dimension>&
68 inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
69 const Dune::FieldVector<DF,IG::dimension>&
70 outside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.outside()->type()).position(0,0);
73 Dune::FieldVector<DF,IG::dimension>
74 inside_global = ig.inside()->geometry().global(inside_local);
75 Dune::FieldVector<DF,IG::dimension>
76 outside_global = ig.outside()->geometry().global(outside_local);
77 inside_global -= outside_global;
78 RF distance = inside_global.two_norm();
81 r_s.accumulate(lfsu_s,0,(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
82 r_n.accumulate(lfsu_n,0,-(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
87 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
89 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
93 typedef typename LFSU::Traits::FiniteElementType::
94 Traits::LocalBasisType::Traits::DomainFieldType DF;
95 typedef typename LFSU::Traits::FiniteElementType::
96 Traits::LocalBasisType::Traits::RangeFieldType RF;
99 const Dune::FieldVector<DF,IG::dimension-1>&
100 face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
103 RF face_volume = ig.geometry().integrationElement(face_local)
104 *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
107 const Dune::FieldVector<DF,IG::dimension>&
108 inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
111 Dune::FieldVector<DF,IG::dimension>
112 inside_global = ig.inside()->geometry().global(inside_local);
113 Dune::FieldVector<DF,IG::dimension>
114 outside_global = ig.geometry().global(face_local);
115 inside_global -= outside_global;
116 RF distance = inside_global.two_norm();
119 typename G::Traits::DomainType x = ig.geometryInInside().global(face_local);
120 typename G::Traits::RangeType y;
121 g.evaluate(*(ig.inside()),x,y);
124 r_s.accumulate(lfsu_s,0,(x_s(lfsu_s,0)-y[0])*face_volume/distance);
Definition: laplacedirichletccfv.hh:40
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: laplacedirichletccfv.hh:39
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538
const IG & ig
Definition: common/constraints.hh:146
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:47
sparsity pattern generator
Definition: pattern.hh:14
Definition: laplacedirichletccfv.hh:25
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: laplacedirichletccfv.hh:88
Definition: laplacedirichletccfv.hh:35
LaplaceDirichletCCFV(const G &g_)
Definition: laplacedirichletccfv.hh:42
Definition: laplacedirichletccfv.hh:36