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.
ConservativeSurfaceTransfer Class Reference

Detailed Description

Definition at line 35 of file ConservativeSurfaceTransfer.H.

Public Member Functions

 ConservativeSurfaceTransfer (meshBase *_source, meshBase *_target)
 
 ~ConservativeSurfaceTransfer () override
 
int transferPointData (const std::vector< int > &arrayIDs=std::vector< int >(), const std::vector< std::string > &newnames=std::vector< std::string >())
 Transfer point data with given ids from source to target. More...
 
int transferPointData (int i, vtkSmartPointer< vtkGenericCell > genCell, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasSource, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasTarget, bool flip)
 
int writeOverlay ()
 
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. More...
 
int transferCellData (int i, vtkSmartPointer< vtkGenericCell > genCell, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasSourceToPoint, std::vector< vtkSmartPointer< vtkDoubleArray >> &dasTarget)
 
int run (const std::vector< std::string > &newnames=std::vector< std::string >())
 Transfer all fields. More...
 
void transfer (int arrayId)
 
int transferPointData (const std::vector< std::string > &arrayNames, const std::vector< std::string > &newnames=std::vector< std::string >())
 Transfer point data with given field names from source to target. More...
 
int transferCellData (const std::vector< std::string > &arrayNames, const std::vector< std::string > &newnames=std::vector< std::string >())
 Transfer cell data with given field names from source to target. More...
 
void setCheckQual (bool x)
 
void setContBool (bool x)
 

Static Public Member Functions

static ConservativeSurfaceTransferCreate (meshBase *_source, meshBase *_target)
 
static std::shared_ptr< ConservativeSurfaceTransferCreateShared (meshBase *_source, meshBase *_target)
 

Protected Attributes

meshBasesource
 
meshBasetarget
 
vtkSmartPointer< vtkStaticCellLocator > srcCellLocator = nullptr
 
vtkSmartPointer< vtkStaticCellLocator > trgCellLocator = nullptr
 
vtkSmartPointer< vtkStaticPointLocator > srcPointLocator = nullptr
 
vtkSmartPointer< vtkStaticPointLocator > trgPointLocator = nullptr
 
bool checkQual
 
bool continuous
 
double c2cTrnsDistTol
 

Private Member Functions

void extractDataFromVTK (vtkDataSet *data, std::vector< double > &coords, std::vector< int > &elems)
 

Private Attributes

int RFC_clear
 
int RFC_read
 
int RFC_write
 
int RFC_overlay
 
int RFC_transfer
 
int RFC_interp
 

Inherits TransferBase.

Constructor & Destructor Documentation

◆ ConservativeSurfaceTransfer()

ConservativeSurfaceTransfer::ConservativeSurfaceTransfer ( meshBase _source,
meshBase _target 
)

Definition at line 45 of file ConservativeSurfaceTransfer.C.

46  {
47  source = _source;
48  target = _target;
49 }
meshBase * target
Definition: TransferBase.H:106
meshBase * source
Definition: TransferBase.H:105

◆ ~ConservativeSurfaceTransfer()

ConservativeSurfaceTransfer::~ConservativeSurfaceTransfer ( )
inlineoverride

Definition at line 39 of file ConservativeSurfaceTransfer.H.

39  {
40  std::cout << "Conservative Surface Transfer destroyed" << std::endl;
41  }

Member Function Documentation

◆ Create()

static ConservativeSurfaceTransfer* ConservativeSurfaceTransfer::Create ( meshBase _source,
meshBase _target 
)
inlinestatic

Definition at line 44 of file ConservativeSurfaceTransfer.H.

Referenced by CreateShared().

45  {
46  return new ConservativeSurfaceTransfer(_source, _target);
47  }
ConservativeSurfaceTransfer(meshBase *_source, meshBase *_target)

◆ CreateShared()

static std::shared_ptr<ConservativeSurfaceTransfer> ConservativeSurfaceTransfer::CreateShared ( meshBase _source,
meshBase _target 
)
inlinestatic

