bilinearForm1D.hh

Go to the documentation of this file.
1 
6 #ifndef bilinearLinFEM1d_hh
7 #define bilinearLinFEM1d_hh
8 
9 #include <memory>
10 #include "basics/typedefs.hh"
12 #include "formula/formula.hh"
13 #include "operator/bilinearForm.hh"
15 
16 
17 namespace linearFEM {
18 
19  using concepts::Real;
20 
21  class Line;
22 
23  // ************************************************************* Laplace1d **
24 
32  class Laplace1d : public concepts::BilinearForm<Real> {
33  public:
34  virtual void operator()(const concepts::Element<Real>& elmX,
35  const concepts::Element<Real>& elmY,
37  virtual Laplace1d* clone() const { return new Laplace1d(); }
38  };
39 
40  // ************************************************************ Identity1d **
41 
49  class Identity1d : public concepts::BilinearForm<Real> {
50  public:
51  virtual void operator()(const concepts::Element<Real>& elmX,
52  const concepts::Element<Real>& elmY,
54  virtual Identity1d* clone() const { return new Identity1d(); }
55  };
56 
57  // *********************************************************** CIdentity1d **
58 
67  class CIdentity1d : public concepts::BilinearForm<Real> {
68  public:
74  const uint gauss_p = 1) :
75  frm_(frm.clone()), quad_(gauss_p) {}
76  virtual void operator()(const concepts::Element<Real>& elmX,
77  const concepts::Element<Real>& elmY,
79  virtual CIdentity1d* clone() const {
80  return new CIdentity1d(*frm_.get(), quad_.n()); }
81  private:
83  std::unique_ptr<const concepts::Formula<Real> > frm_;
86  };
87 
88 } // namespace linearFEM
89 
90 #endif // bilinearLinFEM1d_hh
Discrete equivalent of the Laplacian in 1D for linear FEM.
virtual void operator()(const concepts::Element< Real > &elmX, const concepts::Element< Real > &elmY, concepts::ElementMatrix< Real > &em) const
virtual void operator()(const concepts::Element< Real > &elmX, const concepts::Element< Real > &elmY, concepts::ElementMatrix< Real > &em) const
uint n() const
Returns the number of quadrature points.
Definition: quadrature.hh:109
CIdentity1d(const concepts::Formula< Real > &frm, const uint gauss_p=1)
Constructor.
Discrete equivalent for a reaction term in 1D for linear FEM.
Abstract function class to evaluate a bilinear form.
Definition: bilinearForm.hh:33
Linear FEM in 1D, 2D and 3D.
Definition: spaceTraits.hh:19
const concepts::Quadrature< 4 > quad_
Quadrature rule.
virtual Identity1d * clone() const
Virtual constructor.
std::unique_ptr< const concepts::Formula< Real > > frm_
Formula.
virtual Laplace1d * clone() const
Virtual constructor.
Discrete equivalent for a reaction term in 1D for linear FEM.
virtual CIdentity1d * clone() const
Virtual constructor.
double Real
Type normally used for a floating point number.
Definition: typedefs.hh:36
virtual void operator()(const concepts::Element< Real > &elmX, const concepts::Element< Real > &elmY, concepts::ElementMatrix< Real > &em) const
Page URL: http://wiki.math.ethz.ch/bin/view/Concepts/WebHome
21 August 2020
© 2020 Eidgenössische Technische Hochschule Zürich