29 #ifndef NEMOSYS_CONSERVATIVEVOLUMETRANSFER_H_ 30 #define NEMOSYS_CONSERVATIVEVOLUMETRANSFER_H_ 34 #include "nemosys_export.h" 36 #include <vtkQuadratureSchemeDefinition.h> 38 #include <Eigen/Dense> 39 #include <Eigen/Sparse> 40 #include <Eigen/SparseLU> 50 std::cout <<
"Conservative Volume Transfer destroyed" << std::endl;
59 static std::shared_ptr<ConservativeVolumeTransfer>
61 return std::shared_ptr<ConservativeVolumeTransfer>(
67 const std::vector<int> &arrayIDs,
68 const std::vector<std::string> &newnames = std::vector<std::string>());
71 std::vector<vtkSmartPointer<vtkDoubleArray>> &dasSource,
72 std::vector<vtkSmartPointer<vtkDoubleArray>> &dasTarget,
79 const std::vector<int> &arrayIDs,
80 const std::vector<std::string> &newnames = std::vector<std::string>()) {
85 int i, vtkSmartPointer<vtkGenericCell> genCell,
86 std::vector<vtkSmartPointer<vtkDoubleArray>> &dasSourceToPoint,
87 std::vector<vtkSmartPointer<vtkDoubleArray>> &dasTarget) {
91 int run(
const std::vector<std::string> &newnames =
92 std::vector<std::string>())
override;
96 int transferSurfacePointData(
const std::vector<int> &arrayIDs);
113 int constructSupermesh();
117 int constructMassMatrix();
122 int constructMixedMassMatrix();
126 int constructLoadVector(
int arrayId,
int componentId);
130 int transfer(
int arrayId);
136 int convertSupermeshToUnstructuredGrid();
152 return supermeshGrid;
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};
169 double TET8W[8] = {0.1362178425370874, 0.1137821574629126, 0.1362178425370874,
170 0.1137821574629126, 0.1362178425370874, 0.1137821574629126,
171 0.1362178425370874, 0.1137821574629126};
173 double TET4[16] = {.5854101966249684, .1381966011250105, .1381966011250105,
174 .1381966011250105, .1381966011250105, .5854101966249684,
175 .1381966011250105, .1381966011250105, .1381966011250105,
176 .1381966011250105, .5854101966249684, .1381966011250105,
177 .1381966011250105, .1381966011250105, .1381966011250105,
180 double TET4W[4] = {0.250000000000000, 0.250000000000000, 0.250000000000000,
184 .7784952948213300, .0738349017262234, .0738349017262234,
185 .0738349017262234, .0738349017262234, .7784952948213300,
186 .0738349017262234, .0738349017262234, .0738349017262234,
187 .0738349017262234, .7784952948213300, .0738349017262234,
188 .0738349017262234, .0738349017262234, .0738349017262234,
191 .4062443438840510, .4062443438840510, .0937556561159491,
192 .0937556561159491, .4062443438840510, .0937556561159491,
193 .4062443438840510, .0937556561159491, .4062443438840510,
194 .0937556561159491, .0937556561159491, .4062443438840510,
196 .0937556561159491, .4062443438840510, .4062443438840510,
197 .0937556561159491, .0937556561159491, .4062443438840510,
198 .0937556561159491, .4062443438840510, .0937556561159491,
199 .0937556561159491, .4062443438840510, .4062443438840510,
202 double TET10W[10] = {
203 .0476331348432089, .0476331348432089, .0476331348432089,
206 .1349112434378610, .1349112434378610, .1349112434378610,
207 .1349112434378610, .1349112434378610, .1349112434378610,
223 Eigen::SparseMatrix<double, Eigen::ColMajor, long>
massMatrix;
234 void getLibSupermeshData(vtkDataSet *
data,
long &nnodes,
int &dim,
235 long &nelements,
int &loc,
double *&positions,
236 long *&enlist,
long &initTetId);
240 int interpolateCellDataToPoints(
int cellArrayId);
242 vtkSmartPointer<vtkPolyData> extractSurface(vtkUnstructuredGrid *grid);
246 bool surfaceTransferEnabled =
false;
249 #endif // NEMOSYS_CONSERVATIVEVOLUMETRANSFER_H_
std::vector< std::vector< double > > children_a
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.
Eigen::VectorXd loadVector
static std::shared_ptr< ConservativeVolumeTransfer > CreateShared(meshBase *_source, meshBase *_target)
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.
std::vector< long > parents_b
Eigen::SparseMatrix< double, Eigen::ColMajor, long > mixedMassMatrix
std::vector< long > parents_a
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.
~ConservativeVolumeTransfer() override