From: Marc Fehling Date: Sat, 21 Jan 2023 22:54:16 +0000 (-0700) Subject: Minor changes to PoisseuilleFlow. X-Git-Tag: v9.5.0-rc1~627^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F14707%2Fhead;p=dealii.git Minor changes to PoisseuilleFlow. --- diff --git a/include/deal.II/base/flow_function.h b/include/deal.II/base/flow_function.h index 33b213ca6b..6de3c60247 100644 --- a/include/deal.II/base/flow_function.h +++ b/include/deal.II/base/flow_function.h @@ -174,7 +174,7 @@ namespace Functions std::vector> &values) const override; private: - const double radius; + const double inv_sqr_radius; const double Reynolds; }; diff --git a/source/base/flow_function.cc b/source/base/flow_function.cc index c70b48b86a..08674fa824 100644 --- a/source/base/flow_function.cc +++ b/source/base/flow_function.cc @@ -189,7 +189,7 @@ namespace Functions template PoisseuilleFlow::PoisseuilleFlow(const double r, const double Re) - : radius(r) + : inv_sqr_radius(1 / r / r) , Reynolds(Re) { Assert(Reynolds != 0., ExcMessage("Reynolds number cannot be zero")); @@ -203,8 +203,7 @@ namespace Functions const std::vector> & points, std::vector> &values) const { - unsigned int n = points.size(); - double stretch = 1. / radius; + const unsigned int n = points.size(); Assert(values.size() == dim + 1, ExcDimensionMismatch(values.size(), dim + 1)); @@ -214,13 +213,12 @@ namespace Functions for (unsigned int k = 0; k < n; ++k) { const Point &p = points[k]; - // First, compute the - // square of the distance to - // the x-axis divided by the - // radius. + // First, compute the square of the distance to the x-axis divided by + // the radius. double r2 = 0; for (unsigned int d = 1; d < dim; ++d) - r2 += p(d) * p(d) * stretch * stretch; + r2 += p(d) * p(d); + r2 *= inv_sqr_radius; // x-velocity values[0][k] = 1. - r2; @@ -228,7 +226,7 @@ namespace Functions for (unsigned int d = 1; d < dim; ++d) values[d][k] = 0.; // pressure - values[dim][k] = -2 * (dim - 1) * stretch * stretch * p(0) / Reynolds + + values[dim][k] = -2 * (dim - 1) * inv_sqr_radius * p(0) / Reynolds + this->mean_pressure; } } @@ -241,8 +239,7 @@ namespace Functions const std::vector> & points, std::vector>> &values) const { - unsigned int n = points.size(); - double stretch = 1. / radius; + const unsigned int n = points.size(); Assert(values.size() == dim + 1, ExcDimensionMismatch(values.size(), dim + 1)); @@ -255,12 +252,12 @@ namespace Functions // x-velocity values[0][k][0] = 0.; for (unsigned int d = 1; d < dim; ++d) - values[0][k][d] = -2. * p(d) * stretch * stretch; + values[0][k][d] = -2. * p(d) * inv_sqr_radius; // other velocities for (unsigned int d = 1; d < dim; ++d) values[d][k] = 0.; // pressure - values[dim][k][0] = -2 * (dim - 1) * stretch * stretch / Reynolds; + values[dim][k][0] = -2 * (dim - 1) * inv_sqr_radius / Reynolds; for (unsigned int d = 1; d < dim; ++d) values[dim][k][d] = 0.; } @@ -274,7 +271,7 @@ namespace Functions const std::vector> & points, std::vector> &values) const { - unsigned int n = points.size(); + const unsigned int n = points.size(); (void)n; Assert(values.size() == dim + 1, ExcDimensionMismatch(values.size(), dim + 1));