Reference documentation for deal.II version Git d51799cb54 2020-09-28 09:22:08 +0200
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
flow_function.h
Go to the documentation of this file.
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>
25 
27 
28 namespace Functions
29 {
47  template <int dim>
48  class FlowFunction : public Function<dim>
49  {
50  public:
54  FlowFunction();
55 
59  virtual ~FlowFunction() override = default;
60 
65  void
66  pressure_adjustment(double p);
67 
73  virtual void
74  vector_values(const std::vector<Point<dim>> & points,
75  std::vector<std::vector<double>> &values) const override = 0;
81  virtual void
83  const std::vector<Point<dim>> & points,
84  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override = 0;
93  virtual void
94  vector_laplacians(const std::vector<Point<dim>> & points,
95  std::vector<std::vector<double>> &values) const = 0;
96 
97  virtual void
98  vector_value(const Point<dim> &points,
99  Vector<double> & value) const override;
100  virtual double
101  value(const Point<dim> & points,
102  const unsigned int component) const override;
103  virtual void
104  vector_value_list(const std::vector<Point<dim>> &points,
105  std::vector<Vector<double>> & values) const override;
106  virtual void
108  const std::vector<Point<dim>> & points,
109  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
113  virtual void
114  vector_laplacian_list(const std::vector<Point<dim>> &points,
115  std::vector<Vector<double>> & values) const override;
116 
117  std::size_t
118  memory_consumption() const;
119 
120  protected:
125 
126  private:
131 
135  mutable std::vector<std::vector<double>> aux_values;
136 
140  mutable std::vector<std::vector<Tensor<1, dim>>> aux_gradients;
141  };
142 
150  template <int dim>
151  class PoisseuilleFlow : public FlowFunction<dim>
152  {
153  public:
158  PoisseuilleFlow<dim>(const double r, const double Re);
159 
160  virtual ~PoisseuilleFlow() override = default;
161 
162  virtual void
163  vector_values(const std::vector<Point<dim>> & points,
164  std::vector<std::vector<double>> &values) const override;
165  virtual void
167  const std::vector<Point<dim>> & points,
168  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
169  virtual void
170  vector_laplacians(const std::vector<Point<dim>> & points,
171  std::vector<std::vector<double>> &values) const override;
172 
173  private:
174  const double radius;
175  const double Reynolds;
176  };
177 
178 
191  template <int dim>
192  class StokesCosine : public FlowFunction<dim>
193  {
194  public:
199  StokesCosine(const double viscosity = 1., const double reaction = 0.);
203  void
204  set_parameters(const double viscosity, const double reaction);
205 
206  virtual ~StokesCosine() override = default;
207 
208  virtual void
209  vector_values(const std::vector<Point<dim>> & points,
210  std::vector<std::vector<double>> &values) const override;
211  virtual void
213  const std::vector<Point<dim>> & points,
214  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
215  virtual void
216  vector_laplacians(const std::vector<Point<dim>> & points,
217  std::vector<std::vector<double>> &values) const override;
218 
219  private:
221  double viscosity;
223  double reaction;
224  };
225 
226 
244  {
245  public:
248 
249  virtual void
250  vector_values(const std::vector<Point<2>> & points,
251  std::vector<std::vector<double>> &values) const override;
252  virtual void
254  const std::vector<Point<2>> & points,
255  std::vector<std::vector<Tensor<1, 2>>> &gradients) const override;
256  virtual void
257  vector_laplacians(const std::vector<Point<2>> & points,
258  std::vector<std::vector<double>> &values) const override;
259 
260  private:
262  double
263  Psi(double phi) const;
265  double
266  Psi_1(double phi) const;
268  double
269  Psi_2(double phi) const;
271  double
272  Psi_3(double phi) const;
274  double
275  Psi_4(double phi) const;
277  const double omega;
280  static const double lambda;
282  const double coslo;
284  const double lp;
286  const double lm;
287  };
288 
296  class Kovasznay : public FlowFunction<2>
297  {
298  public:
306  Kovasznay(const double Re, bool Stokes = false);
307 
308  virtual ~Kovasznay() override = default;
309 
310  virtual void
311  vector_values(const std::vector<Point<2>> & points,
312  std::vector<std::vector<double>> &values) const override;
313  virtual void
315  const std::vector<Point<2>> & points,
316  std::vector<std::vector<Tensor<1, 2>>> &gradients) const override;
317  virtual void
318  vector_laplacians(const std::vector<Point<2>> & points,
319  std::vector<std::vector<double>> &values) const override;
320 
322  double
323  lambda() const;
324 
325  private:
326  const double Reynolds;
327  double lbda;
328  double p_average;
329  const bool stokes;
330  };
331 
332 } // namespace Functions
333 
335 
336 #endif
const double Reynolds
std::vector< std::vector< double > > aux_values
const double lm
Auxiliary variable 1-lambda.
const double lp
Auxiliary variable 1+lambda.
virtual void vector_gradient_list(const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim >>> &gradients) const override
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
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
double viscosity
The viscosity.
double reaction
The reaction parameter.
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
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_laplacian_list(const std::vector< Point< dim >> &points, std::vector< Vector< double >> &values) const override
virtual ~FlowFunction() override=default
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
const double coslo
Cosine of lambda times omega.
virtual double value(const Point< dim > &points, const unsigned int component) const override
std::size_t memory_consumption() const
void pressure_adjustment(double p)