| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- | ||
| 2 | // vi: set et ts=4 sw=4 sts=4: | ||
| 3 | // | ||
| 4 | // SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder | ||
| 5 | // SPDX-License-Identifier: GPL-3.0-or-later | ||
| 6 | // | ||
| 7 | /*! | ||
| 8 | * \file | ||
| 9 | * \ingroup EmbeddedCoupling | ||
| 10 | * \brief Helper function to compute points on a circle | ||
| 11 | */ | ||
| 12 | |||
| 13 | #ifndef DUMUX_MULTIDOMAIN_EMBEDDED_CIRCLEPOINTS_HH | ||
| 14 | #define DUMUX_MULTIDOMAIN_EMBEDDED_CIRCLEPOINTS_HH | ||
| 15 | |||
| 16 | #include <vector> | ||
| 17 | #include <cmath> | ||
| 18 | |||
| 19 | #include <dune/common/exceptions.hh> | ||
| 20 | |||
| 21 | #include <dumux/common/math.hh> | ||
| 22 | #include <dumux/geometry/normal.hh> | ||
| 23 | |||
| 24 | namespace Dumux::EmbeddedCoupling { | ||
| 25 | |||
| 26 | /*! | ||
| 27 | * \ingroup EmbeddedCoupling | ||
| 28 | * \param points a vector of points to be filled | ||
| 29 | * \param sincos vector with [sin(a0), cos(a0), sin(a1), cos(a1), ...] for each circumferential sample point ai | ||
| 30 | * \param center the circle center | ||
| 31 | * \param normal the normal to the circle plane | ||
| 32 | * \param radius the circle radius | ||
| 33 | * \note This version allows for a more efficient circle point generator | ||
| 34 | * if the circumferential positions are fixed and can be reused. | ||
| 35 | * In this case, sine and cosine of the corresponding angles can be precomputed. | ||
| 36 | */ | ||
| 37 | template<class GlobalPosition, class Scalar> | ||
| 38 | 46279 | void circlePoints(std::vector<GlobalPosition>& points, | |
| 39 | const std::vector<Scalar>& sincos, | ||
| 40 | const GlobalPosition& center, | ||
| 41 | const GlobalPosition& normal, | ||
| 42 | const Scalar radius) | ||
| 43 | { | ||
| 44 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 46279 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 46279 times.
|
92558 | assert(sincos.size() % 2 == 0 && "Sample angles have to be pairs of sin/cos so size needs to be even."); |
| 45 | static_assert(GlobalPosition::dimension == 3, "Only implemented for world dimension 3"); | ||
| 46 | |||
| 47 | // resize the points vector | ||
| 48 | 46279 | const std::size_t numPoints = sincos.size()/2; | |
| 49 | 46279 | points.resize(numPoints); | |
| 50 | |||
| 51 | // make sure n is a unit vector | ||
| 52 | 46279 | auto n = normal; | |
| 53 | 92558 | n /= n.two_norm(); | |
| 54 | |||
| 55 | // calculate a vector u perpendicular to n | ||
| 56 | 46279 | auto u = unitNormal(n); u*= radius; | |
| 57 | |||
| 58 | // the circle parameterization is p(t) = r*cos(t)*u + r*sin(t)*(n x u) + c | ||
| 59 | 46279 | auto tangent = crossProduct(u, n); | |
| 60 | 92558 | tangent *= radius/tangent.two_norm(); | |
| 61 | |||
| 62 | // insert the vertices | ||
| 63 |
2/2✓ Branch 0 taken 3756998 times.
✓ Branch 1 taken 46279 times.
|
3803277 | for (std::size_t i = 0; i < numPoints; ++i) |
| 64 | { | ||
| 65 | 18784990 | points[i] = GlobalPosition({u[0]*sincos[2*i+1] + tangent[0]*sincos[2*i] + center[0], | |
| 66 | 11270994 | u[1]*sincos[2*i+1] + tangent[1]*sincos[2*i] + center[1], | |
| 67 | 15027992 | u[2]*sincos[2*i+1] + tangent[2]*sincos[2*i] + center[2]}); | |
| 68 | } | ||
| 69 | 46279 | } | |
| 70 | |||
| 71 | /*! | ||
| 72 | * \ingroup EmbeddedCoupling | ||
| 73 | * \brief returns a vector of sin and cos values of a circle parametrization | ||
| 74 | * \param numPoints the number of points | ||
| 75 | */ | ||
| 76 | template<class Scalar = double> | ||
| 77 | 22729 | std::vector<Scalar> circlePointsSinCos(const std::size_t numPoints) | |
| 78 | { | ||
| 79 | 45458 | std::vector<Scalar> sincos(2*numPoints); | |
| 80 | 22729 | Scalar t = 0 + 0.1; // add arbitrary offset | |
| 81 |
2/2✓ Branch 0 taken 319314 times.
✓ Branch 1 taken 22729 times.
|
342043 | for (std::size_t i = 0; i < numPoints; ++i) |
| 82 | { | ||
| 83 | using std::sin; using std::cos; | ||
| 84 |
2/2✓ Branch 0 taken 22729 times.
✓ Branch 1 taken 296585 times.
|
319314 | sincos[2*i] = sin(t); |
| 85 |
2/2✓ Branch 0 taken 22729 times.
✓ Branch 1 taken 296585 times.
|
319314 | sincos[2*i + 1] = cos(t); |
| 86 | 319314 | t += 2*M_PI/numPoints; | |
| 87 |
2/2✓ Branch 0 taken 22729 times.
✓ Branch 1 taken 296585 times.
|
319314 | if(t > 2*M_PI) t -= 2*M_PI; |
| 88 | } | ||
| 89 | 22729 | return sincos; | |
| 90 | } | ||
| 91 | |||
| 92 | /*! | ||
| 93 | * \ingroup EmbeddedCoupling | ||
| 94 | * \brief returns a vector of points on a circle | ||
| 95 | * \param center the circle center | ||
| 96 | * \param normal the normal to the circle plane | ||
| 97 | * \param radius the circle radius | ||
| 98 | * \param numPoints the number of points | ||
| 99 | */ | ||
| 100 | template<class GlobalPosition, class Scalar> | ||
| 101 | 22729 | std::vector<GlobalPosition> circlePoints(const GlobalPosition& center, | |
| 102 | const GlobalPosition& normal, | ||
| 103 | const Scalar radius, | ||
| 104 | const std::size_t numPoints = 20) | ||
| 105 | { | ||
| 106 | // precompute the sin/cos values | ||
| 107 |
1/4✓ Branch 1 taken 22729 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
|
45458 | const auto sincos = circlePointsSinCos<Scalar>(numPoints); |
| 108 | |||
| 109 |
1/2✓ Branch 1 taken 22729 times.
✗ Branch 2 not taken.
|
22729 | std::vector<GlobalPosition> points; |
| 110 |
1/2✓ Branch 1 taken 22729 times.
✗ Branch 2 not taken.
|
22729 | circlePoints(points, sincos, center, normal, radius); |
| 111 | 22729 | return points; | |
| 112 | } | ||
| 113 | |||
| 114 | } // end namespace Dumux::EmbeddedCoupling | ||
| 115 | |||
| 116 | #endif | ||
| 117 |