formula.hh

Go to the documentation of this file.
1 
6 #ifndef formulahp2D_hh
7 #define formulahp2D_hh
8 
9 #include "toolbox/sequence.hh"
10 #include "toolbox/set.hh"
11 #include "toolbox/hashMap.hh"
12 #include "geometry/connectorSet.hh"
14 #include "space/element.hh"
15 #include "space/function.hh"
16 
17 // debugging
18 #include "basics/debug.hh"
19 
20 #define FormulaFromWeightConstr_D 0
21 #define FormulaFromWeightAppl_D 0
22 
23 namespace concepts {
24 
25  template<typename F>
26  class Array;
27 
28  template<typename F>
29  class Vector;
30 
31  class Attribute;
32  class Edge;
33  class Edge2d;
34 
35 } // namespace concepts
36 
37 namespace hp2D {
38 
39  // forward declarations
40  class TraceSpace;
41 
42  template<typename F>
43  class BaseQuad;
44 
45  template<typename F>
46  class Quad;
47 
48  template<typename F>
49  class BaseEdge;
50 
51  using concepts::Real;
52  using concepts::Real2d;
53 
54  // *********************************************************** H1Extension **
55 
66  template<class F>
68  public:
89  const concepts::Set<concepts::Attribute> attributes);
90 
91  virtual ~H1Extension();
92 
94  const Real p, const Real t = 0.0) const;
96  const Real2d& p, const Real t = 0.0) const;
98  const concepts::Real3d& p, const Real t = 0.0) const;
99 
101  virtual H1Extension<F>* clone() const {
102  return new H1Extension<F>(*frm_);
103  }
104  protected:
105  virtual std::ostream& info(std::ostream& os) const;
106  private:
108  std::unique_ptr<const concepts::ElementFormula<F> > frm_;
123 
124  bool compute_(const hp2D::Quad<Real>* quad, const Real2d& p, const Real t,
125  F& val) const;
126 
128  void clear_() const;
129  };
130 
131  // ********************************************** ElementFormulaInterpGrad **
132 
140  template<typename F, uint dim = 2>
142 
143  // ******************************************** ElementFormulaInterpGrad<2> **
144 
145  template<typename F>
147  public concepts::ElementFormula<concepts::Point<F,2> > {
148  public:
150 
151  virtual concepts::Point<F,2> operator()
153  const Real p, const Real t = 0.0) const { return concepts::Point<F,2>(0.0); }
154  virtual concepts::Point<F,2> operator()
156  const concepts::Real2d& p, const Real t = 0.0) const;
157  virtual concepts::Point<F,2> operator()
159  const concepts::Real3d& p, const Real t = 0.0) const { return concepts::Point<F,2>(0.0); }
160 
163  return new ElementFormulaInterpGrad<F,2>(*f_);
164  }
165  protected:
166  virtual std::ostream& info(std::ostream& os) const;
167  private:
169  std::unique_ptr<const concepts::ElementFormula<F> > f_;
173  mutable concepts::Array<Real> x_[2];
175  mutable bool samePoints_;
179  mutable bool zeroElement_;
180  };
181 
182  // ********************************************* ElementFormulaInterpGradN **
183 
192  template<typename F, uint dim = 2>
194 
195  template<typename F>
197  public concepts::ElementFormula<concepts::Point<F,2> > {
198  public:
200 
201  virtual concepts::Point<F,2> operator()
203  const Real p, const Real t = 0.0) const { return concepts::Point<F,2>(0.0); }
204  virtual concepts::Point<F,2> operator()
206  const concepts::Real2d& p, const Real t = 0.0) const;
207  virtual concepts::Point<F,2> operator()
209  const concepts::Real3d& p, const Real t = 0.0) const { return concepts::Point<F,2>(0.0); }
210 
213  return new ElementFormulaInterpGradN<F,2>(*f_);
214  }
215  protected:
216  virtual std::ostream& info(std::ostream& os) const;
217  private:
219  std::unique_ptr<const concepts::ElementFormula<F> > f_;
225  mutable bool zeroElement_;
226  };
227 
228  // ************************************************ ElementFormulaEdgeMean **
229 
235  template<typename F>
237  public:
239  const concepts::Vector<F>& v,
241 
243  const Real p, const Real t = 0.0) const;
245  const concepts::Real2d& p,
246  const Real t = 0.0) const {
247  return (*this)(elm, p[0], t);
248  }
250  const concepts::Real3d& p,
251  const Real t = 0.0) const {
252  return (*this)(elm, p[0], t);
253  }
254 
256  virtual ElementFormulaEdgeMean<F>* clone() const {
257  return new ElementFormulaEdgeMean<F>(spc_, v_, *f_);
258  }
259  protected:
260  virtual std::ostream& info(std::ostream& os) const;
261  private:
267  std::unique_ptr<const concepts::ElementFunction<F> > f_;
268  };
269 
270  // ************************************************ ElementFormulaEdgeJump **
271 
278  template<typename F>
280  public:
282  const concepts::Vector<F>& v,
284 
286  const Real p, const Real t = 0.0) const;
288  const concepts::Real2d& p,
289  const Real t = 0.0) const {
290  return (*this)(elm, p[0], t);
291  }
293  const concepts::Real3d& p,
294  const Real t = 0.0) const {
295  return (*this)(elm, p[0], t);
296  }
297 
299  virtual ElementFormulaEdgeJump<F>* clone() const {
300  return new ElementFormulaEdgeJump<F>(spc_, v_, *f_);
301  }
302  protected:
303  virtual std::ostream& info(std::ostream& os) const;
304  private:
310  std::unique_ptr<const concepts::ElementFunction<F> > f_;
311  };
312 
313  // **************************************** ElementFormulaSignNormalVector **
314 
321  public:
323  const Real p, const Real t = 0.0) const;
325  const concepts::Real2d& p,
326  const Real t = 0.0) const {
327  return (*this)(elm, p[0], t);
328  }
330  const concepts::Real3d& p,
331  const Real t = 0.0) const {
332  return (*this)(elm, p[0], t);
333  }
334 
337  return new ElementFormulaSignNormalVector();
338  }
339  protected:
340  virtual std::ostream& info(std::ostream& os) const;
341  };
342 
343 
344 } // namespace hp2D
345 
346 #endif // formulahp2D_hh
const concepts::Vector< F > & v_
Coefficient vector.
Definition: formula.hh:265
virtual ElementFormulaEdgeJump< F > * clone() const
Virtual copy constructor.
Definition: formula.hh:299
Continuous extension of a formula on edges, e.g.
Definition: formula.hh:67
virtual ElementFormulaInterpGrad< F, 2 > * clone() const
Virtual copy constructor.
Definition: formula.hh:162
std::unique_ptr< const concepts::ElementFunction< F > > f_
Element function.
Definition: formula.hh:267
const hp2D::TraceSpace & spc_
Trace space.
Definition: formula.hh:306
virtual ElementFormulaInterpGradN< F, 2 > * clone() const
Virtual copy constructor.
Definition: formula.hh:212
std::unique_ptr< const concepts::ElementFunction< F > > f_
Element function.
Definition: formula.hh:310
Builds the trace space of an FE space.
Definition: traces.hh:52
ElementFormulaEdgeMean(const hp2D::TraceSpace &spc, const concepts::Vector< F > &v, const concepts::ElementFunction< F > &f)
Element formula on 1D elements based on Edge2d returning the 1.0 if the normal vector is the right on...
Definition: formula.hh:320
const concepts::ElementWithCell< Real > * quad_
Last element.
Definition: formula.hh:221
virtual H1Extension< F > * clone() const
Virtual copy constructor.
Definition: formula.hh:101
concepts::HashMap< F > vtxValues_
Values on the vertices.
Definition: formula.hh:122
virtual F operator()(const concepts::ElementWithCell< Real > &elm, const Real p, const Real t=0.0) const
virtual F operator()(const concepts::ElementWithCell< Real > &elm, const Real p, const Real t=0.0) const
Point< Real, 2 > Real2d
const concepts::Set< concepts::Attribute > attributes_
Attributes of the cells taken into account.
Definition: formula.hh:112
Element formula for a jump value on edges of elements of an other element formula.
Definition: formula.hh:279
virtual std::ostream & info(std::ostream &os) const
Returns information in an output stream.
virtual ElementFormulaEdgeMean< F > * clone() const
Virtual copy constructor.
Definition: formula.hh:256
ElementFormulaInterpGradN(const concepts::ElementFormula< F > &f)
bool zeroElement_
Element with all zero values.
Definition: formula.hh:225
2D hp-FEM for H1-conforming elements.
const concepts::ElementWithCell< Real > * quad_
Last element.
Definition: formula.hh:171
std::unique_ptr< const concepts::ElementFormula< F > > frm_
The formula on the edges.
Definition: formula.hh:108
Element formula for the gradient of an scalar element formula.
Definition: formula.hh:141
const concepts::Vector< F > & v_
Coefficient vector.
Definition: formula.hh:308
virtual std::ostream & info(std::ostream &os) const
Returns information in an output stream.
An array of objects.
Definition: bilinearForm.hh:23
H1Extension(const concepts::ElementFormula< F > &frm, const concepts::Set< concepts::Attribute > attributes)
Constructor.
virtual std::ostream & info(std::ostream &os) const
Returns information in an output stream.
Interface for a formula defined element by element.
virtual std::ostream & info(std::ostream &os) const
Returns information in an output stream.
Element formula for the gradient of an scalar element formula.
Definition: formula.hh:193
virtual ElementFormulaSignNormalVector * clone() const
Virtual copy constructor.
Definition: formula.hh:336
H1Extension(const concepts::ElementFormula< F > &frm)
Constructor for extension from formula on edges.
virtual ~H1Extension()
concepts::Array< F > coeff_
Coefficients of the formula after projection.
Definition: formula.hh:223
bool compute_(const hp2D::Quad< Real > *quad, const Real2d &p, const Real t, F &val) const
bool zeroElement_
Element with all zero values.
Definition: formula.hh:179
virtual Real operator()(const concepts::ElementWithCell< Real > &elm, const Real p, const Real t=0.0) const
virtual std::ostream & info(std::ostream &os) const
Returns information in an output stream.
void clear_() const
Clear the temporary information (elmEdges_, cellEdges_).
ElementFormulaEdgeJump(const hp2D::TraceSpace &spc, const concepts::Vector< F > &v, const concepts::ElementFunction< F > &f)
Element formula for a mean value on edges of elements of an other element formula.
Definition: formula.hh:236
std::unique_ptr< const concepts::ElementFormula< F > > f_
Element formula.
Definition: formula.hh:169
concepts::HashMap< const concepts::Edge * > edges_
The last used edges.
Definition: formula.hh:120
concepts::Array< concepts::ElementWithCell< Real > * > elmEdges_
The last used elements on the edges of elm_.
Definition: formula.hh:116
virtual F operator()(const concepts::ElementWithCell< Real > &elm, const Real p, const Real t=0.0) const
std::unique_ptr< const concepts::ElementFormula< F > > f_
Element formula.
Definition: formula.hh:219
const concepts::ElementWithCell< Real > * elm_
The last used element.
Definition: formula.hh:114
concepts::Array< F > values_
Values of the element formula in the current quadrilateral.
Definition: formula.hh:177
virtual std::ostream & info(std::ostream &os) const
Returns information in an output stream.
bool samePoints_
Flag, if same number of points in both directions.
Definition: formula.hh:175
concepts::Sequence< concepts::Edge2d * > cellEdges_
List of the generated cells on the edges.
Definition: formula.hh:118
double Real
Type normally used for a floating point number.
Definition: typedefs.hh:36
const hp2D::TraceSpace & spc_
Trace space.
Definition: formula.hh:263
ElementFormulaInterpGrad(const concepts::ElementFormula< F > &f)
Basic namespace for Concepts-2.
Definition: pml_formula.h:16
Page URL: http://wiki.math.ethz.ch/bin/view/Concepts/WebHome
21 August 2020
© 2020 Eidgenössische Technische Hochschule Zürich