Definition at line 50 of file ConservativeSurfaceTransfer.H.

References Create(), and TransferBase::transferPointData().

Referenced by NEM::DRV::TransferDriver::CreateTransferObject().

50  {
51  return std::shared_ptr<ConservativeSurfaceTransfer>(
52  ConservativeSurfaceTransfer::Create(_source, _target));
53  }
static ConservativeSurfaceTransfer * Create(meshBase *_source, meshBase *_target)

◆ extractDataFromVTK()

void ConservativeSurfaceTransfer::extractDataFromVTK ( vtkDataSet *  data,
std::vector< double > &  coords,
std::vector< int > &  elems 
)
private

Definition at line 242 of file ConservativeSurfaceTransfer.C.

Referenced by transferPointData(), and writeOverlay().

243  {
244  for (int ptId = 0; ptId < data->GetNumberOfPoints(); ++ptId) {
245  double x[3];
246 
247  data->GetPoint(ptId, x);
248 
249  coords.push_back(x[0]);
250  coords.push_back(x[1]);
251  coords.push_back(x[2]);
252  }
253 
254  for (int cellId = 0; cellId < data->GetNumberOfCells(); ++cellId) {
255  if (data->GetCellType(cellId) != VTK_TRIANGLE)
256  continue;
257 
258  vtkCell *cell = data->GetCell(cellId);
259 
260  elems.push_back(cell->GetPointId(0) + 1);
261  elems.push_back(cell->GetPointId(1) + 1);
262  elems.push_back(cell->GetPointId(2) + 1);
263  }
264 }
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).

◆ run()

int ConservativeSurfaceTransfer::run ( const std::vector< std::string > &  newnames = std::vector< std::string >())
inlinevirtual
Parameters
newnamesoptional array of names to be applied to transferred fields
Returns
0 upon completion

Implements TransferBase.

Definition at line 83 of file ConservativeSurfaceTransfer.H.

References data, and elems.

84  {
85  return 0;
86  }

◆ setCheckQual()

void TransferBase::setCheckQual ( bool  x)
inlineinherited

Definition at line 100 of file TransferBase.H.

100 { checkQual = x; }

◆ setContBool()

void TransferBase::setContBool ( bool  x)
inlineinherited

Definition at line 102 of file TransferBase.H.

102 { continuous = x; }

◆ transfer()

void ConservativeSurfaceTransfer::transfer ( int  arrayId)

◆ transferCellData() [1/3]

int ConservativeSurfaceTransfer::transferCellData ( const std::vector< int > &  arrayIDs,
const std::vector< std::string > &  newnames = std::vector< std::string >() 
)
inlinevirtual
Parameters
arrayIDsarray of array ids to specify which fields to transfer
newnamesoptional array of names to be applied to transferred fields
Returns
0 upon completion

Implements TransferBase.

Definition at line 70 of file ConservativeSurfaceTransfer.H.

72  {
73  return 0;
74  }

◆ transferCellData() [2/3]

int ConservativeSurfaceTransfer::transferCellData ( int  i,
vtkSmartPointer< vtkGenericCell >  genCell,
std::vector< vtkSmartPointer< vtkDoubleArray >> &  dasSourceToPoint,
std::vector< vtkSmartPointer< vtkDoubleArray >> &  dasTarget 
)
inline

Definition at line 76 of file ConservativeSurfaceTransfer.H.

79  {
80  return 0;
81  }

◆ transferCellData() [3/3]

int TransferBase::transferCellData ( const std::vector< std::string > &  arrayNames,
const std::vector< std::string > &  newnames = std::vector<std::string>() 
)
inherited
Parameters
arrayIDsarray of array names to specify which fields to transfer
newnamesoptional array of names to be applied to transferred fields
Returns
0 upon completion

Definition at line 47 of file TransferBase.C.

