GCC Code Coverage Report

Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/assembly/jacobianpattern.hh
Date: 2024-05-04 19:09:25
Exec Total Coverage
Lines: 85 85 100.0%
Functions: 93 97 95.9%
Branches: 118 170 69.4%
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 Assembly
10 * \brief Helper function to generate Jacobian pattern for different discretization methods
11 */
12 #ifndef DUMUX_JACOBIAN_PATTERN_HH
13 #define DUMUX_JACOBIAN_PATTERN_HH
14
15 #include <type_traits>
16 #include <dune/istl/matrixindexset.hh>
17 #include <dumux/discretization/method.hh>
18
19 namespace Dumux {
20
21 /*!
22 * \ingroup Assembly
23 * \brief Helper function to generate Jacobian pattern for cell-centered methods
24 */
25 template<bool isImplicit, class GridGeometry,
26 typename std::enable_if_t<( (GridGeometry::discMethod == DiscretizationMethods::cctpfa)
27 || (GridGeometry::discMethod == DiscretizationMethods::ccmpfa) ), int> = 0>
28 441 Dune::MatrixIndexSet getJacobianPattern(const GridGeometry& gridGeometry)
29 {
30
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 441 const auto numDofs = gridGeometry.numDofs();
31
1/2
✓ Branch 1 taken 279 times.
✗ Branch 2 not taken.
 441 Dune::MatrixIndexSet pattern;
32
1/2
✓ Branch 1 taken 279 times.
✗ Branch 2 not taken.
 441 pattern.resize(numDofs, numDofs);
33
34 // matrix pattern for implicit Jacobians
35 if (isImplicit)
36 {
37
2/2
✓ Branch 0 taken 1047913 times.
✓ Branch 1 taken 277 times.
 2376738 for (unsigned int globalI = 0; globalI < numDofs; ++globalI)
38 {
39
2/4
✓ Branch 1 taken 1047913 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1047913 times.
✗ Branch 5 not taken.
 4752614 pattern.add(globalI, globalI);
40
4/4
✓ Branch 0 taken 4253808 times.
✓ Branch 1 taken 1047913 times.
✓ Branch 2 taken 4253808 times.
✓ Branch 3 taken 1047913 times.
 34643157 for (const auto& dataJ : gridGeometry.connectivityMap()[globalI])
41
1/2
✓ Branch 1 taken 4253808 times.
✗ Branch 2 not taken.
 11413266 pattern.add(dataJ.globalJ, globalI);
42 }
43 }
44
45 // matrix pattern for explicit Jacobians -> diagonal matrix
46 else
47 {
48
2/2
✓ Branch 0 taken 5000 times.
✓ Branch 1 taken 2 times.
 18102 for (unsigned int globalI = 0; globalI < numDofs; ++globalI)
49
1/2
✓ Branch 1 taken 5000 times.
✗ Branch 2 not taken.
 18092 pattern.add(globalI, globalI);
50 }
51
52 441 return pattern;
53 }
54
55 /*!
56 * \ingroup Assembly
57 * \brief Helper function to generate Jacobian pattern for the staggered method
58 */
59 template<bool isImplicit, class GridGeometry,
60 typename std::enable_if_t<( (GridGeometry::discMethod == DiscretizationMethods::staggered) ), int> = 0>
61 102 auto getJacobianPattern(const GridGeometry& gridGeometry)
62 {
63 // resize the jacobian and the residual
64 102 const auto numDofs = gridGeometry.numDofs();
65 102 Dune::MatrixIndexSet pattern(numDofs, numDofs);
66
67 102 const auto& connectivityMap = gridGeometry.connectivityMap();
68
69 102 auto fvGeometry = localView(gridGeometry);
70 // evaluate the actual pattern
71 551756 for (const auto& element : elements(gridGeometry.gridView()))
72 {
73 if(gridGeometry.isCellCenter())
74 {
75 // the global index of the element at hand
76 static constexpr auto cellCenterIdx = GridGeometry::cellCenterIdx();
77 179184 const auto ccGlobalI = gridGeometry.elementMapper().index(element);
78 89592 pattern.add(ccGlobalI, ccGlobalI);
79 1063452 for (auto&& ccGlobalJ : connectivityMap(cellCenterIdx, cellCenterIdx, ccGlobalI))
80 352542 pattern.add(ccGlobalI, ccGlobalJ);
81 }
82 else
83 {
84 static constexpr auto faceIdx = GridGeometry::faceIdx();
85 89592 fvGeometry.bindElement(element);
86
87 // loop over sub control faces
88 537552 for (auto&& scvf : scvfs(fvGeometry))
89 {
90 358368 const auto faceGlobalI = scvf.dofIndex();
91 358368 pattern.add(faceGlobalI, faceGlobalI);
92 6543216 for (auto&& faceGlobalJ : connectivityMap(faceIdx, faceIdx, scvf.index()))
93 2554872 pattern.add(faceGlobalI, faceGlobalJ);
94 }
95 }
96 }
97
98 102 return pattern;
99 }
100
101 /*!
102 * \ingroup Assembly
103 * \brief Helper function to generate Jacobian pattern for finite element scheme
104 */
105 template<class FEBasis>
106 17 Dune::MatrixIndexSet getFEJacobianPattern(const FEBasis& feBasis)
107 {
108 17 const auto numDofs = feBasis.size();
109
110
1/2
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
 17 Dune::MatrixIndexSet pattern;
111
1/2
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
 17 pattern.resize(numDofs, numDofs);
112
113
1/2
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
 18 auto localView = feBasis.localView();
114 // matrix pattern for implicit Jacobians
115
4/8
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 13 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 76913 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 76900 times.
✗ Branch 11 not taken.
 154463 for (const auto& element : elements(feBasis.gridView()))
116 {
117
1/2
✓ Branch 1 taken 76900 times.
✗ Branch 2 not taken.
 77213 localView.bind(element);
118
119 77213 const auto& finiteElement = localView.tree().finiteElement();
120
2/4
✓ Branch 1 taken 76900 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 76900 times.
✗ Branch 5 not taken.
 77213 const auto numLocalDofs = finiteElement.localBasis().size();
121
2/2
✓ Branch 0 taken 192900 times.
✓ Branch 1 taken 76900 times.
 272352 for (std::size_t i = 0; i < numLocalDofs; i++)
122 {
123 195139 const auto dofIdxI = localView.index(i);
124
2/2
✓ Branch 0 taken 660900 times.
✓ Branch 1 taken 192900 times.
 873956 for (std::size_t j = 0; j < numLocalDofs; j++)
125 {
126
1/2
✓ Branch 1 taken 660900 times.
✗ Branch 2 not taken.
 678817 const auto dofIdxJ = localView.index(j);
127
3/6
✓ Branch 1 taken 660900 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 660900 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 660900 times.
✗ Branch 8 not taken.
 2036451 pattern.add(dofIdxI, dofIdxJ);
128 }
129 }
130 }
131
132 17 return pattern;
133 }
134
135 /*!
136 * \ingroup Assembly
137 * \brief Helper function to generate Jacobian pattern for finite element scheme
138 * \note This interface is for compatibility with the other schemes. The pattern
139 * in fem is the same independent of the time discretization scheme.
140 */
141 template<bool isImplicit, class GridGeometry,
142 typename std::enable_if_t<(GridGeometry::discMethod == DiscretizationMethods::fem), int> = 0>
143 Dune::MatrixIndexSet getJacobianPattern(const GridGeometry& gridGeometry)
144 { return getFEJacobianPattern(gridGeometry.feBasis()); }
145
146 /*!
147 * \ingroup Assembly
148 * \brief Helper function to generate Jacobian pattern for the face-centered methods
149 */
150 template<bool isImplicit, class GridGeometry,
151 typename std::enable_if_t<(GridGeometry::discMethod == DiscretizationMethods::fcstaggered), int> = 0>
152 60 Dune::MatrixIndexSet getJacobianPattern(const GridGeometry& gridGeometry)
153 {
154 // resize the jacobian and the residual
155 60 const auto numDofs = gridGeometry.numDofs();
156 60 Dune::MatrixIndexSet pattern(numDofs, numDofs);
157
158
1/2
✓ Branch 1 taken 55 times.
✗ Branch 2 not taken.
 60 const auto& connectivityMap = gridGeometry.connectivityMap();
159
1/2
✓ Branch 1 taken 55 times.
✗ Branch 2 not taken.
 60 auto fvGeometry = localView(gridGeometry);
160
161 // set the pattern
162
14/19
✓ Branch 1 taken 59 times.
✓ Branch 2 taken 101 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 159 times.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 330591 times.
✓ Branch 7 taken 6 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 1250 times.
✓ Branch 10 taken 998 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 3536 times.
✓ Branch 13 taken 6 times.
✓ Branch 14 taken 3536 times.
✓ Branch 15 taken 6 times.
✓ Branch 17 taken 3536 times.
✗ Branch 18 not taken.
✓ Branch 20 taken 3536 times.
✗ Branch 21 not taken.
 664935 for (const auto& element : elements(gridGeometry.gridView()))
163 {
164
1/2
✓ Branch 1 taken 5778 times.
✗ Branch 2 not taken.
 334174 fvGeometry.bind(element);
165
5/5
✓ Branch 0 taken 12072 times.
✓ Branch 1 taken 1329442 times.
✓ Branch 2 taken 331156 times.
✓ Branch 3 taken 1326424 times.
✓ Branch 4 taken 331156 times.
 1672670 for (const auto& scv : scvs(fvGeometry))
166 {
167 1338496 const auto globalI = scv.dofIndex();
168
1/2
✓ Branch 1 taken 1338496 times.
✗ Branch 2 not taken.
 1338496 pattern.add(globalI, globalI);
169
170
4/4
✓ Branch 0 taken 9335360 times.
✓ Branch 1 taken 1338496 times.
✓ Branch 2 taken 9335360 times.
✓ Branch 3 taken 1338496 times.
 14689344 for (const auto& scvIdxJ : connectivityMap[scv.index()])
171 {
172
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 9252480 times.
 9335360 const auto globalJ = fvGeometry.scv(scvIdxJ).dofIndex();
173
1/2
✓ Branch 1 taken 9335360 times.
✗ Branch 2 not taken.
 9335360 pattern.add(globalI, globalJ);
174
175
2/2
✓ Branch 0 taken 2720 times.
✓ Branch 1 taken 9332640 times.
 9335360 if (gridGeometry.isPeriodic())
176 {
177
3/4
✓ Branch 1 taken 136 times.
✓ Branch 2 taken 2584 times.
✓ Branch 3 taken 136 times.
✗ Branch 4 not taken.
 2720 if (gridGeometry.dofOnPeriodicBoundary(globalI) && globalI != globalJ)
178 {
179
1/2
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
 136 const auto globalIP = gridGeometry.periodicallyMappedDof(globalI);
180
1/2
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
 136 pattern.add(globalIP, globalI);
181
1/2
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
 136 pattern.add(globalI, globalIP);
182
183
2/2
✓ Branch 0 taken 68 times.
✓ Branch 1 taken 68 times.
 136 if (globalI > globalIP)
184
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 68 pattern.add(globalIP, globalJ);
185 }
186 }
187 }
188 }
189 }
190
191 60 return pattern;
192 }
193
194 /*!
195 * \ingroup Assembly
196 * \brief Compute the Jacobian matrix sparsity pattern for control volume finite element schemes
197 */
198 template<bool isImplicit, class GridGeometry,
199 typename std::enable_if_t<DiscretizationMethods::isCVFE<typename GridGeometry::DiscretizationMethod>, int> = 0>
200 260 Dune::MatrixIndexSet getJacobianPattern(const GridGeometry& gridGeometry)
201 {
202 // resize the jacobian and the residual
203
1/3
✗ Branch 0 not taken.
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
 260 const auto numDofs = gridGeometry.numDofs();
204 260 Dune::MatrixIndexSet pattern(numDofs, numDofs);
205
206 // matrix pattern for implicit Jacobians
207 if constexpr (isImplicit)
208 {
209
1/4
✓ Branch 1 taken 163 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 305 auto fvGeometry = localView(gridGeometry);
210
18/21
✓ Branch 1 taken 165 times.
✓ Branch 2 taken 52439 times.
✓ Branch 3 taken 28 times.
✓ Branch 4 taken 52603 times.
✓ Branch 5 taken 29 times.
✓ Branch 6 taken 226975 times.
✓ Branch 7 taken 35504 times.
✓ Branch 8 taken 104293 times.
✓ Branch 9 taken 37 times.
✓ Branch 10 taken 51301 times.
✓ Branch 11 taken 104293 times.
✓ Branch 12 taken 50386 times.
✓ Branch 13 taken 27 times.
✓ Branch 14 taken 154679 times.
✓ Branch 15 taken 27 times.
✓ Branch 17 taken 50386 times.
✗ Branch 18 not taken.
✓ Branch 20 taken 50386 times.
✗ Branch 21 not taken.
✓ Branch 25 taken 19 times.
✗ Branch 26 not taken.
 846988 for (const auto& element : elements(gridGeometry.gridView()))
211 {
212
1/2
✓ Branch 1 taken 251393 times.
✗ Branch 2 not taken.
 585488 fvGeometry.bindElement(element);
213
4/4
✓ Branch 0 taken 1247041 times.
✓ Branch 1 taken 329698 times.
✓ Branch 2 taken 1173811 times.
✓ Branch 3 taken 293083 times.
 3389398 for (const auto& scv : scvs(fvGeometry))
214 {
215 2255037 const auto globalI = scv.dofIndex();
216
4/4
✓ Branch 0 taken 5278499 times.
✓ Branch 1 taken 1247041 times.
✓ Branch 2 taken 5132039 times.
✓ Branch 3 taken 1173811 times.
 14897997 for (const auto& scvJ : scvs(fvGeometry))
217 {
218 9470491 const auto globalJ = scvJ.dofIndex();
219
1/2
✓ Branch 1 taken 5278499 times.
✗ Branch 2 not taken.
 9470491 pattern.add(globalI, globalJ);
220
221
2/2
✓ Branch 0 taken 160000 times.
✓ Branch 1 taken 4835711 times.
 9470491 if (gridGeometry.isPeriodic())
222 {
223
4/4
✓ Branch 1 taken 1600 times.
✓ Branch 2 taken 158400 times.
✓ Branch 3 taken 1200 times.
✓ Branch 4 taken 400 times.
 160000 if (gridGeometry.dofOnPeriodicBoundary(globalI) && globalI != globalJ)
224 {
225
1/2
✓ Branch 1 taken 1200 times.
✗ Branch 2 not taken.
 1200 const auto globalIP = gridGeometry.periodicallyMappedDof(globalI);
226
1/2
✓ Branch 1 taken 1200 times.
✗ Branch 2 not taken.
 1200 pattern.add(globalIP, globalI);
227
1/2
✓ Branch 1 taken 1200 times.
✗ Branch 2 not taken.
 1200 pattern.add(globalI, globalIP);
228
229
2/2
✓ Branch 0 taken 600 times.
✓ Branch 1 taken 600 times.
 1200 if (globalI > globalIP)
230
1/2
✓ Branch 1 taken 600 times.
✗ Branch 2 not taken.
 600 pattern.add(globalIP, globalJ);
231 }
232 }
233 }
234 }
235 }
236 }
237
238 // matrix pattern for explicit Jacobians -> diagonal matrix
239 else
240 {
241
2/2
✓ Branch 1 taken 5202 times.
✓ Branch 2 taken 2 times.
 6408 for (unsigned int globalI = 0; globalI < numDofs; ++globalI)
242
1/2
✓ Branch 1 taken 5202 times.
✗ Branch 2 not taken.
 6404 pattern.add(globalI, globalI);
243 }
244
245 260 return pattern;
246 }
247
248
249 } // namespace Dumux
250
251 #endif
252