dune-pdelab  2.4.1
laplacedirichletccfv.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_LAPLACEDIRICHLETCCFV_HH
3 #define DUNE_PDELAB_LAPLACEDIRICHLETCCFV_HH
4 #warning This file is deprecated and will be removed after the Dune-PDELab 2.4 release! Use the ConvectionDiffusionCCFV local operator from dune/pdelab/localoperator/convectiondiffusionccfv.hh instead!
5 
6 #include<dune/common/exceptions.hh>
7 #include<dune/common/fvector.hh>
8 #include<dune/geometry/referenceelements.hh>
9 
11 
12 #include"pattern.hh"
13 #include"flags.hh"
14 
15 
16 namespace Dune {
17  namespace PDELab {
18 
19  // a local operator for solving the Laplace equation with Dirichlet boundary conditions
20  // - \Delta u = 0 in \Omega,
21  // u = g on \partial\Omega
22  // with cell centered finite volumes on axiparallel cube grids
23  // G : grid function for Dirichlet boundary conditions
24  template<typename G>
25  class LaplaceDirichletCCFV : public NumericalJacobianApplySkeleton<LaplaceDirichletCCFV<G> >,
26  public NumericalJacobianApplyBoundary<LaplaceDirichletCCFV<G> >,
27  public NumericalJacobianSkeleton<LaplaceDirichletCCFV<G> >,
28  public NumericalJacobianBoundary<LaplaceDirichletCCFV<G> >,
29  public FullSkeletonPattern,
30  public FullVolumePattern,
32  {
33  public:
34  // pattern assembly flags
35  enum { doPatternVolume = true };
36  enum { doPatternSkeleton = true };
37 
38  // residual assembly flags
39  enum { doAlphaSkeleton = true };
40  enum { doAlphaBoundary = true };
41 
42  DUNE_DEPRECATED_MSG("Deprecated in Dune-PDELab 2.4, use the local operator ConvectionDiffusionCCFV instead!")
43  LaplaceDirichletCCFV (const G& g_) : g(g_) {}
44 
45  // skeleton integral depending on test and ansatz functions
46  // each face is only visited ONCE!
47  template<typename IG, typename LFSU, typename X, typename LFSV, typename R>
48  void alpha_skeleton (const IG& ig,
49  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
50  const LFSU& lfsu_n, const X& x_n, const LFSV& lfsv_n,
51  R& r_s, R& r_n) const
52  {
53  // domain and range field type
54  typedef typename LFSU::Traits::FiniteElementType::
55  Traits::LocalBasisType::Traits::DomainFieldType DF;
56  typedef typename LFSU::Traits::FiniteElementType::
57  Traits::LocalBasisType::Traits::RangeFieldType RF;
58 
59  // center in face's reference element
60  const Dune::FieldVector<DF,IG::dimension-1>&
61  face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
62 
63  // face volume for integration
64  RF face_volume = ig.geometry().integrationElement(face_local)
65  *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
66 
67  // cell centers in references elements
68  const Dune::FieldVector<DF,IG::dimension>&
69  inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
70  const Dune::FieldVector<DF,IG::dimension>&
71  outside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.outside()->type()).position(0,0);
72 
73  // distance between the two cell centers
74  Dune::FieldVector<DF,IG::dimension>
75  inside_global = ig.inside()->geometry().global(inside_local);
76  Dune::FieldVector<DF,IG::dimension>
77  outside_global = ig.outside()->geometry().global(outside_local);
78  inside_global -= outside_global;
79  RF distance = inside_global.two_norm();
80 
81  // contribution to residual on inside element, other residual is computed by symmetric call
82  r_s.accumulate(lfsu_s,0,(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
83  r_n.accumulate(lfsu_n,0,-(x_s(lfsu_s,0)-x_n(lfsu_n,0))*face_volume/distance);
84  }
85 
86  // skeleton integral depending on test and ansatz functions
87  // We put the Dirchlet evaluation also in the alpha term to save som e geometry evaluations
88  template<typename IG, typename LFSU, typename X, typename LFSV, typename R>
89  void alpha_boundary (const IG& ig,
90  const LFSU& lfsu_s, const X& x_s, const LFSV& lfsv_s,
91  R& r_s) const
92  {
93  // domain and range field type
94  typedef typename LFSU::Traits::FiniteElementType::
95  Traits::LocalBasisType::Traits::DomainFieldType DF;
96  typedef typename LFSU::Traits::FiniteElementType::
97  Traits::LocalBasisType::Traits::RangeFieldType RF;
98 
99  // center in face's reference element
100  const Dune::FieldVector<DF,IG::dimension-1>&
101  face_local = Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).position(0,0);
102 
103  // face volume for integration
104  RF face_volume = ig.geometry().integrationElement(face_local)
105  *Dune::ReferenceElements<DF,IG::dimension-1>::general(ig.geometry().type()).volume();
106 
107  // cell center in reference element
108  const Dune::FieldVector<DF,IG::dimension>&
109  inside_local = Dune::ReferenceElements<DF,IG::dimension>::general(ig.inside()->type()).position(0,0);
110 
111  // distance between cell center and face center
112  Dune::FieldVector<DF,IG::dimension>
113  inside_global = ig.inside()->geometry().global(inside_local);
114  Dune::FieldVector<DF,IG::dimension>
115  outside_global = ig.geometry().global(face_local);
116  inside_global -= outside_global;
117  RF distance = inside_global.two_norm();
118 
119  // evaluate boundary condition function
120  typename G::Traits::DomainType x = ig.geometryInInside().global(face_local);
121  typename G::Traits::RangeType y;
122  g.evaluate(*(ig.inside()),x,y);
123 
124  // contribution to residual on inside element
125  r_s.accumulate(lfsu_s,0,(x_s(lfsu_s,0)-y[0])*face_volume/distance);
126  }
127 
128  private:
129  const G& g;
130  };
131 
133  } // namespace PDELab
134 } // namespace Dune
135 
136 #endif
sparsity pattern generator
Definition: pattern.hh:13
Implements linear and nonlinear versions of jacobian_apply_skeleton() based on alpha_skeleton() ...
Definition: numericaljacobianapply.hh:180
const IG & ig
Definition: constraints.hh:148
Implement jacobian_boundary() based on alpha_boundary()
Definition: numericaljacobian.hh:250
Definition: adaptivity.hh:27
sparsity pattern generator
Definition: pattern.hh:29
Implements linear and nonlinear versions of jacobian_apply_boundary() based on alpha_boundary() ...
Definition: numericaljacobianapply.hh:285
Definition: laplacedirichletccfv.hh:40
Default flags for all local operators.
Definition: flags.hh:18
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:89
Definition: laplacedirichletccfv.hh:25
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:48
Definition: laplacedirichletccfv.hh:36
Definition: laplacedirichletccfv.hh:35
Definition: laplacedirichletccfv.hh:39
Implement jacobian_skeleton() based on alpha_skeleton()
Definition: numericaljacobian.hh:156