dune-localfunctions  2.9.0
p0localbasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 // SPDX-FileCopyrightInfo: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
4 // SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5 #ifndef DUNE_P0LOCALBASIS_HH
6 #define DUNE_P0LOCALBASIS_HH
7 
8 #include <numeric>
9 
10 #include <dune/common/fmatrix.hh>
11 
13 
14 namespace Dune
15 {
28  template<class D, class R, int d>
30  {
31  public:
33  typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,
34  Dune::FieldMatrix<R,1,d> > Traits;
35 
37  unsigned int size () const
38  {
39  return 1;
40  }
41 
43  inline void evaluateFunction (const typename Traits::DomainType&,
44  std::vector<typename Traits::RangeType>& out) const
45  {
46  out.resize(1);
47  out[0] = 1;
48  }
49 
51  inline void
52  evaluateJacobian (const typename Traits::DomainType&, // position
53  std::vector<typename Traits::JacobianType>& out) const // return value
54  {
55  out.resize(1);
56  for (int i=0; i<d; i++)
57  out[0][0][i] = 0;
58  }
59 
65  void partial(const std::array<unsigned int,d>& order,
66  const typename Traits::DomainType& in,
67  std::vector<typename Traits::RangeType>& out) const
68  {
69  auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
70  if (totalOrder == 0) {
71  evaluateFunction(in, out);
72  } else {
73  out.resize(1);
74  out[0] = 0;
75  }
76  }
77 
79  unsigned int order () const
80  {
81  return 0;
82  }
83  };
84 
85 }
86 
87 #endif
Definition: bdfmcube.hh:18
Type traits for LocalBasisVirtualInterface.
Definition: common/localbasis.hh:34
D DomainType
domain type
Definition: common/localbasis.hh:42
Constant shape function.
Definition: p0localbasis.hh:30
unsigned int order() const
Polynomial order of the shape functions.
Definition: p0localbasis.hh:79
void evaluateJacobian(const typename Traits::DomainType &, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: p0localbasis.hh:52
void evaluateFunction(const typename Traits::DomainType &, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p0localbasis.hh:43
void partial(const std::array< unsigned int, d > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of any order of all shape functions.
Definition: p0localbasis.hh:65
unsigned int size() const
number of shape functions
Definition: p0localbasis.hh:37
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
export type traits for function signature
Definition: p0localbasis.hh:34