dune-pdelab 2.7-git
Loading...
Searching...
No Matches
linearacousticsparameter.hh
Go to the documentation of this file.
1// -*- tab-width: 8; indent-tabs-mode: nil -*-
2#ifndef DUNE_PDELAB_LOCALOPERATOR_LINEARACOUSTICSPARAMETER_HH
3#define DUNE_PDELAB_LOCALOPERATOR_LINEARACOUSTICSPARAMETER_HH
4
5#include<vector>
6
7#include<dune/common/exceptions.hh>
8#include<dune/common/fvector.hh>
9#include<dune/geometry/type.hh>
10#include<dune/geometry/referenceelements.hh>
11#include<dune/geometry/quadraturerules.hh>
15
16namespace Dune {
17 namespace PDELab {
18
25 template<typename GV, typename RF>
27 {
29 typedef GV GridViewType;
30
32 enum {
34 dimDomain = GV::dimension
35 };
36
38 typedef typename GV::Grid::ctype DomainFieldType;
39
41 typedef Dune::FieldVector<DomainFieldType,dimDomain> DomainType;
42
44 typedef Dune::FieldVector<DomainFieldType,dimDomain-1> IntersectionDomainType;
45
47 typedef RF RangeFieldType;
48
50 typedef Dune::FieldVector<RF,GV::dimension> RangeType;
51
53 typedef Dune::FieldVector<RF,GV::dimension+1> StateType;
54
56 typedef typename GV::Traits::template Codim<0>::Entity ElementType;
57 typedef typename GV::Intersection IntersectionType;
58 };
59
60 template<typename T>
62 : public Dune::PDELab::GridFunctionBase<Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
63 typename T::Traits::RangeFieldType,
64 T::Traits::dimDomain+1,Dune::FieldVector<typename T::Traits::RangeFieldType,T::Traits::dimDomain+1> >
65 ,LinearAcousticsInitialValueAdapter<T> >
66 {
67 public:
68 typedef Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
69 typename T::Traits::RangeFieldType,
70 T::Traits::dimDomain+1,Dune::FieldVector<typename T::Traits::RangeFieldType,T::Traits::dimDomain+1> > Traits;
71
73 LinearAcousticsInitialValueAdapter (const typename Traits::GridViewType& g_, const T& t_) : g(g_), t(t_) {}
74
76 inline void evaluate (const typename Traits::ElementType& e,
77 const typename Traits::DomainType& x,
78 typename Traits::RangeType& y) const
79 {
80 y = t.u0(e,x);
81 }
82
83 inline const typename Traits::GridViewType& getGridView () const
84 {
85 return g;
86 }
87
88 private:
89 typename Traits::GridViewType g;
90 const T& t;
91 };
92
93 template<typename GV, typename RF>
95 {
96 public:
98
100 : pi(3.141592653589793238462643), time(0.0)
101 {
102 }
103
106 c (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
107 {
108 return 340.0;
109 }
110
112 typename Traits::StateType
113 g (const typename Traits::IntersectionType& is, const typename Traits::IntersectionDomainType& x, const typename Traits::StateType& s) const
114 {
115 typename Traits::DomainType xglobal = is.geometry().global(x);
116 if (xglobal[0]<1e-6)
117 {
118 typename Traits::StateType u(0.0);
119 u[0] = s[0];
120 u[1] = 1.224*(1+0.5*sin(2*pi*1500.0*time));
121 return u;
122 }
123 if (xglobal[0]>1.0-1e-6)
124 {
125 typename Traits::StateType u(0.0);
126 return u;
127 }
128 typename Traits::StateType u(0.0);
129 u[0] = s[0];
130 u[2] = 0.0;
131 return u;
132 }
133
135 typename Traits::StateType
136 q (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
137 {
138 typename Traits::StateType rhs(0.0);
139 return rhs;
140 }
141
143 typename Traits::StateType
144 u0 (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
145 {
146 typename Traits::StateType u(0.0);
147 return u;
148 }
149
151 void setTime (RF t)
152 {
153 time = t;
154 }
155
156 private:
157 double pi;
158 RF time;
159 };
160
161
162
163 }
164}
165
166#endif // DUNE_PDELAB_LOCALOPERATOR_LINEARACOUSTICSPARAMETER_HH
const std::string s
Definition: function.hh:843
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
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
Traits class for linear acoustics parameters.
Definition: linearacousticsparameter.hh:27
GV::Traits::template Codim< 0 >::Entity ElementType
grid types
Definition: linearacousticsparameter.hh:56
Dune::FieldVector< DomainFieldType, dimDomain > DomainType
domain type
Definition: linearacousticsparameter.hh:41
GV::Intersection IntersectionType
Definition: linearacousticsparameter.hh:57
Dune::FieldVector< RF, GV::dimension > RangeType
range type
Definition: linearacousticsparameter.hh:50
GV::Grid::ctype DomainFieldType
Export type for domain field.
Definition: linearacousticsparameter.hh:38
@ dimDomain
dimension of the domain
Definition: linearacousticsparameter.hh:34
RF RangeFieldType
Export type for range field.
Definition: linearacousticsparameter.hh:47
Dune::FieldVector< RF, GV::dimension+1 > StateType
range type
Definition: linearacousticsparameter.hh:53
Dune::FieldVector< DomainFieldType, dimDomain-1 > IntersectionDomainType
domain type
Definition: linearacousticsparameter.hh:44
GV GridViewType
the grid view
Definition: linearacousticsparameter.hh:29
Definition: linearacousticsparameter.hh:66
const Traits::GridViewType & getGridView() const
Definition: linearacousticsparameter.hh:83
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Evaluate the GridFunction at given position.
Definition: linearacousticsparameter.hh:76
Dune::PDELab::GridFunctionTraits< typename T::Traits::GridViewType, typename T::Traits::RangeFieldType, T::Traits::dimDomain+1, Dune::FieldVector< typename T::Traits::RangeFieldType, T::Traits::dimDomain+1 > > Traits
Definition: linearacousticsparameter.hh:70
LinearAcousticsInitialValueAdapter(const typename Traits::GridViewType &g_, const T &t_)
constructor
Definition: linearacousticsparameter.hh:73
Definition: linearacousticsparameter.hh:95
void setTime(RF t)
set time for subsequent evaluation
Definition: linearacousticsparameter.hh:151
Traits::StateType g(const typename Traits::IntersectionType &is, const typename Traits::IntersectionDomainType &x, const typename Traits::StateType &s) const
Dirichlet boundary condition value.
Definition: linearacousticsparameter.hh:113
LinearAcousticsParameterTraits< GV, RF > Traits
Definition: linearacousticsparameter.hh:97
LinearAcousticsModelProblem()
Definition: linearacousticsparameter.hh:99
Traits::StateType q(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
right hand side
Definition: linearacousticsparameter.hh:136
Traits::StateType u0(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
initial value
Definition: linearacousticsparameter.hh:144
Traits::RangeFieldType c(const typename Traits::ElementType &e, const typename Traits::DomainType &x) const
speed of sound
Definition: linearacousticsparameter.hh:106