2 #ifndef DUNE_PDELAB_DIFFUSIONCCFV_HH
3 #define DUNE_PDELAB_DIFFUSIONCCFV_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>
30 template<
typename K,
typename A0,
typename F,
typename B,
typename J,
typename G>
55 DiffusionCCFV (
const K& k_,
const A0& a0_,
const F& f_,
const B& b_,
const J& j_,
const G& g_)
56 : k(k_), a0(a0_), f(f_), b(b_), j(j_), g(g_)
61 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
62 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
65 typedef typename LFSU::Traits::FiniteElementType::
66 Traits::LocalBasisType::Traits::DomainFieldType DF;
69 const int dim = EG::Geometry::dimension;
72 const Dune::FieldVector<DF,dim>&
73 inside_local = Dune::ReferenceElements<DF,dim>::general(eg.entity().type()).position(0,0);
76 typename A0::Traits::RangeType a0value;
77 a0.evaluate(eg.entity(),inside_local,a0value);
79 r.accumulate(lfsu,0,a0value*x(lfsu,0)*eg.geometry().volume());
82 typename F::Traits::RangeType fvalue;
83 f.evaluate(eg.entity(),inside_local,fvalue);
85 r.accumulate(lfsu,0,-fvalue*eg.geometry().volume());
90 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
92 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
93 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
97 typedef typename LFSU::Traits::FiniteElementType::
98 Traits::LocalBasisType::Traits::DomainFieldType DF;
99 typedef typename LFSU::Traits::FiniteElementType::
100 Traits::LocalBasisType::Traits::RangeFieldType RF;
103 const Dune::FieldVector<DF,IG::dimension-1>&
104 face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
107 RF face_volume = ig.geometry().integrationElement(face_local)
108 *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
111 const Dune::FieldVector<DF,IG::dimension>&
112 inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
113 const Dune::FieldVector<DF,IG::dimension>&
114 outside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.outside()->type()).position(0,0);
117 typename K::Traits::RangeType k_inside, k_outside;
118 k.evaluate(*(ig.inside()),inside_local,k_inside);
119 k.evaluate(*(ig.outside()),outside_local,k_outside);
120 typename K::Traits::RangeType k_avg = 2.0/(1.0/(k_inside+1E-30) + 1.0/(k_outside+1E-30));
123 Dune::FieldVector<DF,IG::dimension>
124 inside_global = ig.inside()->geometry().global(inside_local);
125 Dune::FieldVector<DF,IG::dimension>
126 outside_global = ig.outside()->geometry().global(outside_local);
129 inside_global -= outside_global;
130 RF distance = inside_global.two_norm();
133 r_s.accumulate(lfsu_s,0,k_avg*(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
134 r_n.accumulate(lfsu_n,0,-k_avg*(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
139 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
141 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
145 typedef typename LFSU::Traits::FiniteElementType::
146 Traits::LocalBasisType::Traits::DomainFieldType DF;
147 typedef typename LFSU::Traits::FiniteElementType::
148 Traits::LocalBasisType::Traits::RangeFieldType RF;
151 const Dune::FieldVector<DF,IG::dimension-1>&
152 face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
155 RF face_volume = ig.geometry().integrationElement(face_local)
156 *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
159 const Dune::FieldVector<DF,IG::dimension>&
160 inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
163 typename B::Traits::RangeType bctype;
164 b.evaluate(ig,face_local,bctype);
170 Dune::FieldVector<DF,IG::dimension>
171 inside_global = ig.inside()->geometry().global(inside_local);
172 Dune::FieldVector<DF,IG::dimension>
173 outside_global = ig.geometry().global(face_local);
174 inside_global -= outside_global;
175 RF distance = inside_global.two_norm();
178 typename K::Traits::RangeType k_inside;
179 k.evaluate(*(ig.inside()),inside_local,k_inside);
182 typename G::Traits::DomainType x = ig.geometryInInside().global(face_local);
183 typename G::Traits::RangeType y;
184 g.evaluate(*(ig.inside()),x,y);
187 r_s.accumulate(lfsu_s,0,k_inside*(x_s(lfsu_s,0)-y[0])*face_volume/distance);
195 typename J::Traits::DomainType x = ig.geometryInInside().global(face_local);
196 typename J::Traits::RangeType jvalue;
197 j.evaluate(*(ig.inside()),x,jvalue);
200 r_s.accumulate(lfsu_s,0,jvalue*face_volume);
DiffusionCCFV(const K &k_, const A0 &a0_, const F &f_, const B &b_, const J &j_, const G &g_)
Definition: diffusionccfv.hh:55
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: diffusionccfv.hh:91
sparsity pattern generator
Definition: pattern.hh:30
Implement jacobian_apply_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:445
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: diffusionccfv.hh:140
Default flags for all local operators.
Definition: flags.hh:18
Definition: diffusionccfv.hh:31
Implement jacobian_boundary() based on alpha_boundary()
Definition: defaultimp.hh:251
Implement jacobian_apply_boundary() based on alpha_boundary()
Definition: defaultimp.hh:538
Definition: diffusionccfv.hh:48
static const int dim
Definition: adaptivity.hh:82
Definition: diffusionccfv.hh:50
static bool isDirichlet(Type i)
Test for Dirichlet boundary condition.
Definition: diffusionparam.hh:25
Implement jacobian_volume() based on alpha_volume()
Definition: defaultimp.hh:32
const IG & ig
Definition: common/constraints.hh:146
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: diffusionccfv.hh:62
Definition: diffusionccfv.hh:44
Definition: diffusionccfv.hh:53
Definition: diffusionccfv.hh:45
sparsity pattern generator
Definition: pattern.hh:14
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
const EG & eg
Definition: common/constraints.hh:277
Definition: diffusionccfv.hh:49
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: defaultimp.hh:157
Definition: diffusionccfv.hh:51
Definition: diffusionccfv.hh:52