Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
flow_function.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2007 - 2023 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
23#include <deal.II/base/mutex.h>
24#include <deal.II/base/point.h>
25
27
28namespace Functions
29{
47 template <int dim>
48 class FlowFunction : public Function<dim>
49 {
50 public:
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
120 virtual std::size_t
121 memory_consumption() const override;
122
123 protected:
128
129 private:
134
138 mutable std::vector<std::vector<double>> aux_values;
139
143 mutable std::vector<std::vector<Tensor<1, dim>>> aux_gradients;
144 };
145
153 template <int dim>
154 class PoisseuilleFlow : public FlowFunction<dim>
155 {
156 public:
161 PoisseuilleFlow(const double r, const double Re);
162
163 virtual ~PoisseuilleFlow() override = default;
164
165 virtual void
166 vector_values(const std::vector<Point<dim>> & points,
167 std::vector<std::vector<double>> &values) const override;
168 virtual void
170 const std::vector<Point<dim>> & points,
171 std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
172 virtual void
173 vector_laplacians(const std::vector<Point<dim>> & points,
174 std::vector<std::vector<double>> &values) const override;
175
176 private:
177 const double inv_sqr_radius;
178 const double Reynolds;
179 };
180
181
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
247 {
248 public:
251
252 virtual void
253 vector_values(const std::vector<Point<2>> & points,
254 std::vector<std::vector<double>> &values) const override;
255 virtual void
257 const std::vector<Point<2>> & points,
258 std::vector<std::vector<Tensor<1, 2>>> &gradients) const override;
259 virtual void
260 vector_laplacians(const std::vector<Point<2>> & points,
261 std::vector<std::vector<double>> &values) const override;
262
263 private:
265 double
266 Psi(double phi) const;
268 double
269 Psi_1(double phi) const;
271 double
272 Psi_2(double phi) const;
274 double
275 Psi_3(double phi) const;
277 double
278 Psi_4(double phi) const;
280 const double omega;
283 static const double lambda;
285 const double coslo;
287 const double lp;
289 const double lm;
290 };
291
299 class Kovasznay : public FlowFunction<2>
300 {
301 public:
309 Kovasznay(const double Re, bool Stokes = false);
310
311 virtual ~Kovasznay() override = default;
312
313 virtual void
314 vector_values(const std::vector<Point<2>> & points,
315 std::vector<std::vector<double>> &values) const override;
316 virtual void
318 const std::vector<Point<2>> & points,
319 std::vector<std::vector<Tensor<1, 2>>> &gradients) const override;
320 virtual void
321 vector_laplacians(const std::vector<Point<2>> & points,
322 std::vector<std::vector<double>> &values) const override;
323
325 double
326 lambda() const;
327
328 private:
329 const double Reynolds;
330 double lbda;
331 double p_average;
332 const bool stokes;
333 };
334
335} // namespace Functions
336
338
339#endif
virtual void vector_value_list(const std::vector< Point< dim > > &points, std::vector< Vector< double > > &values) const override
std::vector< std::vector< Tensor< 1, dim > > > aux_gradients
void pressure_adjustment(double p)
virtual void vector_gradient_list(const std::vector< Point< dim > > &points, std::vector< std::vector< Tensor< 1, dim > > > &gradients) const override
virtual void vector_laplacian_list(const std::vector< Point< dim > > &points, std::vector< Vector< double > > &values) const override
virtual void vector_laplacians(const std::vector< Point< dim > > &points, std::vector< std::vector< double > > &values) const =0
virtual void vector_gradients(const std::vector< Point< dim > > &points, std::vector< std::vector< Tensor< 1, dim > > > &gradients) const override=0
virtual void vector_values(const std::vector< Point< dim > > &points, std::vector< std::vector< double > > &values) const override=0
virtual std::size_t memory_consumption() const override
virtual ~FlowFunction() override=default
virtual void vector_value(const Point< dim > &points, Vector< double > &value) const override
std::vector< std::vector< double > > aux_values
virtual double value(const Point< dim > &points, const unsigned int component) const override
virtual ~Kovasznay() override=default
virtual void vector_values(const std::vector< Point< 2 > > &points, std::vector< std::vector< double > > &values) const override
virtual void vector_gradients(const std::vector< Point< 2 > > &points, std::vector< std::vector< Tensor< 1, 2 > > > &gradients) const override
virtual void vector_laplacians(const std::vector< Point< 2 > > &points, std::vector< std::vector< double > > &values) const override
double lambda() const
The value of lambda.
virtual void vector_gradients(const std::vector< Point< dim > > &points, std::vector< std::vector< Tensor< 1, dim > > > &gradients) const override
virtual void vector_values(const std::vector< Point< dim > > &points, std::vector< std::vector< double > > &values) const override
virtual void vector_laplacians(const std::vector< Point< dim > > &points, std::vector< std::vector< double > > &values) const override
virtual ~PoisseuilleFlow() override=default
double reaction
The reaction parameter.
virtual void vector_values(const std::vector< Point< dim > > &points, std::vector< std::vector< double > > &values) const override
virtual ~StokesCosine() override=default
double viscosity
The viscosity.
void set_parameters(const double viscosity, const double reaction)
virtual void vector_gradients(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 override
virtual void vector_gradients(const std::vector< Point< 2 > > &points, std::vector< std::vector< Tensor< 1, 2 > > > &gradients) const override
virtual void vector_values(const std::vector< Point< 2 > > &points, std::vector< std::vector< double > > &values) const override
double Psi_1(double phi) const
The derivative of Psi()
double Psi(double phi) const
The auxiliary function Psi.
const double lp
Auxiliary variable 1+lambda.
const double lm
Auxiliary variable 1-lambda.
const double omega
The angle of the reentrant corner, set to 3*pi/2.
virtual void vector_laplacians(const std::vector< Point< 2 > > &points, std::vector< std::vector< double > > &values) const override
double Psi_3(double phi) const
The 3rd derivative of Psi()
StokesLSingularity()
Constructor setting up some data.
double Psi_4(double phi) const
The 4th derivative of Psi()
double Psi_2(double phi) const
The 2nd derivative of Psi()
const double coslo
Cosine of lambda times omega.
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473