Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
flow_function.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2007 - 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 #ifndef dealii_flow_function_h
17 #define dealii_flow_function_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/function.h>
23 #include <deal.II/base/point.h>
24 #include <deal.II/base/thread_management.h>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
28 namespace Functions
29 {
48  template <int dim>
49  class FlowFunction : public Function<dim>
50  {
51  public:
55  FlowFunction();
56 
60  virtual ~FlowFunction() override = default;
61 
66  void
67  pressure_adjustment(double p);
68 
74  virtual void
75  vector_values(const std::vector<Point<dim>> & points,
76  std::vector<std::vector<double>> &values) const override = 0;
82  virtual void
84  const std::vector<Point<dim>> & points,
85  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override = 0;
94  virtual void
95  vector_laplacians(const std::vector<Point<dim>> & points,
96  std::vector<std::vector<double>> &values) const = 0;
97 
98  virtual void
99  vector_value(const Point<dim> &points,
100  Vector<double> & value) const override;
101  virtual double
102  value(const Point<dim> & points,
103  const unsigned int component) const override;
104  virtual void
105  vector_value_list(const std::vector<Point<dim>> &points,
106  std::vector<Vector<double>> & values) const override;
107  virtual void
108  vector_gradient_list(
109  const std::vector<Point<dim>> & points,
110  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
114  virtual void
115  vector_laplacian_list(const std::vector<Point<dim>> &points,
116  std::vector<Vector<double>> & values) const override;
117 
118  std::size_t
119  memory_consumption() const;
120 
121  protected:
126 
127  private:
132 
136  mutable std::vector<std::vector<double>> aux_values;
137 
141  mutable std::vector<std::vector<Tensor<1, dim>>> aux_gradients;
142  };
143 
152  template <int dim>
153  class PoisseuilleFlow : public FlowFunction<dim>
154  {
155  public:
160  PoisseuilleFlow<dim>(const double r, const double Re);
161 
162  virtual ~PoisseuilleFlow() override = default;
163 
164  virtual void
165  vector_values(const std::vector<Point<dim>> & points,
166  std::vector<std::vector<double>> &values) const override;
167  virtual void
169  const std::vector<Point<dim>> & points,
170  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
171  virtual void
172  vector_laplacians(const std::vector<Point<dim>> & points,
173  std::vector<std::vector<double>> &values) const override;
174 
175  private:
176  const double radius;
177  const double Reynolds;
178  };
179 
180 
194  template <int dim>
195  class StokesCosine : public FlowFunction<dim>
196  {
197  public:
202  StokesCosine(const double viscosity = 1., const double reaction = 0.);
206  void
207  set_parameters(const double viscosity, const double reaction);
208 
209  virtual ~StokesCosine() override = default;
210 
211  virtual void
212  vector_values(const std::vector<Point<dim>> & points,
213  std::vector<std::vector<double>> &values) const override;
214  virtual void
216  const std::vector<Point<dim>> & points,
217  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
218  virtual void
219  vector_laplacians(const std::vector<Point<dim>> & points,
220  std::vector<std::vector<double>> &values) const override;
221 
222  private:
224  double viscosity;
226  double reaction;
227  };
228 
229 
248  {
249  public:
252 
253  virtual void
254  vector_values(const std::vector<Point<2>> & points,
255  std::vector<std::vector<double>> &values) const override;
256  virtual void
257  vector_gradients(
258  const std::vector<Point<2>> & points,
259  std::vector<std::vector<Tensor<1, 2>>> &gradients) const override;
260  virtual void
261  vector_laplacians(const std::vector<Point<2>> & points,
262  std::vector<std::vector<double>> &values) const override;
263 
264  private:
266  double
267  Psi(double phi) const;
269  double
270  Psi_1(double phi) const;
272  double
273  Psi_2(double phi) const;
275  double
276  Psi_3(double phi) const;
278  double
279  Psi_4(double phi) const;
281  const double omega;
284  static const double lambda;
286  const double coslo;
288  const double lp;
290  const double lm;
291  };
292 
301  class Kovasznay : public FlowFunction<2>
302  {
303  public:
311  Kovasznay(const double Re, bool Stokes = false);
312 
313  virtual ~Kovasznay() override = default;
314 
315  virtual void
316  vector_values(const std::vector<Point<2>> & points,
317  std::vector<std::vector<double>> &values) const override;
318  virtual void
319  vector_gradients(
320  const std::vector<Point<2>> & points,
321  std::vector<std::vector<Tensor<1, 2>>> &gradients) const override;
322  virtual void
323  vector_laplacians(const std::vector<Point<2>> & points,
324  std::vector<std::vector<double>> &values) const override;
325 
327  double
328  lambda() const;
329 
330  private:
331  const double Reynolds;
332  double lbda;
333  double p_average;
334  const bool stokes;
335  };
336 
337 } // namespace Functions
338 
339 DEAL_II_NAMESPACE_CLOSE
340 
341 #endif
std::vector< std::vector< double > > aux_values
virtual void vector_laplacians(const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const override
const double lm
Auxiliary variable 1-lambda.
virtual void vector_laplacians(const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const override
virtual void vector_gradients(const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim >>> &gradients) const override
const double lp
Auxiliary variable 1+lambda.
virtual void vector_laplacians(const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const =0
virtual void vector_value_list(const std::vector< Point< dim >> &points, std::vector< Vector< double >> &values) const override
static const double lambda
StokesCosine(const double viscosity=1., const double reaction=0.)
virtual void vector_values(const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const override=0
virtual void vector_gradients(const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim >>> &gradients) const override=0
void set_parameters(const double viscosity, const double reaction)
double viscosity
The viscosity.
double Psi(double phi) const
The auxiliary function Psi.
double lambda() const
The value of lambda.
virtual void vector_values(const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const override
double reaction
The reaction parameter.
std::vector< std::vector< Tensor< 1, dim > > > aux_gradients
virtual void vector_value(const Point< dim > &points, Vector< double > &value) const override
const double omega
The angle of the reentrant corner, set to 3*pi/2.
virtual void vector_values(const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const override
virtual void vector_laplacian_list(const std::vector< Point< dim >> &points, std::vector< Vector< double >> &values) const override
double Psi_1(double phi) const
The derivative of Psi()
virtual ~FlowFunction() override=default
const double coslo
Cosine of lambda times omega.
double Psi_4(double phi) const
The 4th derivative of Psi()
virtual double value(const Point< dim > &points, const unsigned int component) const override
double Psi_3(double phi) const
The 3rd derivative of Psi()
StokesLSingularity()
Constructor setting up some data.
virtual void vector_gradients(const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim >>> &gradients) const override
Kovasznay(const double Re, bool Stokes=false)
double Psi_2(double phi) const
The 2nd derivative of Psi()
void pressure_adjustment(double p)