NEMoSys  0.63.0
A modular, extensible resource with robust automated mesh generation, mesh quality analysis, adaptive mesh refinement, and data transfer between arbitrary meshes.
ConservativeVolumeTransfer.H
Go to the documentation of this file.
1 /*******************************************************************************
2 * Promesh *
3 * Copyright (C) 2022, IllinoisRocstar LLC. All rights reserved. *
4 * *
5 * Promesh is the property of IllinoisRocstar LLC. *
6 * *
7 * IllinoisRocstar LLC *
8 * Champaign, IL *
9 * www.illinoisrocstar.com *
10 * promesh@illinoisrocstar.com *
11 *******************************************************************************/
12 /*******************************************************************************
13 * This file is part of Promesh *
14 * *
15 * This version of Promesh is free software: you can redistribute it and/or *
16 * modify it under the terms of the GNU Lesser General Public License as *
17 * published by the Free Software Foundation, either version 3 of the License, *
18 * or (at your option) any later version. *
19 * *
20 * Promesh is distributed in the hope that it will be useful, but WITHOUT ANY *
21 * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS *
22 * FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more *
23 * details. *
24 * *
25 * You should have received a copy of the GNU Lesser General Public License *
26 * along with this program. If not, see <https://www.gnu.org/licenses/>. *
27 * *
28 *******************************************************************************/
29 #ifndef NEMOSYS_CONSERVATIVEVOLUMETRANSFER_H_
30 #define NEMOSYS_CONSERVATIVEVOLUMETRANSFER_H_
31 
33 #include "Transfer/TransferBase.H"
34 #include "nemosys_export.h"
35 
36 #include <vtkQuadratureSchemeDefinition.h>
37 
38 #include <Eigen/Dense>
39 #include <Eigen/Sparse>
40 #include <Eigen/SparseLU>
41 
42 #include <cstring>
43 
44 class NEMOSYS_EXPORT ConservativeVolumeTransfer : public TransferBase {
45 
46 public:
47  ConservativeVolumeTransfer(meshBase *_source, meshBase *_target);
48 
50  std::cout << "Conservative Volume Transfer destroyed" << std::endl;
51  }
52 
53 public:
55  meshBase *_target) {
56  return new ConservativeVolumeTransfer(_source, _target);
57  }
58 
59  static std::shared_ptr<ConservativeVolumeTransfer>
60  CreateShared(meshBase *_source, meshBase *_target) {
61  return std::shared_ptr<ConservativeVolumeTransfer>(
62  ConservativeVolumeTransfer::Create(_source, _target));
63  }
64 
65 public:
67  const std::vector<int> &arrayIDs,
68  const std::vector<std::string> &newnames = std::vector<std::string>());
69 
70  int transferPointData(int i, vtkSmartPointer<vtkGenericCell> genCell,
71  std::vector<vtkSmartPointer<vtkDoubleArray>> &dasSource,
72  std::vector<vtkSmartPointer<vtkDoubleArray>> &dasTarget,
73  bool flip) {
74  return 0;
75  }
76 
77 public:
79  const std::vector<int> &arrayIDs,
80  const std::vector<std::string> &newnames = std::vector<std::string>()) {
81  return 0;
82  }
83 
85  int i, vtkSmartPointer<vtkGenericCell> genCell,
86  std::vector<vtkSmartPointer<vtkDoubleArray>> &dasSourceToPoint,
87  std::vector<vtkSmartPointer<vtkDoubleArray>> &dasTarget) {
88  return 0;
89  }
90 
91  int run(const std::vector<std::string> &newnames =
92  std::vector<std::string>()) override;
93 
94  int initialize();
95 
96  int transferSurfacePointData(const std::vector<int> &arrayIDs);
97 
98 public:
99  /**
100  * @brief Enable surface transfer. This is necessary for meshes with a shared
101  * boundary.
102  */
103  void enableSurfaceTransfer() { surfaceTransferEnabled = true; }
104 
105  /**
106  * @brief Disable surface transfer.
107  */
108  void disableSurfaceTransfer() { surfaceTransferEnabled = false; }
109 
110  /**
111  * @brief Construct supermesh.
112  */
113  int constructSupermesh();
114  /**
115  * @brief Construct mass matrix corresponding to target.
116  */
117  int constructMassMatrix();
118  /**
119  * @brief Construct mixed mass matrix corresponding to the source and target.
120  * Requires a constructed supermesh.
121  */
122  int constructMixedMassMatrix();
123  /**
124  * @brief Construct load vector corresponding to a given array and component.
125  */
126  int constructLoadVector(int arrayId, int componentId);
127  /**
128  * @brief Transfer array given its id.
129  */
130  int transfer(int arrayId);
131 
132  /**
133  * @brief Convertes supermesh to an unstructured grid.
134  * @see getSupermeshGrid
135  */
136  int convertSupermeshToUnstructuredGrid();
137 
138  /**
139  * @brief Get the source grid in the transfer.
140  */
141  vtkSmartPointer<vtkUnstructuredGrid> getSourceGrid() { return sourceGrid; }
142  /**
143  * @brief Get the target grid in the transfer.
144  */
145  vtkSmartPointer<vtkUnstructuredGrid> getTargetGrid() { return targetGrid; }
146 
147  /**
148  * @brief Get the supermesh grid in the transfer. Call
149  * convertSupermeshToUnstructuredGrid before this query.
150  */
151  vtkSmartPointer<vtkUnstructuredGrid> getSupermeshGrid() {
152  return supermeshGrid;
153  }
154 
155 private:
156  double TET8[32] = {
157  0.3281633025163817, 0.3281633025163817, 0.3281633025163817,
158  0.01551009245085488, 0.1080472498984286, 0.1080472498984286,
159  0.1080472498984286, 0.6758582503047142, 0.3281633025163817,
160  0.3281633025163817, 0.01551009245085488, 0.3281633025163817,
161  0.1080472498984286, 0.1080472498984286, 0.6758582503047142,
162  0.1080472498984286, 0.3281633025163817, 0.01551009245085488,
163  0.3281633025163817, 0.3281633025163817, 0.1080472498984286,
164  0.6758582503047142, 0.1080472498984286, 0.1080472498984286,
165  0.01551009245085488, 0.3281633025163817, 0.3281633025163817,
166  0.3281633025163817, 0.6758582503047142, 0.1080472498984286,
167  0.1080472498984286, 0.1080472498984286};
168 
169  double TET8W[8] = {0.1362178425370874, 0.1137821574629126, 0.1362178425370874,
170  0.1137821574629126, 0.1362178425370874, 0.1137821574629126,
171  0.1362178425370874, 0.1137821574629126};
172 
173  double TET4[16] = {.5854101966249684, .1381966011250105, .1381966011250105,
174  .1381966011250105, .1381966011250105, .5854101966249684,
175  .1381966011250105, .1381966011250105, .1381966011250105,
176  .1381966011250105, .5854101966249684, .1381966011250105,
177  .1381966011250105, .1381966011250105, .1381966011250105,
178  .5854101966249684};
179 
180  double TET4W[4] = {0.250000000000000, 0.250000000000000, 0.250000000000000,
181  0.250000000000000};
182 
183  double TET10[40] = {
184  .7784952948213300, .0738349017262234, .0738349017262234,
185  .0738349017262234, .0738349017262234, .7784952948213300,
186  .0738349017262234, .0738349017262234, .0738349017262234,
187  .0738349017262234, .7784952948213300, .0738349017262234,
188  .0738349017262234, .0738349017262234, .0738349017262234,
189  .7784952948213300,
190 
191  .4062443438840510, .4062443438840510, .0937556561159491,
192  .0937556561159491, .4062443438840510, .0937556561159491,
193  .4062443438840510, .0937556561159491, .4062443438840510,
194  .0937556561159491, .0937556561159491, .4062443438840510,
195 
196  .0937556561159491, .4062443438840510, .4062443438840510,
197  .0937556561159491, .0937556561159491, .4062443438840510,
198  .0937556561159491, .4062443438840510, .0937556561159491,
199  .0937556561159491, .4062443438840510, .4062443438840510,
200  };
201 
202  double TET10W[10] = {
203  .0476331348432089, .0476331348432089, .0476331348432089,
204  .0476331348432089,
205 
206  .1349112434378610, .1349112434378610, .1349112434378610,
207  .1349112434378610, .1349112434378610, .1349112434378610,
208  };
209 
212 
215 
216  vtkSmartPointer<vtkUnstructuredGrid> sourceGrid;
217  vtkSmartPointer<vtkUnstructuredGrid> targetGrid;
218  vtkSmartPointer<vtkUnstructuredGrid> supermeshGrid;
219 
220  vtkSmartPointer<vtkQuadratureSchemeDefinition> quadrature;
221 
222  // Eigen requires indexing via a signed type (can't use size_t)
223  Eigen::SparseMatrix<double, Eigen::ColMajor, long> massMatrix;
224  Eigen::SparseMatrix<double, Eigen::ColMajor, long> mixedMassMatrix;
225  Eigen::VectorXd loadVector;
226 
227  std::vector<std::vector<double>> children_a;
228  std::vector<std::vector<double>> children_b;
229 
230  std::vector<std::vector<double>> tets_c;
231  std::vector<long> parents_a;
232  std::vector<long> parents_b;
233 
234  void getLibSupermeshData(vtkDataSet *data, long &nnodes, int &dim,
235  long &nelements, int &loc, double *&positions,
236  long *&enlist, long &initTetId);
237 
238  // Given cell array id, interpolate to points in source mesh and return
239  // corresponding point array id in source.
240  int interpolateCellDataToPoints(int cellArrayId);
241 
242  vtkSmartPointer<vtkPolyData> extractSurface(vtkUnstructuredGrid *grid);
243 
244  void volumeCheck();
245 
246  bool surfaceTransferEnabled = false;
247 };
248 
249 #endif // NEMOSYS_CONSERVATIVEVOLUMETRANSFER_H_
std::vector< std::vector< double > > children_a
double TET4W[]
Definition: Cubature.C:77
data_type data
Edge/face with sorted point ids (a, b, c, ...) is located at some index i in data[b], with data[b][i].first == [a, c] (for edges, third point id treated as -1).
vtkSmartPointer< vtkUnstructuredGrid > targetGrid
std::vector< std::vector< double > > children_b
static ConservativeVolumeTransfer * Create(meshBase *_source, meshBase *_target)
std::vector< std::vector< double > > tets_c
A brief description of meshBase.
Definition: meshBase.H:64
static std::shared_ptr< ConservativeVolumeTransfer > CreateShared(meshBase *_source, meshBase *_target)
double TET4[]
Definition: Cubature.C:68
int transferPointData(int i, vtkSmartPointer< vtkGenericCell > genCell, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasSource, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasTarget, bool flip)
void disableSurfaceTransfer()
Disable surface transfer.
virtual int run(const std::vector< std::string > &newnames=std::vector< std::string >())=0
Transfer all fields.
Eigen::SparseMatrix< double, Eigen::ColMajor, long > massMatrix
vtkSmartPointer< vtkUnstructuredGrid > getTargetGrid()
Get the target grid in the transfer.
int transferCellData(int i, vtkSmartPointer< vtkGenericCell > genCell, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasSourceToPoint, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasTarget)
vtkSmartPointer< vtkUnstructuredGrid > supermeshGrid
virtual int transferPointData(const std::vector< int > &arrayIDs, const std::vector< std::string > &newnames=std::vector< std::string >())=0
Transfer point data with given ids from source to target.
Eigen::SparseMatrix< double, Eigen::ColMajor, long > mixedMassMatrix
vtkSmartPointer< vtkUnstructuredGrid > sourceGrid
int transferCellData(const std::vector< int > &arrayIDs, const std::vector< std::string > &newnames=std::vector< std::string >())
Transfer cell data with given ids from source to target.
vtkSmartPointer< vtkUnstructuredGrid > getSupermeshGrid()
Get the supermesh grid in the transfer.
void enableSurfaceTransfer()
Enable surface transfer.
vtkSmartPointer< vtkQuadratureSchemeDefinition > quadrature
vtkSmartPointer< vtkUnstructuredGrid > getSourceGrid()
Get the source grid in the transfer.