References TransferBase::getArrayIDs(), meshBase::getDataSet(), TransferBase::source, and TransferBase::transferCellData().

49 {
50  vtkCellData* cellData = source->getDataSet()->GetCellData();
51  std::vector<int> arrayIDs = getArrayIDs(arrayNames, cellData);
52 
53  return transferCellData(arrayIDs, newnames);
54 }
virtual int transferCellData(const std::vector< int > &arrayIDs, const std::vector< std::string > &newnames=std::vector< std::string >())=0
Transfer cell data with given ids from source to target.
vtkSmartPointer< vtkDataSet > getDataSet() const
get this meshes&#39; dataSet
Definition: meshBase.H:308
std::vector< int > getArrayIDs(const std::vector< std::string > &arrayNames, vtkFieldData *fieldData)
given array names and field data, return vector of corresponding array ids in the field data ...
Definition: TransferBase.C:56
meshBase * source
Definition: TransferBase.H:105

◆ transferPointData() [1/3]

int ConservativeSurfaceTransfer::transferPointData ( const std::vector< int > &  arrayIDs = std::vector< int >(),
const std::vector< std::string > &  newnames = std::vector< std::string >() 
)
virtual
Parameters
arrayIDsarray of array ids to specify which fields to transfer
newnamesoptional array of names to be applied to transferred fields
Returns
0 upon completion

Implements TransferBase.

Definition at line 51 of file ConservativeSurfaceTransfer.C.

References data, extractDataFromVTK(), meshBase::getDataSet(), NEM::MSH::New(), RFC_clear, RFC_interp, RFC_overlay, RFC_read, RFC_transfer, RFC_write, TransferBase::source, and TransferBase::target.

