GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/test/multidomain/embedded/2d3d/1p_1p/problem_matrix.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 33 36 91.7%
Functions: 3 5 60.0%
Branches: 64 114 56.1%

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 EmbeddedTests
10 * \brief The matrix problem.
11 */
12
13 #ifndef DUMUX_MATRIX_ROBLEM_HH
14 #define DUMUX_MATRIX_ROBLEM_HH
15
16 #include <dune/geometry/quadraturerules.hh>
17 #include <dune/localfunctions/lagrange/pqkfactory.hh>
18
19 #include <dumux/common/properties.hh>
20 #include <dumux/common/boundarytypes.hh>
21 #include <dumux/common/math.hh>
22 #include <dumux/common/parameters.hh>
23 #include <dumux/common/numeqvector.hh>
24
25 #include <dumux/porousmediumflow/problem.hh>
26 #include <dumux/porousmediumflow/1p/incompressiblelocalresidual.hh>
27
28 namespace Dumux {
29
30 /*!
31 * \ingroup EmbeddedTests
32 * \brief The matrix problem.
33 */
34 template <class TypeTag>
35 class MatrixProblem : public PorousMediumFlowProblem<TypeTag>
36 {
37 using ParentType = PorousMediumFlowProblem<TypeTag>;
38 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
39 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
40 using GridView = typename GridGeometry::GridView;
41 using FVElementGeometry = typename GridGeometry::LocalView;
42 using SubControlVolume = typename GridGeometry::SubControlVolume;
43 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
44 using NumEqVector = Dumux::NumEqVector<PrimaryVariables>;
45 using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
46 using GridVariables = GetPropType<TypeTag, Properties::GridVariables>;
47 using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
48 using PointSource = GetPropType<TypeTag, Properties::PointSource>;
49 using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
50
51 enum {
52 // world dimension
53 dim = GridView::dimension,
54 dimWorld = GridView::dimensionworld
55 };
56
57 using Element = typename GridView::template Codim<0>::Entity;
58 using GlobalPosition = Dune::FieldVector<Scalar, dimWorld>;
59
60 using CouplingManager = GetPropType<TypeTag, Properties::CouplingManager>;
61
62 public:
63 2 MatrixProblem(std::shared_ptr<const GridGeometry> gridGeometry,
64 std::shared_ptr<typename ParentType::SpatialParams> spatialParams,
65 std::shared_ptr<CouplingManager> couplingManager,
66 const std::string& paramGroup = "Matrix")
67 : ParentType(gridGeometry, spatialParams, paramGroup)
68
3/12
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
8 , couplingManager_(couplingManager)
69 {
70 // read parameters from input file
71
9/26
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 2 times.
✗ Branch 14 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✓ Branch 19 taken 2 times.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✓ Branch 23 taken 2 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
4 name_ = getParam<std::string>("Vtk.OutputName") + "_" + getParamFromGroup<std::string>(this->paramGroup(), "Problem.Name");
72 2 }
73
74 /*!
75 * \name Problem parameters
76 */
77 // \{
78
79 /*!
80 * \brief The problem name.
81 *
82 * This is used as a prefix for files generated by the simulation.
83 */
84 const std::string& name() const
85
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 { return name_; }
86
87 // \}
88
89 /*!
90 * \name Boundary conditions
91 */
92 // \{
93
94 /*!
95 * \brief Specifies which kind of boundary condition should be
96 * used for which equation on a given boundary segment.
97 *
98 * \param globalPos The position for which the bc type should be evaluated
99 */
100 BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
101 {
102
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 9600 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 8000 times.
17600 BoundaryTypes values;
103
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 9600 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 8000 times.
17600 values.setAllNeumann();
104 return values;
105 }
106
107 // \}
108
109 /*!
110 * \name Volume terms
111 */
112 // \{
113
114 /*!
115 * \brief Applies a vector of point sources which are possibly solution dependent.
116 *
117 * \param pointSources A vector of Dumux::PointSource s that contain
118 source values for all phases and space positions.
119 *
120 * For this method, the \a values method of the point source
121 * has to return the absolute mass rate in kg/s. Positive values mean
122 * that mass is created, negative ones mean that it vanishes.
123 */
124 void addPointSources(std::vector<PointSource>& pointSources) const
125
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 { pointSources = this->couplingManager().bulkPointSources(); }
126
127 /*!
128 * \brief Evaluates the point sources (added by addPointSources)
129 * for all phases within a given sub-control volume.
130 *
131 * This is the method for the case where the point source is
132 * solution dependent and requires some quantities that
133 * are specific to the fully-implicit method.
134 *
135 * \param source A single point source
136 * \param element The finite element
137 * \param fvGeometry The finite-volume geometry
138 * \param elemVolVars All volume variables for the element
139 * \param scv The sub-control volume within the element
140 *
141 * For this method, the \a values() method of the point sources returns
142 * the absolute rate mass generated or annihilated in kg/s. Positive values mean
143 * that mass is created, negative ones mean that it vanishes.
144 */
145 template<class ElementVolumeVariables>
146 324928 void pointSource(PointSource& source,
147 const Element &element,
148 const FVElementGeometry& fvGeometry,
149 const ElementVolumeVariables& elemVolVars,
150 const SubControlVolume &scv) const
151 {
152 // compute source at every integration point
153 974784 const Scalar pressure3D = this->couplingManager().bulkPriVars(source.id())[Indices::pressureIdx];
154 974784 const Scalar pressure1D = this->couplingManager().lowDimPriVars(source.id())[Indices::pressureIdx];
155
156 // calculate the source
157 974784 const Scalar meanDistance = this->couplingManager().averageDistance(source.id());
158
5/6
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 324923 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 3 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
324928 static const Scalar matrixPerm = getParamFromGroup<Scalar>("Matrix", "SpatialParams.Permeability");
159
5/6
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 324924 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
324928 static const Scalar rho = getParam<Scalar>("Component.LiquidDensity");
160
5/6
✓ Branch 0 taken 7 times.
✓ Branch 1 taken 324921 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 5 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
324928 static const Scalar mu = getParam<Scalar>("Component.LiquidKinematicViscosity")*rho;
161 324928 const Scalar sourceValue = rho*(pressure1D - pressure3D)/meanDistance*matrixPerm/mu;
162 324928 source = sourceValue*source.quadratureWeight()*source.integrationElement();
163 324928 }
164
165 /*!
166 * \brief Evaluates the initial value for a control volume.
167 *
168 * \param globalPos The position for which the initial condition should be evaluated
169 *
170 * For this method, the \a values parameter stores primary
171 * variables.
172 */
173 PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
174 16000 { return PrimaryVariables({1e5}); }
175
176 //! Called after every time step
177 //! Output the total global exchange term
178 2 void computeSourceIntegral(const SolutionVector& sol, const GridVariables& gridVars)
179 {
180
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 NumEqVector source(0.0);
181
2/4
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
4 auto fvGeometry = localView(this->gridGeometry());
182
2/4
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
4 auto elemVolVars = localView(gridVars.curGridVolVars());
183
4/8
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 8002 times.
✗ Branch 10 not taken.
16006 for (const auto& element : elements(this->gridGeometry().gridView()))
184 {
185 8000 fvGeometry.bindElement(element);
186 8000 elemVolVars.bindElement(element, fvGeometry, sol);
187
188
5/6
✓ Branch 2 taken 8000 times.
✓ Branch 3 taken 8000 times.
✓ Branch 4 taken 8000 times.
✓ Branch 5 taken 8000 times.
✓ Branch 7 taken 8000 times.
✗ Branch 8 not taken.
32000 for (auto&& scv : scvs(fvGeometry))
189 {
190
1/2
✓ Branch 1 taken 8000 times.
✗ Branch 2 not taken.
8000 auto pointSources = this->scvPointSources(element, fvGeometry, elemVolVars, scv);
191 8000 pointSources *= scv.volume()*elemVolVars[scv].extrusionFactor();
192 8000 source += pointSources;
193 }
194 }
195
196
2/4
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
4 std::cout << "Global integrated source (3D): " << source << '\n';
197 2 }
198
199 //! Get the coupling manager
200 const CouplingManager& couplingManager() const
201
14/16
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 324923 times.
✓ Branch 2 taken 5 times.
✓ Branch 3 taken 324923 times.
✓ Branch 4 taken 5 times.
✓ Branch 5 taken 324923 times.
✓ Branch 6 taken 5 times.
✓ Branch 7 taken 324923 times.
✓ Branch 8 taken 5 times.
✓ Branch 9 taken 324923 times.
✓ Branch 10 taken 5 times.
✓ Branch 11 taken 324923 times.
✓ Branch 13 taken 2 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
1949572 { return *couplingManager_; }
202
203 private:
204
205 static constexpr Scalar eps_ = 1.5e-7;
206 std::string name_;
207
208 std::shared_ptr<CouplingManager> couplingManager_;
209 };
210
211 } // end namespace Dumux
212
213 #endif
214