GCC Code Coverage Report

Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/freeflow/navierstokes/momentum/velocityreconstruction.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 43 43 100.0%
Functions: 44 72 61.1%
Branches: 49 78 62.8%
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 NavierStokesModel
10 * \copydoc Dumux::StaggeredVelocityReconstruction
11 */
12 #ifndef DUMUX_NAVIERSTOKES_STAGGERED_VELOCITYRECONSTRUCTION_HH
13 #define DUMUX_NAVIERSTOKES_STAGGERED_VELOCITYRECONSTRUCTION_HH
14
15 #include <numeric>
16 #include <algorithm>
17 #include <dune/common/reservedvector.hh>
18 #include <dumux/discretization/method.hh>
19
20 namespace Dumux {
21
22 /*!
23 * \ingroup NavierStokesModel
24 * \brief Helper class for reconstructing the velocity.
25 */
26 struct StaggeredVelocityReconstruction
27 {
28 //! Return the velocity vector at the center of the primal grid.
29 template<class VelocityHelper, class FVElementGeometry>
30 922347 static auto cellCenterVelocity(const VelocityHelper& getFaceVelocity,
31 const FVElementGeometry& fvGeometry)
32 {
33 static_assert(FVElementGeometry::GridGeometry::discMethod == DiscretizationMethods::cctpfa);
34 using VelocityVector = typename FVElementGeometry::GridGeometry::GlobalCoordinate;
35 922347 VelocityVector result(0.0);
36
37 3693024 const auto directionIndex = [&](const auto& vector)
38 {
39
14/36
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 36 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 128 times.
✗ Branch 13 not taken.
✓ Branch 14 taken 128 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 3728 times.
✓ Branch 17 taken 7456 times.
✓ Branch 18 taken 3728 times.
✓ Branch 19 taken 7456 times.
✓ Branch 20 taken 1844548 times.
✓ Branch 21 taken 1844676 times.
✓ Branch 22 taken 1844548 times.
✓ Branch 23 taken 1844548 times.
✓ Branch 24 taken 1844548 times.
✗ Branch 25 not taken.
✓ Branch 26 taken 1844548 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 33 taken 3692824 times.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
 25861864 return std::find_if(vector.begin(), vector.end(), [](const auto& x) { return std::abs(x) > 1e-8; } ) - vector.begin();
40 };
41
42
4/4
✓ Branch 0 taken 3692988 times.
✓ Branch 1 taken 922347 times.
✓ Branch 2 taken 3692988 times.
✓ Branch 3 taken 922347 times.
 9230670 for (const auto& scvf : scvfs(fvGeometry))
43 {
44 7385976 const auto dirIdx = directionIndex(scvf.unitOuterNormal());
45 3693096 result[dirIdx] += 0.5*getFaceVelocity(fvGeometry, scvf)[dirIdx];
46 }
47 922347 return result;
48 }
49
50 //! Return the velocity vector at dof position given an scv
51 template<class SubControlVolume, class FVElementGeometry, class VelocityHelper>
52 90380 static auto faceVelocity(const SubControlVolume& scv,
53 const FVElementGeometry& fvGeometry,
54 const VelocityHelper& getNormalVelocityDofValue)
55 {
56 90380 int axis = scv.dofAxis();
57 90380 const int dim = FVElementGeometry::GridGeometry::GridView::dimension;
58 using Scalar = typename SubControlVolume::Traits::Scalar;
59 using GlobalPosition = typename FVElementGeometry::GridGeometry::GlobalCoordinate;
60 using VelocityVector = typename FVElementGeometry::GridGeometry::GlobalCoordinate;
61 90380 VelocityVector faceVelocityVector(0.0);
62
63 // per dimension, we have at max two velocities from which we'll compute an average
64
4/4
✓ Branch 0 taken 90380 times.
✓ Branch 1 taken 90380 times.
✓ Branch 2 taken 90380 times.
✓ Branch 3 taken 90380 times.
 271140 std::array<Dune::ReservedVector<Scalar, 2>, dim> normalVelocities;
65
4/6
✓ Branch 0 taken 576 times.
✓ Branch 1 taken 89804 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 576 times.
✓ Branch 5 taken 576 times.
✗ Branch 6 not taken.
 90380 if (scv.boundary() && !fvGeometry.gridGeometry().dofOnPeriodicBoundary(scv.dofIndex()))
66 {
67 // iterate through the inner lateral velocities,
68
4/4
✓ Branch 1 taken 2304 times.
✓ Branch 2 taken 576 times.
✓ Branch 3 taken 2304 times.
✓ Branch 4 taken 576 times.
 2880 for (const auto& scvf : scvfs(fvGeometry, scv))
69 {
70
2/2
✓ Branch 0 taken 1152 times.
✓ Branch 1 taken 1152 times.
 2304 if (scvf.isFrontal())
71 continue;
72
73 // at a lateral velocity, find the inner and outer normal velocities
74 1152 const auto& orthogonalScvf = fvGeometry.lateralOrthogonalScvf(scvf);
75
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 1152 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1152 times.
 2304 const auto& lateralScv = fvGeometry.scv(orthogonalScvf.insideScvIdx());
76 1152 auto lateralAxis = lateralScv.dofAxis();
77 2304 normalVelocities[lateralAxis].push_back( getNormalVelocityDofValue(lateralScv) ) ;
78 }
79 }
80 else
81 {
82 // Find the location of interpolation, if periodic, from both sides
83 89804 const GlobalPosition selfPosition = scv.dofPosition();
84 89804 GlobalPosition outsideSelfPosition = selfPosition;
85
3/4
✗ Branch 0 not taken.
✓ Branch 1 taken 89804 times.
✓ Branch 3 taken 12 times.
✓ Branch 4 taken 89792 times.
 89804 if (fvGeometry.gridGeometry().dofOnPeriodicBoundary(scv.dofIndex()))
86 12 outsideSelfPosition = fvGeometry.outsidePeriodicScv(scv).dofPosition();
87
88 // iterate through the lateral velocities to get to the normal locations
89
4/4
✓ Branch 1 taken 269412 times.
✓ Branch 2 taken 89804 times.
✓ Branch 3 taken 269412 times.
✓ Branch 4 taken 89804 times.
 359216 for (const auto& scvf : scvfs(fvGeometry, scv))
90 {
91
2/2
✓ Branch 0 taken 179608 times.
✓ Branch 1 taken 89804 times.
 269412 if (scvf.isFrontal())
92 continue;
93
94 // at a lateral velocity, find the inner and outer normal velocities
95 179608 const auto& orthogonalScvf = fvGeometry.lateralOrthogonalScvf(scvf);
96
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 179608 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 179608 times.
 359216 const auto& insideLateralScv = fvGeometry.scv(orthogonalScvf.insideScvIdx());
97 359216 const auto& outsideLateralScv = fvGeometry.scv(orthogonalScvf.outsideScvIdx());
98 179608 const auto& lateralAxis = insideLateralScv.dofAxis();
99
100 // Find the inside normal velocities
101 179608 const auto& insideNormalVelocity = getNormalVelocityDofValue(insideLateralScv);
102 359216 const auto& insideNormalPosition = insideLateralScv.dofPosition()[axis];
103
104 // Find the outside normal velocities
105 179608 const auto& outsideNormalVelocity = getNormalVelocityDofValue(outsideLateralScv);
106 359216 const auto& outsideNormalPosition = outsideLateralScv.dofPosition()[axis];
107
108 // Linear interpolation at the face plane and add to normal velocity collection
109 359216 const auto& innerDistance = std::abs(insideNormalPosition - selfPosition[axis]);
110 359216 const auto& totalDistance = std::abs(outsideNormalPosition - outsideSelfPosition[axis]) + innerDistance;
111 179608 const auto& velDiff = outsideNormalVelocity - insideNormalVelocity;
112
113 538824 normalVelocities[lateralAxis].push_back(insideNormalVelocity + (velDiff * innerDistance / totalDistance));
114 }
115 }
116
117
2/2
✓ Branch 0 taken 180760 times.
✓ Branch 1 taken 90380 times.
 271140 for (int i = 0; i < faceVelocityVector.size(); i++)
118 {
119
2/2
✓ Branch 0 taken 90380 times.
✓ Branch 1 taken 90380 times.
 180760 if (i == axis)
120 90452 faceVelocityVector[i] = getNormalVelocityDofValue(scv);
121 else
122 723040 faceVelocityVector[i] = std::accumulate(normalVelocities[i].begin(), normalVelocities[i].end(), 0.0) / normalVelocities[i].size();
123 }
124
125 90380 return faceVelocityVector;
126 }
127
128 };
129
130 } // end namespace Dumux
131
132 #endif // DUMUX_NAVIERSTOKES_STAGGERED_VELOCITYRECONSTRUCTION_HH
133