53  {
54 
55  COM_init(NULL, NULL);
56  COM_set_profiling(1);
57  // COM_set_verbose(10);
58 
59  MPI_Comm comm = MPI_COMM_WORLD;
60 
61  COM_LOAD_MODULE_STATIC_DYNAMIC(SurfX, "RFC");
62 
63  RFC_clear = COM_get_function_handle("RFC.clear_overlay");
64  RFC_read = COM_get_function_handle("RFC.read_overlay");
65  RFC_write = COM_get_function_handle("RFC.write_overlay");
66  RFC_overlay = COM_get_function_handle("RFC.overlay");
67  RFC_transfer = COM_get_function_handle("RFC.least_squares_transfer");
68  RFC_interp = COM_get_function_handle("RFC.interpolate");
69 
70  vtkDataSet *sourceDataSet = source->getDataSet();
71  vtkDataSet *targetDataSet = target->getDataSet();
72 
73  std::vector<double> sourceCoords;
74  std::vector<int> sourceElems;
75  extractDataFromVTK(sourceDataSet, sourceCoords, sourceElems);
76 
77  // === SOURCE WINDOW ===
78 
79  COM_new_window("sourceWindow");
80 
81  // set source nodal coordinates (nc)
82  COM_set_size("sourceWindow.nc", 1, sourceCoords.size() / 3);
83  COM_set_array("sourceWindow.nc", 1, &sourceCoords[0], 3);
84  // set source element node ids (:t3:)
85  COM_set_size("sourceWindow.:t3:", 1, sourceElems.size() / 3);
86  COM_set_array("sourceWindow.:t3:", 1, &sourceElems[0], 3);
87 
88  // initialize source soln array - this is reused
89  COM_new_dataitem("sourceWindow.soln", 'n', COM_DOUBLE, 1, "");
90 
91  COM_window_init_done("sourceWindow");
92 
93  // === TARGET WINDOW ===
94 
95  std::vector<double> targetCoords;
96  std::vector<int> targetElems;
97  extractDataFromVTK(target->getDataSet(), targetCoords, targetElems);
98 
99  COM_new_window("targetWindow");
100  // set target nodal coordinates (nc)
101  COM_set_size("targetWindow.nc", 1, targetCoords.size() / 3);
102  COM_set_array("targetWindow.nc", 1, &targetCoords[0], 3);
103  // set target element node ids (:t3:)
104  COM_set_size("targetWindow.:t3:", 1, targetElems.size() / 3);
105  COM_set_array("targetWindow.:t3:", 1, &targetElems[0], 3);
106 
107  // initialize target soln array - this is reused
108  std::vector<double> targetData(targetCoords.size(), -1.);
109  COM_new_dataitem("targetWindow.soln", 'n', COM_DOUBLE, 1, "m/s");
110  COM_resize_array("targetWindow.soln");
111 
112  COM_window_init_done("targetWindow");
113 
114  // === OVERLAY ===
115  int srcMesh = COM_get_dataitem_handle("sourceWindow.mesh");
116  int tgtMesh = COM_get_dataitem_handle("targetWindow.mesh");
117  COM_call_function(RFC_overlay, &srcMesh, &tgtMesh);
118 
119  // === TRANSFER ===
120  // set source data and transfer
121  int srcData = COM_get_dataitem_handle("sourceWindow.soln");
122  int tgtData = COM_get_dataitem_handle("targetWindow.soln");
123 
124  // source data buffer
125  std::vector<double> sourceData(sourceCoords.size());
126  // transfer specified arrays
127  for (const int &arrayId : arrayIds) {
128  vtkSmartPointer<vtkDataArray> sourceArray =
129  sourceDataSet->GetPointData()->GetArray(arrayId);
130 
131  auto targetPointValues = vtkSmartPointer<vtkDoubleArray>::New();
132  targetPointValues->SetName(sourceArray->GetName());
133  targetPointValues->SetNumberOfComponents(
134  sourceArray->GetNumberOfComponents());
135 
136  for (int componentId = 0;
137  componentId < sourceArray->GetNumberOfComponents(); ++componentId) {
138  auto sourcePointValues = vtkDoubleArray::SafeDownCast(
139  sourceDataSet->GetPointData()->GetArray(arrayId));
140 
141  for (int srcNodeId = 0; srcNodeId < sourceCoords.size() / 3;
142  ++srcNodeId) {
143  sourceData[srcNodeId] =
144  sourcePointValues->GetComponent(srcNodeId, componentId);
145  }
146 
147  COM_set_array("sourceWindow.soln", 1, &sourceData[0]);
148 
149  // transfer
150  COM_call_function(RFC_transfer, &srcData, &tgtData);
151 
152  void *addr;
153  double *data;
154 
155  // window var, pane number, addr
156  COM_get_array("targetWindow.soln", 1, &addr);
157  data = (double *)addr;
158 
159  for (int trgNodeId = 0; trgNodeId < targetCoords.size() / 3;
160  ++trgNodeId) {
161  targetPointValues->InsertComponent(trgNodeId, componentId,
162  data[trgNodeId]);
163  }
164  }
165 
166  targetDataSet->GetPointData()->AddArray(targetPointValues);
167  }
168 
169  COM_finalize();
170 
171  return 0;
172 }
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).
geoMeshBase * New(MeshType meshType)
Create a new mesh object.
void extractDataFromVTK(vtkDataSet *data, std::vector< double > &coords, std::vector< int > &elems)
meshBase * target
Definition: TransferBase.H:106
vtkSmartPointer< vtkDataSet > getDataSet() const
get this meshes&#39; dataSet
Definition: meshBase.H:308
meshBase * source
Definition: TransferBase.H:105

◆ transferPointData() [2/3]

int ConservativeSurfaceTransfer::transferPointData ( int  i,
vtkSmartPointer< vtkGenericCell >  genCell,
std::vector< vtkSmartPointer< vtkDoubleArray >> &  dasSource,
std::vector< vtkSmartPointer< vtkDoubleArray >> &  dasTarget,
bool  flip 
)
inline

Definition at line 60 of file ConservativeSurfaceTransfer.H.

