dune-pdelab 2.7-git
Loading...
Searching...
No Matches
permeability_adapter.hh
Go to the documentation of this file.
1#ifndef DUNE_PDELAB_LOCALOPERATOR_PERMEABILITY_ADAPTER_HH
2#define DUNE_PDELAB_LOCALOPERATOR_PERMEABILITY_ADAPTER_HH
3
5
10template<typename T>
12 : public Dune::PDELab::GridFunctionBase<Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
13 typename T::Traits::RangeFieldType,
14 1,Dune::FieldVector<typename T::Traits::RangeFieldType,1> >
15 ,PermeabilityAdapter<T> >
16{
17public:
18 typedef Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
19 typename T::Traits::RangeFieldType,
20 1,Dune::FieldVector<typename T::Traits::RangeFieldType,1> > Traits;
21
23 PermeabilityAdapter (const typename Traits::GridViewType& g_, T& t_)
24 : g(g_), t(t_)
25 {}
26
28 inline void evaluate (const typename Traits::ElementType& e,
29 const typename Traits::DomainType& x,
30 typename Traits::RangeType& y) const
31 {
32 using namespace std;
33 y = log(t.A(e,x)[0][0]);
34 }
35
36 inline const typename Traits::GridViewType& getGridView () const
37 {
38 return g;
39 }
40
41 inline void setTime(double time_)
42 {
43 t.setTime(time_);
44 }
45
46private:
47 const typename Traits::GridViewType& g;
48 T& t;
49};
50
55template<typename T>
57 : public Dune::PDELab::GridFunctionBase<Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
58 typename T::Traits::RangeFieldType,
59 T::Traits::dimDomain,Dune::FieldVector<typename T::Traits::RangeFieldType,T::Traits::dimDomain> >
60 ,DiagonalPermeabilityAdapter<T> >
61{
62public:
63 typedef Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
64 typename T::Traits::RangeFieldType,
65 T::Traits::dimDomain,Dune::FieldVector<typename T::Traits::RangeFieldType,T::Traits::dimDomain> > Traits;
66
69 : g(g_), t(t_)
70 {}
71
73 inline void evaluate (const typename Traits::ElementType& e,
74 const typename Traits::DomainType& x,
75 typename Traits::RangeType& y) const
76 {
77 using namespace std;
78 for (int i=0; i<T::Traits::dimDomain; i++)
79 y[i] = log(t.A(e,x)[i][i]);
80 }
81
82 inline const typename Traits::GridViewType& getGridView () const
83 {
84 return g;
85 }
86
87 inline void setTime(double time_)
88 {
89 t.setTime(time_);
90 }
91
92private:
93 const typename Traits::GridViewType& g;
94 T& t;
95};
96
97#endif // DUNE_PDELAB_LOCALOPERATOR_PERMEABILITY_ADAPTER_HH
STL namespace.
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:50
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:119
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:116
traits class holding the function signature, same as in local function
Definition: function.hh:183
leaf of a function tree
Definition: function.hh:302
Definition: permeability_adapter.hh:16
void setTime(double time_)
Definition: permeability_adapter.hh:41
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: permeability_adapter.hh:28
PermeabilityAdapter(const typename Traits::GridViewType &g_, T &t_)
constructor
Definition: permeability_adapter.hh:23
const Traits::GridViewType & getGridView() const
Definition: permeability_adapter.hh:36
Dune::PDELab::GridFunctionTraits< typename T::Traits::GridViewType, typename T::Traits::RangeFieldType, 1, Dune::FieldVector< typename T::Traits::RangeFieldType, 1 > > Traits
Definition: permeability_adapter.hh:20
Definition: permeability_adapter.hh:61
void setTime(double time_)
Definition: permeability_adapter.hh:87
DiagonalPermeabilityAdapter(const typename Traits::GridViewType &g_, T &t_)
constructor
Definition: permeability_adapter.hh:68
Dune::PDELab::GridFunctionTraits< typename T::Traits::GridViewType, typename T::Traits::RangeFieldType, T::Traits::dimDomain, Dune::FieldVector< typename T::Traits::RangeFieldType, T::Traits::dimDomain > > Traits
Definition: permeability_adapter.hh:65
const Traits::GridViewType & getGridView() const
Definition: permeability_adapter.hh:82
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: permeability_adapter.hh:73