Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
polynomials_rannacher_turek.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2015 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_polynomials_rannacher_turek_h
18 #define dealii_polynomials_rannacher_turek_h
19 
20 #include <deal.II/base/point.h>
21 #include <deal.II/base/tensor.h>
22 
23 #include <vector>
24 
25 DEAL_II_NAMESPACE_OPEN
26 
27 
40 template <int dim>
42 {
43 public:
47  static const unsigned int dimension = dim;
48 
53 
57  double
58  compute_value(const unsigned int i, const Point<dim> &p) const;
59 
65  template <int order>
67  compute_derivative(const unsigned int i, const Point<dim> &p) const;
68 
73  compute_grad(const unsigned int i, const Point<dim> &p) const;
74 
79  compute_grad_grad(const unsigned int i, const Point<dim> &p) const;
80 
87  void
88  compute(const Point<dim> & unit_point,
89  std::vector<double> & values,
90  std::vector<Tensor<1, dim>> &grads,
91  std::vector<Tensor<2, dim>> &grad_grads,
92  std::vector<Tensor<3, dim>> &third_derivatives,
93  std::vector<Tensor<4, dim>> &fourth_derivatives) const;
94 };
95 
96 
97 namespace internal
98 {
99  namespace PolynomialsRannacherTurekImplementation
100  {
101  template <int order, int dim>
102  inline Tensor<order, dim>
103  compute_derivative(const unsigned int, const Point<dim> &)
104  {
105  Assert(dim == 2, ExcNotImplemented());
106  return Tensor<order, dim>();
107  }
108 
109 
110  template <int order>
111  inline Tensor<order, 2>
112  compute_derivative(const unsigned int i, const Point<2> &p)
113  {
114  const unsigned int dim = 2;
115 
116  Tensor<order, dim> derivative;
117  switch (order)
118  {
119  case 1:
120  {
121  Tensor<1, dim> &grad =
122  *reinterpret_cast<Tensor<1, dim> *>(&derivative);
123  if (i == 0)
124  {
125  grad[0] = -2.5 + 3 * p(0);
126  grad[1] = 1.5 - 3 * p(1);
127  }
128  else if (i == 1)
129  {
130  grad[0] = -0.5 + 3.0 * p(0);
131  grad[1] = 1.5 - 3.0 * p(1);
132  }
133  else if (i == 2)
134  {
135  grad[0] = 1.5 - 3.0 * p(0);
136  grad[1] = -2.5 + 3.0 * p(1);
137  }
138  else if (i == 3)
139  {
140  grad[0] = 1.5 - 3.0 * p(0);
141  grad[1] = -0.5 + 3.0 * p(1);
142  }
143  else
144  {
145  Assert(false, ExcNotImplemented());
146  }
147  return derivative;
148  }
149  case 2:
150  {
151  Tensor<2, dim> &grad_grad =
152  *reinterpret_cast<Tensor<2, dim> *>(&derivative);
153  if (i == 0)
154  {
155  grad_grad[0][0] = 3;
156  grad_grad[0][1] = 0;
157  grad_grad[1][0] = 0;
158  grad_grad[1][1] = -3;
159  }
160  else if (i == 1)
161  {
162  grad_grad[0][0] = 3;
163  grad_grad[0][1] = 0;
164  grad_grad[1][0] = 0;
165  grad_grad[1][1] = -3;
166  }
167  else if (i == 2)
168  {
169  grad_grad[0][0] = -3;
170  grad_grad[0][1] = 0;
171  grad_grad[1][0] = 0;
172  grad_grad[1][1] = 3;
173  }
174  else if (i == 3)
175  {
176  grad_grad[0][0] = -3;
177  grad_grad[0][1] = 0;
178  grad_grad[1][0] = 0;
179  grad_grad[1][1] = 3;
180  }
181  return derivative;
182  }
183  default:
184  {
185  // higher derivatives are all zero
186  return Tensor<order, dim>();
187  }
188  }
189  }
190  } // namespace PolynomialsRannacherTurekImplementation
191 } // namespace internal
192 
193 
194 
195 // template functions
196 template <int dim>
197 template <int order>
200  const Point<dim> & p) const
201 {
202  return internal::PolynomialsRannacherTurekImplementation::compute_derivative<
203  order>(i, p);
204 }
205 
206 
207 DEAL_II_NAMESPACE_CLOSE
208 
209 #endif
Tensor< 1, dim > compute_grad(const unsigned int i, const Point< dim > &p) const
static const unsigned int dimension
double compute_value(const unsigned int i, const Point< dim > &p) const
#define Assert(cond, exc)
Definition: exceptions.h:1407
Tensor< order, dim > compute_derivative(const unsigned int i, const Point< dim > &p) const
Tensor< 2, dim > compute_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: mpi.h:90
static ::ExceptionBase & ExcNotImplemented()
void compute(const Point< dim > &unit_point, std::vector< double > &values, std::vector< Tensor< 1, dim >> &grads, std::vector< Tensor< 2, dim >> &grad_grads, std::vector< Tensor< 3, dim >> &third_derivatives, std::vector< Tensor< 4, dim >> &fourth_derivatives) const