63  {
64  return 0;
65  }

◆ transferPointData() [3/3]

int TransferBase::transferPointData ( const std::vector< std::string > &  arrayNames,
const std::vector< std::string > &  newnames = std::vector<std::string>() 
)
inherited
Parameters
arrayIDsarray of array names to specify which fields to transfer
newnamesoptional array of names to be applied to transferred fields
Returns
0 upon completion

Definition at line 36 of file TransferBase.C.

References TransferBase::getArrayIDs(), meshBase::getDataSet(), TransferBase::source, and TransferBase::transferPointData().

38 {
39  vtkPointData* pointData = source->getDataSet()->GetPointData();
40  std::vector<int> arrayIDs = getArrayIDs(arrayNames, pointData);
41 
42  return transferPointData(arrayIDs, newnames);
43 }
vtkSmartPointer< vtkDataSet > getDataSet() const
get this meshes&#39; dataSet
Definition: meshBase.H:308
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< int > getArrayIDs(const std::vector< std::string > &arrayNames, vtkFieldData *fieldData)
given array names and field data, return vector of corresponding array ids in the field data ...
Definition: TransferBase.C:56
meshBase * source
Definition: TransferBase.H:105

◆ writeOverlay()

int ConservativeSurfaceTransfer::writeOverlay ( )

Definition at line 174 of file ConservativeSurfaceTransfer.C.

References extractDataFromVTK(), meshBase::getDataSet(), RFC_overlay, RFC_write, TransferBase::source, and TransferBase::target.

175 {
176  // TODO Find a way to resolve scoping with COM
177  // this initialization should not need to happen in more than one function
178  COM_init(NULL, NULL);
179  COM_set_profiling(1);
180  // COM_set_verbose(10);
181 
182  MPI_Comm comm = MPI_COMM_WORLD;
183 
184  COM_LOAD_MODULE_STATIC_DYNAMIC(SurfX, "RFC");
185 
186  RFC_write = COM_get_function_handle("RFC.write_overlay");
187  RFC_overlay = COM_get_function_handle("RFC.overlay");
188 
189  vtkDataSet *sourceDataSet = source->getDataSet();
190  vtkDataSet *targetDataSet = target->getDataSet();
191 
192  std::vector<double> sourceCoords;
193  std::vector<int> sourceElems;
194  extractDataFromVTK(sourceDataSet, sourceCoords, sourceElems);
195 
196  // === SOURCE WINDOW ===
197  COM_new_window("sourceWindow");
198 
199  // set source nodal coordinates (nc)
200  COM_set_size("sourceWindow.nc", 1, sourceCoords.size() / 3);
201  COM_set_array("sourceWindow.nc", 1, &sourceCoords[0], 3);
202  // set source element node ids (:t3:)
203  COM_set_size("sourceWindow.:t3:", 1, sourceElems.size() / 3);
204  COM_set_array("sourceWindow.:t3:", 1, &sourceElems[0], 3);
205 
206  // initialize source soln array - this is reused
207  COM_new_dataitem("sourceWindow.soln", 'n', COM_DOUBLE, 1, "");
208 
209  COM_window_init_done("sourceWindow");
210 
211  // === TARGET WINDOW ===
212  std::vector<double> targetCoords;
213  std::vector<int> targetElems;
214  extractDataFromVTK(target->getDataSet(), targetCoords, targetElems);
215 
216  COM_new_window("targetWindow");
217  // set target nodal coordinates (nc)
218  COM_set_size("targetWindow.nc", 1, targetCoords.size() / 3);
219  COM_set_array("targetWindow.nc", 1, &targetCoords[0], 3);
220  // set target element node ids (:t3:)
221  COM_set_size("targetWindow.:t3:", 1, targetElems.size() / 3);
222  COM_set_array("targetWindow.:t3:", 1, &targetElems[0], 3);
223 
224  // initialize target soln array - this is reused
225  std::vector<double> targetData(targetCoords.size(), -1.);
226  COM_new_dataitem("targetWindow.soln", 'n', COM_DOUBLE, 1, "m/s");
227  COM_resize_array("targetWindow.soln");
228 
229  COM_window_init_done("targetWindow");
230 
231  int srcMesh = COM_get_dataitem_handle("sourceWindow.mesh");
232  int tgtMesh = COM_get_dataitem_handle("targetWindow.mesh");
233  COM_call_function(RFC_overlay, &srcMesh, &tgtMesh);
234 
235  COM_call_function(RFC_write, &srcMesh, &tgtMesh, "srcMesh", "tgtMesh", "CGNS");
236 
237  COM_finalize();
238 
239  return 0;
240 }
void extractDataFromVTK(vtkDataSet *data, std::vector< double > &coords, std::vector< int > &elems)
meshBase * target
Definition: TransferBase.H:106
vtkSmartPointer< vtkDataSet > getDataSet() const
get this meshes&#39; dataSet
Definition: meshBase.H:308
meshBase * source
Definition: TransferBase.H:105

Member Data Documentation

◆ c2cTrnsDistTol

double TransferBase::c2cTrnsDistTol
protectedinherited

Definition at line 116 of file TransferBase.H.

Referenced by FETransfer::transferCellData().

◆ checkQual

bool TransferBase::checkQual
protectedinherited

Definition at line 114 of file TransferBase.H.

Referenced by FETransfer::transferPointData().

◆ continuous

bool TransferBase::continuous
protectedinherited

Definition at line 115 of file TransferBase.H.

Referenced by FETransfer::transferCellData().

◆ RFC_clear

int ConservativeSurfaceTransfer::RFC_clear
private

Definition at line 96 of file ConservativeSurfaceTransfer.H.

Referenced by transferPointData().

◆ RFC_interp

int ConservativeSurfaceTransfer::RFC_interp
private

Definition at line 96 of file ConservativeSurfaceTransfer.H.

Referenced by transferPointData().

◆ RFC_overlay

int ConservativeSurfaceTransfer::RFC_overlay
private

Definition at line 96 of file ConservativeSurfaceTransfer.H.

Referenced by transferPointData(), and writeOverlay().

◆ RFC_read

int ConservativeSurfaceTransfer::RFC_read
private

Definition at line 96 of file ConservativeSurfaceTransfer.H.

Referenced by transferPointData().

◆ RFC_transfer

int ConservativeSurfaceTransfer::RFC_transfer
private

Definition at line 96 of file ConservativeSurfaceTransfer.H.

Referenced by transferPointData().

◆ RFC_write

int ConservativeSurfaceTransfer::RFC_write
private

Definition at line 96 of file ConservativeSurfaceTransfer.H.

Referenced by transferPointData(), and writeOverlay().

◆ source

◆ srcCellLocator

vtkSmartPointer<vtkStaticCellLocator> TransferBase::srcCellLocator = nullptr
protectedinherited

Definition at line 108 of file TransferBase.H.

Referenced by FETransfer::FETransfer(), and FETransfer::getClosestSourceCell().

◆ srcPointLocator

vtkSmartPointer<vtkStaticPointLocator> TransferBase::srcPointLocator = nullptr
protectedinherited

Definition at line 111 of file TransferBase.H.

Referenced by FETransfer::FETransfer().

◆ target

◆ trgCellLocator

vtkSmartPointer<vtkStaticCellLocator> TransferBase::trgCellLocator = nullptr
protectedinherited

Definition at line 109 of file TransferBase.H.

Referenced by FETransfer::FETransfer(), and FETransfer::getClosestSourceCell().

◆ trgPointLocator

vtkSmartPointer<vtkStaticPointLocator> TransferBase::trgPointLocator = nullptr
protectedinherited

Definition at line 112 of file TransferBase.H.

Referenced by FETransfer::FETransfer().


The documentation for this class was generated from the following files: