VTK  9.2.5
vtkCellArray.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkCellArray.h
5 
6  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7  All rights reserved.
8  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9 
10  This software is distributed WITHOUT ANY WARRANTY; without even
11  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12  PURPOSE. See the above copyright notice for more information.
13 
14 =========================================================================*/
250 #ifndef vtkCellArray_h
251 #define vtkCellArray_h
252 
253 #include "vtkCommonDataModelModule.h" // For export macro
254 #include "vtkObject.h"
255 
256 #include "vtkAOSDataArrayTemplate.h" // Needed for inline methods
257 #include "vtkCell.h" // Needed for inline methods
258 #include "vtkDataArrayRange.h" // Needed for inline methods
259 #include "vtkFeatures.h" // for VTK_USE_MEMKIND
260 #include "vtkSmartPointer.h" // For vtkSmartPointer
261 #include "vtkTypeInt32Array.h" // Needed for inline methods
262 #include "vtkTypeInt64Array.h" // Needed for inline methods
263 #include "vtkTypeList.h" // Needed for ArrayList definition
264 
265 #include <cassert> // for assert
266 #include <initializer_list> // for API
267 #include <type_traits> // for std::is_same
268 #include <utility> // for std::forward
269 
290 #define VTK_CELL_ARRAY_V2
291 
293 class vtkIdTypeArray;
294 
295 class VTKCOMMONDATAMODEL_EXPORT vtkCellArray : public vtkObject
296 {
297 public:
298  using ArrayType32 = vtkTypeInt32Array;
299  using ArrayType64 = vtkTypeInt64Array;
300 
302 
306  static vtkCellArray* New();
307  vtkTypeMacro(vtkCellArray, vtkObject);
308  void PrintSelf(ostream& os, vtkIndent indent) override;
309  void PrintDebug(ostream& os);
311 
320  using StorageArrayList = vtkTypeList::Create<ArrayType32, ArrayType64>;
321 
333 
342  vtkTypeBool Allocate(vtkIdType sz, vtkIdType vtkNotUsed(ext) = 1000)
343  {
344  return this->AllocateExact(sz, sz) ? 1 : 0;
345  }
346 
356  bool AllocateEstimate(vtkIdType numCells, vtkIdType maxCellSize)
357  {
358  return this->AllocateExact(numCells, numCells * maxCellSize);
359  }
360 
370  bool AllocateExact(vtkIdType numCells, vtkIdType connectivitySize);
371 
382  {
383  return this->AllocateExact(other->GetNumberOfCells(), other->GetNumberOfConnectivityIds());
384  }
385 
395  bool ResizeExact(vtkIdType numCells, vtkIdType connectivitySize);
396 
400  void Initialize();
401 
405  void Reset();
406 
412  void Squeeze();
413 
424  bool IsValid();
425 
430  {
431  if (this->Storage.Is64Bit())
432  {
433  return this->Storage.GetArrays64().Offsets->GetNumberOfValues() - 1;
434  }
435  else
436  {
437  return this->Storage.GetArrays32().Offsets->GetNumberOfValues() - 1;
438  }
439  }
440 
446  {
447  if (this->Storage.Is64Bit())
448  {
449  return this->Storage.GetArrays64().Offsets->GetNumberOfValues();
450  }
451  else
452  {
453  return this->Storage.GetArrays32().Offsets->GetNumberOfValues();
454  }
455  }
456 
461  {
462  if (this->Storage.Is64Bit())
463  {
464  return this->Storage.GetArrays64().Offsets->GetValue(cellId);
465  }
466  else
467  {
468  return this->Storage.GetArrays32().Offsets->GetValue(cellId);
469  }
470  }
471 
479  {
480  if (this->Storage.Is64Bit())
481  {
482  return this->Storage.GetArrays64().Connectivity->GetNumberOfValues();
483  }
484  else
485  {
486  return this->Storage.GetArrays32().Connectivity->GetNumberOfValues();
487  }
488  }
489 
496 
497 #ifndef __VTK_WRAP__ // The wrappers have issues with some of these templates
508  void SetData(vtkTypeInt32Array* offsets, vtkTypeInt32Array* connectivity);
509  void SetData(vtkTypeInt64Array* offsets, vtkTypeInt64Array* connectivity);
510  void SetData(vtkIdTypeArray* offsets, vtkIdTypeArray* connectivity);
513  void SetData(
516 #endif // __VTK_WRAP__
517 
530  bool SetData(vtkDataArray* offsets, vtkDataArray* connectivity);
531 
545  bool SetData(vtkIdType cellSize, vtkDataArray* connectivity);
546 
551  bool IsStorage64Bit() const { return this->Storage.Is64Bit(); }
552 
559  bool IsStorageShareable() const
560  {
561  if (this->Storage.Is64Bit())
562  {
564  }
565  else
566  {
568  }
569  }
570 
625  {
626  if (this->Storage.Is64Bit())
627  {
628  return this->GetOffsetsArray64();
629  }
630  else
631  {
632  return this->GetOffsetsArray32();
633  }
634  }
646  {
647  if (this->Storage.Is64Bit())
648  {
649  return this->GetConnectivityArray64();
650  }
651  else
652  {
653  return this->GetConnectivityArray32();
654  }
655  }
669 
679  void InitTraversal();
680 
695  int GetNextCell(vtkIdType& npts, vtkIdType const*& pts) VTK_SIZEHINT(pts, npts);
696 
707  int GetNextCell(vtkIdList* pts);
708 
719  void GetCellAtId(vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints)
720  VTK_SIZEHINT(cellPoints, cellSize) VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
721 
731  void GetCellAtId(
732  vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints, vtkIdList* ptIds)
733  VTK_SIZEHINT(cellPoints, cellSize) VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
734 
740  void GetCellAtId(vtkIdType cellId, vtkIdList* pts)
741  VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
742 
746  vtkIdType GetCellSize(const vtkIdType cellId) const;
747 
751  vtkIdType InsertNextCell(vtkCell* cell);
752 
757  vtkIdType InsertNextCell(vtkIdType npts, const vtkIdType* pts) VTK_SIZEHINT(pts, npts);
758 
763  vtkIdType InsertNextCell(vtkIdList* pts);
764 
772  vtkIdType InsertNextCell(const std::initializer_list<vtkIdType>& cell)
773  {
774  return this->InsertNextCell(static_cast<vtkIdType>(cell.size()), cell.begin());
775  }
776 
783  vtkIdType InsertNextCell(int npts);
784 
789  void InsertCellPoint(vtkIdType id);
790 
795  void UpdateCellCount(int npts);
796 
811  void ReverseCellAtId(vtkIdType cellId) VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
812 
821  void ReplaceCellAtId(vtkIdType cellId, vtkIdList* list);
822  void ReplaceCellAtId(vtkIdType cellId, vtkIdType cellSize, const vtkIdType* cellPoints)
823  VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells()) VTK_SIZEHINT(cellPoints, cellSize);
833  void ReplaceCellAtId(vtkIdType cellId, const std::initializer_list<vtkIdType>& cell)
834  {
835  return this->ReplaceCellAtId(cellId, static_cast<vtkIdType>(cell.size()), cell.begin());
836  }
837 
843 
848 
853 
857  void Append(vtkCellArray* src, vtkIdType pointOffset = 0);
858 
870 
898  void AppendLegacyFormat(const vtkIdType* data, vtkIdType len, vtkIdType ptOffset = 0)
899  VTK_SIZEHINT(data, len);
910  unsigned long GetActualMemorySize() const;
911 
912  // The following code is used to support
913 
914  // The wrappers get understandably confused by some of the template code below
915 #ifndef __VTK_WRAP__
916 
917  // Holds connectivity and offset arrays of the given ArrayType.
918  template <typename ArrayT>
919  struct VisitState
920  {
921  using ArrayType = ArrayT;
922  using ValueType = typename ArrayType::ValueType;
923  using CellRangeType = decltype(vtk::DataArrayValueRange<1>(std::declval<ArrayType>()));
924 
925  // We can't just use is_same here, since binary compatible representations
926  // (e.g. int and long) are distinct types. Instead, ensure that ValueType
927  // is a signed integer the same size as vtkIdType.
928  // If this value is true, ValueType pointers may be safely converted to
929  // vtkIdType pointers via reinterpret cast.
930  static constexpr bool ValueTypeIsSameAsIdType = std::is_integral<ValueType>::value &&
931  std::is_signed<ValueType>::value && (sizeof(ValueType) == sizeof(vtkIdType));
932 
933  ArrayType* GetOffsets() { return this->Offsets; }
934  const ArrayType* GetOffsets() const { return this->Offsets; }
935 
936  ArrayType* GetConnectivity() { return this->Connectivity; }
937  const ArrayType* GetConnectivity() const { return this->Connectivity; }
938 
940 
942 
944 
946 
948 
949  friend class vtkCellArray;
950 
951  protected:
953  {
954  this->Connectivity = vtkSmartPointer<ArrayType>::New();
955  this->Offsets = vtkSmartPointer<ArrayType>::New();
956  this->Offsets->InsertNextValue(0);
958  {
959  this->IsInMemkind = true;
960  }
961  }
962  ~VisitState() = default;
963  void* operator new(size_t nSize)
964  {
965  void* r;
966 #ifdef VTK_USE_MEMKIND
968 #else
969  r = malloc(nSize);
970 #endif
971  return r;
972  }
973  void operator delete(void* p)
974  {
975 #ifdef VTK_USE_MEMKIND
976  VisitState* a = static_cast<VisitState*>(p);
977  if (a->IsInMemkind)
978  {
980  }
981  else
982  {
983  free(p);
984  }
985 #else
986  free(p);
987 #endif
988  }
989 
992 
993  private:
994  VisitState(const VisitState&) = delete;
995  VisitState& operator=(const VisitState&) = delete;
996  bool IsInMemkind = false;
997  };
998 
999 private: // Helpers that allow Visit to return a value:
1000  template <typename Functor, typename... Args>
1001  using GetReturnType = decltype(
1002  std::declval<Functor>()(std::declval<VisitState<ArrayType32>&>(), std::declval<Args>()...));
1003 
1004  template <typename Functor, typename... Args>
1005  struct ReturnsVoid : std::is_same<GetReturnType<Functor, Args...>, void>
1006  {
1007  };
1008 
1009 public:
1079  template <typename Functor, typename... Args,
1080  typename = typename std::enable_if<ReturnsVoid<Functor, Args...>::value>::type>
1081  void Visit(Functor&& functor, Args&&... args)
1082  {
1083  if (this->Storage.Is64Bit())
1084  {
1085  // If you get an error on the next line, a call to Visit(functor, Args...)
1086  // is being called with arguments that do not match the functor's call
1087  // signature. See the Visit documentation for details.
1088  functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1089  }
1090  else
1091  {
1092  // If you get an error on the next line, a call to Visit(functor, Args...)
1093  // is being called with arguments that do not match the functor's call
1094  // signature. See the Visit documentation for details.
1095  functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1096  }
1097  }
1098 
1099  template <typename Functor, typename... Args,
1100  typename = typename std::enable_if<ReturnsVoid<Functor, Args...>::value>::type>
1101  void Visit(Functor&& functor, Args&&... args) const
1102  {
1103  if (this->Storage.Is64Bit())
1104  {
1105  // If you get an error on the next line, a call to Visit(functor, Args...)
1106  // is being called with arguments that do not match the functor's call
1107  // signature. See the Visit documentation for details.
1108  functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1109  }
1110  else
1111  {
1112  // If you get an error on the next line, a call to Visit(functor, Args...)
1113  // is being called with arguments that do not match the functor's call
1114  // signature. See the Visit documentation for details.
1115  functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1116  }
1117  }
1118 
1119  template <typename Functor, typename... Args,
1120  typename = typename std::enable_if<!ReturnsVoid<Functor, Args...>::value>::type>
1121  GetReturnType<Functor, Args...> Visit(Functor&& functor, Args&&... args)
1122  {
1123  if (this->Storage.Is64Bit())
1124  {
1125  // If you get an error on the next line, a call to Visit(functor, Args...)
1126  // is being called with arguments that do not match the functor's call
1127  // signature. See the Visit documentation for details.
1128  return functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1129  }
1130  else
1131  {
1132  // If you get an error on the next line, a call to Visit(functor, Args...)
1133  // is being called with arguments that do not match the functor's call
1134  // signature. See the Visit documentation for details.
1135  return functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1136  }
1137  }
1138  template <typename Functor, typename... Args,
1139  typename = typename std::enable_if<!ReturnsVoid<Functor, Args...>::value>::type>
1140  GetReturnType<Functor, Args...> Visit(Functor&& functor, Args&&... args) const
1141  {
1142  if (this->Storage.Is64Bit())
1143  {
1144  // If you get an error on the next line, a call to Visit(functor, Args...)
1145  // is being called with arguments that do not match the functor's call
1146  // signature. See the Visit documentation for details.
1147  return functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1148  }
1149  else
1150  {
1151  // If you get an error on the next line, a call to Visit(functor, Args...)
1152  // is being called with arguments that do not match the functor's call
1153  // signature. See the Visit documentation for details.
1154  return functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1155  }
1156  }
1157 
1160 #endif // __VTK_WRAP__
1161 
1162  //=================== Begin Legacy Methods ===================================
1163  // These should be deprecated at some point as they are confusing or very slow
1164 
1172 
1184  vtkIdType EstimateSize(vtkIdType numCells, int maxPtsPerCell);
1185 
1195 
1203 
1213  void GetCell(vtkIdType loc, vtkIdType& npts, const vtkIdType*& pts)
1214  VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries()) VTK_SIZEHINT(pts, npts);
1215 
1222  void GetCell(vtkIdType loc, vtkIdList* pts)
1223  VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries());
1224 
1231  vtkIdType GetInsertLocation(int npts);
1232 
1240  vtkIdType GetTraversalLocation();
1241  vtkIdType GetTraversalLocation(vtkIdType npts);
1242  void SetTraversalLocation(vtkIdType loc);
1252  void ReverseCell(vtkIdType loc) VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries());
1253 
1265  void ReplaceCell(vtkIdType loc, int npts, const vtkIdType pts[])
1266  VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries()) VTK_SIZEHINT(pts, npts);
1267 
1282  void SetCells(vtkIdType ncells, vtkIdTypeArray* cells);
1283 
1295 
1296  //=================== End Legacy Methods =====================================
1297 
1298  friend class vtkCellArrayIterator;
1299 
1300 protected:
1302  ~vtkCellArray() override;
1303 
1304  // Encapsulates storage of the internal arrays as a discriminated union
1305  // between 32-bit and 64-bit storage.
1306  struct Storage
1307  {
1308  // Union type that switches 32 and 64 bit array storage
1309  union ArraySwitch {
1310  ArraySwitch() = default; // handled by Storage
1311  ~ArraySwitch() = default; // handle by Storage
1314  };
1315 
1317  {
1318 #ifdef VTK_USE_MEMKIND
1319  this->Arrays =
1321 #else
1322  this->Arrays = new ArraySwitch;
1323 #endif
1324 
1325  // Default to the compile-time setting:
1326 #ifdef VTK_USE_64BIT_IDS
1327 
1328  this->Arrays->Int64 = new VisitState<ArrayType64>;
1329  this->StorageIs64Bit = true;
1330 
1331 #else // VTK_USE_64BIT_IDS
1332 
1333  this->Arrays->Int32 = new VisitState<ArrayType32>;
1334  this->StorageIs64Bit = false;
1335 
1336 #endif // VTK_USE_64BIT_IDS
1337 #ifdef VTK_USE_MEMKIND
1339  {
1340  this->IsInMemkind = true;
1341  }
1342 #else
1343  (void)this->IsInMemkind; // comp warning workaround
1344 #endif
1345  }
1346 
1348  {
1349  if (this->StorageIs64Bit)
1350  {
1351  this->Arrays->Int64->~VisitState();
1352  delete this->Arrays->Int64;
1353  }
1354  else
1355  {
1356  this->Arrays->Int32->~VisitState();
1357  delete this->Arrays->Int32;
1358  }
1359 #ifdef VTK_USE_MEMKIND
1360  if (this->IsInMemkind)
1361  {
1363  }
1364  else
1365  {
1366  free(this->Arrays);
1367  }
1368 #else
1369  delete this->Arrays;
1370 #endif
1371  }
1372 
1373  // Switch the internal arrays to be 32-bit. Any old data is lost. Returns
1374  // true if the storage changes.
1376  {
1377  if (!this->StorageIs64Bit)
1378  {
1379  return false;
1380  }
1381 
1382  this->Arrays->Int64->~VisitState();
1383  delete this->Arrays->Int64;
1384  this->Arrays->Int32 = new VisitState<ArrayType32>;
1385  this->StorageIs64Bit = false;
1386 
1387  return true;
1388  }
1389 
1390  // Switch the internal arrays to be 64-bit. Any old data is lost. Returns
1391  // true if the storage changes.
1393  {
1394  if (this->StorageIs64Bit)
1395  {
1396  return false;
1397  }
1398 
1399  this->Arrays->Int32->~VisitState();
1400  delete this->Arrays->Int32;
1401  this->Arrays->Int64 = new VisitState<ArrayType64>;
1402  this->StorageIs64Bit = true;
1403 
1404  return true;
1405  }
1406 
1407  // Returns true if the storage is currently configured to be 64 bit.
1408  bool Is64Bit() const { return this->StorageIs64Bit; }
1409 
1410  // Get the VisitState for 32-bit arrays
1412  {
1413  assert(!this->StorageIs64Bit);
1414  return *this->Arrays->Int32;
1415  }
1416 
1418  {
1419  assert(!this->StorageIs64Bit);
1420  return *this->Arrays->Int32;
1421  }
1422 
1423  // Get the VisitState for 64-bit arrays
1425  {
1426  assert(this->StorageIs64Bit);
1427  return *this->Arrays->Int64;
1428  }
1429 
1431  {
1432  assert(this->StorageIs64Bit);
1433  return *this->Arrays->Int64;
1434  }
1435 
1436  private:
1437  // Access restricted to ensure proper union construction/destruction thru
1438  // API.
1439  ArraySwitch* Arrays;
1440  bool StorageIs64Bit;
1441  bool IsInMemkind = false;
1442  };
1443 
1446  vtkIdType TraversalCellId{ 0 };
1447 
1449 
1450 private:
1451  vtkCellArray(const vtkCellArray&) = delete;
1452  void operator=(const vtkCellArray&) = delete;
1453 };
1454 
1455 template <typename ArrayT>
1457 {
1458  return this->Offsets->GetNumberOfValues() - 1;
1459 }
1460 
1461 template <typename ArrayT>
1463 {
1464  return static_cast<vtkIdType>(this->Offsets->GetValue(cellId));
1465 }
1466 
1467 template <typename ArrayT>
1469 {
1470  return static_cast<vtkIdType>(this->Offsets->GetValue(cellId + 1));
1471 }
1472 
1473 template <typename ArrayT>
1475 {
1476  return this->GetEndOffset(cellId) - this->GetBeginOffset(cellId);
1477 }
1478 
1479 template <typename ArrayT>
1482 {
1483  return vtk::DataArrayValueRange<1>(
1484  this->GetConnectivity(), this->GetBeginOffset(cellId), this->GetEndOffset(cellId));
1485 }
1486 
1488 {
1489 
1491 {
1492  // Insert full cell
1493  template <typename CellStateT>
1494  vtkIdType operator()(CellStateT& state, const vtkIdType npts, const vtkIdType pts[])
1495  {
1496  using ValueType = typename CellStateT::ValueType;
1497  auto* conn = state.GetConnectivity();
1498  auto* offsets = state.GetOffsets();
1499 
1500  const vtkIdType cellId = offsets->GetNumberOfValues() - 1;
1501 
1502  offsets->InsertNextValue(static_cast<ValueType>(conn->GetNumberOfValues() + npts));
1503 
1504  for (vtkIdType i = 0; i < npts; ++i)
1505  {
1506  conn->InsertNextValue(static_cast<ValueType>(pts[i]));
1507  }
1508 
1509  return cellId;
1510  }
1511 
1512  // Just update offset table (for incremental API)
1513  template <typename CellStateT>
1514  vtkIdType operator()(CellStateT& state, const vtkIdType npts)
1515  {
1516  using ValueType = typename CellStateT::ValueType;
1517  auto* conn = state.GetConnectivity();
1518  auto* offsets = state.GetOffsets();
1519 
1520  const vtkIdType cellId = offsets->GetNumberOfValues() - 1;
1521 
1522  offsets->InsertNextValue(static_cast<ValueType>(conn->GetNumberOfValues() + npts));
1523 
1524  return cellId;
1525  }
1526 };
1527 
1528 // for incremental API:
1530 {
1531  template <typename CellStateT>
1532  void operator()(CellStateT& state, const vtkIdType npts)
1533  {
1534  using ValueType = typename CellStateT::ValueType;
1535 
1536  auto* offsets = state.GetOffsets();
1537  const ValueType cellBegin = offsets->GetValue(offsets->GetMaxId() - 1);
1538  offsets->SetValue(offsets->GetMaxId(), static_cast<ValueType>(cellBegin + npts));
1539  }
1540 };
1541 
1543 {
1544  template <typename CellStateT>
1545  vtkIdType operator()(CellStateT& state, const vtkIdType cellId)
1546  {
1547  return state.GetCellSize(cellId);
1548  }
1549 };
1550 
1552 {
1553  template <typename CellStateT>
1554  void operator()(CellStateT& state, const vtkIdType cellId, vtkIdList* ids)
1555  {
1556  using ValueType = typename CellStateT::ValueType;
1557 
1558  const auto cellPts = state.GetCellRange(cellId);
1559 
1560  ids->SetNumberOfIds(cellPts.size());
1561  vtkIdType* idPtr = ids->GetPointer(0);
1562 
1563  for (ValueType ptId : cellPts)
1564  {
1565  *idPtr++ = static_cast<vtkIdType>(ptId);
1566  }
1567  }
1568 
1569  // SFINAE helper to check if a VisitState's connectivity array's memory
1570  // can be used as a vtkIdType*.
1571  template <typename CellStateT>
1573  {
1574  private:
1575  using ValueType = typename CellStateT::ValueType;
1576  using ArrayType = typename CellStateT::ArrayType;
1578  static constexpr bool ValueTypeCompat = CellStateT::ValueTypeIsSameAsIdType;
1579  static constexpr bool ArrayTypeCompat = std::is_base_of<AOSArrayType, ArrayType>::value;
1580 
1581  public:
1582  static constexpr bool value = ValueTypeCompat && ArrayTypeCompat;
1583  };
1584 
1585  template <typename CellStateT>
1587  CellStateT& state, const vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints,
1588  vtkIdList* vtkNotUsed(temp))
1589  {
1590  const vtkIdType beginOffset = state.GetBeginOffset(cellId);
1591  const vtkIdType endOffset = state.GetEndOffset(cellId);
1592  cellSize = endOffset - beginOffset;
1593  // This is safe, see CanShareConnPtr helper above.
1594  cellPoints = reinterpret_cast<vtkIdType*>(state.GetConnectivity()->GetPointer(beginOffset));
1595  }
1596 
1597  template <typename CellStateT>
1599  CellStateT& state, const vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints,
1600  vtkIdList* temp)
1601  {
1602  using ValueType = typename CellStateT::ValueType;
1603 
1604  const auto cellPts = state.GetCellRange(cellId);
1605  cellSize = cellPts.size();
1606 
1607  // ValueType differs from vtkIdType, so we have to copy into a temporary
1608  // buffer:
1609  temp->SetNumberOfIds(cellSize);
1610  vtkIdType* tempPtr = temp->GetPointer(0);
1611  for (ValueType ptId : cellPts)
1612  {
1613  *tempPtr++ = static_cast<vtkIdType>(ptId);
1614  }
1615 
1616  cellPoints = temp->GetPointer(0);
1617  }
1618 };
1619 
1621 {
1622  template <typename CellStateT>
1623  void operator()(CellStateT& state)
1624  {
1625  state.GetOffsets()->Reset();
1626  state.GetConnectivity()->Reset();
1627  state.GetOffsets()->InsertNextValue(0);
1628  }
1629 };
1630 
1631 } // end namespace vtkCellArray_detail
1632 
1633 //----------------------------------------------------------------------------
1635 {
1636  this->TraversalCellId = 0;
1637 }
1638 
1639 //----------------------------------------------------------------------------
1640 inline int vtkCellArray::GetNextCell(vtkIdType& npts, vtkIdType const*& pts) VTK_SIZEHINT(pts, npts)
1641 {
1642  if (this->TraversalCellId < this->GetNumberOfCells())
1643  {
1644  this->GetCellAtId(this->TraversalCellId, npts, pts);
1645  ++this->TraversalCellId;
1646  return 1;
1647  }
1648 
1649  npts = 0;
1650  pts = nullptr;
1651  return 0;
1652 }
1653 
1654 //----------------------------------------------------------------------------
1656 {
1657  if (this->TraversalCellId < this->GetNumberOfCells())
1658  {
1659  this->GetCellAtId(this->TraversalCellId, pts);
1660  ++this->TraversalCellId;
1661  return 1;
1662  }
1663 
1664  pts->Reset();
1665  return 0;
1666 }
1667 //----------------------------------------------------------------------------
1669 {
1670  return this->Visit(vtkCellArray_detail::GetCellSizeImpl{}, cellId);
1671 }
1672 
1673 //----------------------------------------------------------------------------
1674 inline void vtkCellArray::GetCellAtId(vtkIdType cellId, vtkIdType& cellSize,
1675  vtkIdType const*& cellPoints) VTK_SIZEHINT(cellPoints, cellSize)
1676 {
1677  this->Visit(vtkCellArray_detail::GetCellAtIdImpl{}, cellId, cellSize, cellPoints, this->TempCell);
1678 }
1679 
1680 //----------------------------------------------------------------------------
1681 inline void vtkCellArray::GetCellAtId(vtkIdType cellId, vtkIdType& cellSize,
1682  vtkIdType const*& cellPoints, vtkIdList* ptIds) VTK_SIZEHINT(cellPoints, cellSize)
1683 {
1684  this->Visit(vtkCellArray_detail::GetCellAtIdImpl{}, cellId, cellSize, cellPoints, ptIds);
1685 }
1686 
1687 //----------------------------------------------------------------------------
1689 {
1690  this->Visit(vtkCellArray_detail::GetCellAtIdImpl{}, cellId, pts);
1691 }
1692 
1693 //----------------------------------------------------------------------------
1695  VTK_SIZEHINT(pts, npts)
1696 {
1697  return this->Visit(vtkCellArray_detail::InsertNextCellImpl{}, npts, pts);
1698 }
1699 
1700 //----------------------------------------------------------------------------
1702 {
1703  return this->Visit(vtkCellArray_detail::InsertNextCellImpl{}, npts);
1704 }
1705 
1706 //----------------------------------------------------------------------------
1708 {
1709  if (this->Storage.Is64Bit())
1710  {
1711  using ValueType = typename ArrayType64::ValueType;
1712  this->Storage.GetArrays64().Connectivity->InsertNextValue(static_cast<ValueType>(id));
1713  }
1714  else
1715  {
1716  using ValueType = typename ArrayType32::ValueType;
1717  this->Storage.GetArrays32().Connectivity->InsertNextValue(static_cast<ValueType>(id));
1718  }
1719 }
1720 
1721 //----------------------------------------------------------------------------
1722 inline void vtkCellArray::UpdateCellCount(int npts)
1723 {
1725 }
1726 
1727 //----------------------------------------------------------------------------
1729 {
1730  return this->Visit(
1732 }
1733 
1734 //----------------------------------------------------------------------------
1736 {
1737  vtkIdList* pts = cell->GetPointIds();
1738  return this->Visit(
1740 }
1741 
1742 //----------------------------------------------------------------------------
1743 inline void vtkCellArray::Reset()
1744 {
1746 }
1747 
1748 #endif // vtkCellArray.h
Encapsulate traversal logic for vtkCellArray.
object to represent cell connectivity
Definition: vtkCellArray.h:296
void SetData(vtkAOSDataArrayTemplate< int > *offsets, vtkAOSDataArrayTemplate< int > *connectivity)
Set the internal data arrays to the supplied offsets and connectivity arrays.
int GetNextCell(vtkIdType &npts, vtkIdType const *&pts)
vtkIdType GetNumberOfConnectivityEntries()
Return the size of the array that would be returned from ExportLegacyFormat().
vtkTypeBool Allocate(vtkIdType sz, vtkIdType vtkNotUsed(ext)=1000)
Allocate memory.
Definition: vtkCellArray.h:342
void UseDefaultStorage()
Initialize internal data structures to use 32- or 64-bit storage.
Storage Storage
vtkIdType GetOffset(vtkIdType cellId)
Get the offset (into the connectivity) for a specified cell id.
Definition: vtkCellArray.h:460
bool AllocateCopy(vtkCellArray *other)
Pre-allocate memory in internal data structures to match the used size of the input vtkCellArray.
Definition: vtkCellArray.h:381
void AppendLegacyFormat(vtkIdTypeArray *data, vtkIdType ptOffset=0)
Append an array of data with the legacy vtkCellArray layout, e.g.
bool IsStorageShareable() const
Definition: vtkCellArray.h:559
ArrayType32 * GetOffsetsArray32()
Return the array used to store cell offsets.
Definition: vtkCellArray.h:635
bool IsValid()
Check that internal storage is consistent and in a valid state.
void AppendLegacyFormat(const vtkIdType *data, vtkIdType len, vtkIdType ptOffset=0)
Append an array of data with the legacy vtkCellArray layout, e.g.
bool SetData(vtkIdType cellSize, vtkDataArray *connectivity)
Sets the internal arrays to the supported connectivity array with an offsets array automatically gene...
vtkIdType GetTraversalCellId()
Get/Set the current cellId for traversal.
void Visit(Functor &&functor, Args &&... args)
vtkTypeInt32Array ArrayType32
Definition: vtkCellArray.h:298
void Reset()
Reuse list.
ArrayType64 * GetOffsetsArray64()
Return the array used to store cell offsets.
Definition: vtkCellArray.h:636
void SetData(vtkAOSDataArrayTemplate< long long > *offsets, vtkAOSDataArrayTemplate< long long > *connectivity)
Set the internal data arrays to the supplied offsets and connectivity arrays.
GetReturnType< Functor, Args... > Visit(Functor &&functor, Args &&... args)
bool AllocateExact(vtkIdType numCells, vtkIdType connectivitySize)
Pre-allocate memory in internal data structures.
bool ResizeExact(vtkIdType numCells, vtkIdType connectivitySize)
ResizeExact() resizes the internal structures to hold numCells total cell offsets and connectivitySiz...
vtkIdType EstimateSize(vtkIdType numCells, int maxPtsPerCell)
Utility routines help manage memory of cell array.
bool AllocateEstimate(vtkIdType numCells, vtkIdType maxCellSize)
Pre-allocate memory in internal data structures.
Definition: vtkCellArray.h:356
vtkIdType GetNumberOfOffsets() const
Get the number of elements in the offsets array.
Definition: vtkCellArray.h:445
vtkTypeInt64Array ArrayType64
Definition: vtkCellArray.h:299
vtkIdType TraversalCellId
void InitTraversal()
void Use64BitStorage()
Initialize internal data structures to use 32- or 64-bit storage.
bool ConvertToDefaultStorage()
Convert internal data structures to use 32- or 64-bit storage.
bool CanConvertToDefaultStorage() const
Check if the existing data can safely be converted to use 32- or 64- bit storage.
void GetCell(vtkIdType loc, vtkIdType &npts, const vtkIdType *&pts)
Internal method used to retrieve a cell given a legacy offset location.
bool CanConvertTo32BitStorage() const
Check if the existing data can safely be converted to use 32- or 64- bit storage.
void InsertCellPoint(vtkIdType id)
Used in conjunction with InsertNextCell(npts) to add another point to the list of cells.
void GetCellAtId(vtkIdType cellId, vtkIdType &cellSize, vtkIdType const *&cellPoints)
Return the point ids for the cell at cellId.
void ReplaceCellAtId(vtkIdType cellId, vtkIdList *list)
Replaces the point ids for the specified cell with the supplied list.
void ReverseCellAtId(vtkIdType cellId)
Reverses the order of the point ids for the specified cell.
bool SetData(vtkDataArray *offsets, vtkDataArray *connectivity)
Sets the internal arrays to the supplied offsets and connectivity arrays.
void Squeeze()
Reclaim any extra memory while preserving data.
void SetData(vtkIdTypeArray *offsets, vtkIdTypeArray *connectivity)
Set the internal data arrays to the supplied offsets and connectivity arrays.
static vtkCellArray * New()
Standard methods for instantiation, type information, and printing.
ArrayType32 * GetConnectivityArray32()
Return the array used to store the point ids that define the cells' connectivity.
Definition: vtkCellArray.h:656
ArrayType64 * GetConnectivityArray64()
Return the array used to store the point ids that define the cells' connectivity.
Definition: vtkCellArray.h:657
bool ConvertTo32BitStorage()
Convert internal data structures to use 32- or 64-bit storage.
vtkIdType IsHomogeneous()
Check if all cells have the same number of vertices.
void Visit(Functor &&functor, Args &&... args) const
int GetMaxCellSize()
Returns the size of the largest cell.
bool ConvertToSmallestStorage()
Convert internal data structures to use 32- or 64-bit storage.
void ShallowCopy(vtkCellArray *ca)
Shallow copy ca into this cell array.
void ExportLegacyFormat(vtkIdTypeArray *data)
Fill data with the old-style vtkCellArray data layout, e.g.
bool ConvertTo64BitStorage()
Convert internal data structures to use 32- or 64-bit storage.
void SetData(vtkAOSDataArrayTemplate< long > *offsets, vtkAOSDataArrayTemplate< long > *connectivity)
Set the internal data arrays to the supplied offsets and connectivity arrays.
void ImportLegacyFormat(const vtkIdType *data, vtkIdType len)
Import an array of data with the legacy vtkCellArray layout, e.g.
void ImportLegacyFormat(vtkIdTypeArray *data)
Import an array of data with the legacy vtkCellArray layout, e.g.
void Use32BitStorage()
Initialize internal data structures to use 32- or 64-bit storage.
void Initialize()
Free any memory and reset to an empty state.
void DeepCopy(vtkCellArray *ca)
Perform a deep copy (no reference counting) of the given cell array.
typename vtkTypeList::Unique< vtkTypeList::Create< vtkAOSDataArrayTemplate< int >, vtkAOSDataArrayTemplate< long >, vtkAOSDataArrayTemplate< long long > >>::Result InputArrayList
List of possible ArrayTypes that are compatible with internal storage.
Definition: vtkCellArray.h:332
vtkDataArray * GetConnectivityArray()
Return the array used to store the point ids that define the cells' connectivity.
Definition: vtkCellArray.h:645
void SetTraversalCellId(vtkIdType cellId)
Get/Set the current cellId for traversal.
void SetData(vtkTypeInt64Array *offsets, vtkTypeInt64Array *connectivity)
Set the internal data arrays to the supplied offsets and connectivity arrays.
vtkIdType GetNumberOfCells() const
Get the number of cells in the array.
Definition: vtkCellArray.h:429
vtkIdType InsertNextCell(vtkCell *cell)
Insert a cell object.
virtual void SetNumberOfCells(vtkIdType)
Set the number of cells in the array.
vtkIdType GetNumberOfConnectivityIds() const
Get the size of the connectivity array that stores the point ids.
Definition: vtkCellArray.h:478
vtkNew< vtkIdTypeArray > LegacyData
void PrintSelf(ostream &os, vtkIndent indent) override
Standard methods for instantiation, type information, and printing.
void PrintDebug(ostream &os)
Standard methods for instantiation, type information, and printing.
vtkCellArrayIterator * NewIterator()
NewIterator returns a new instance of vtkCellArrayIterator that is initialized to point at the first ...
void UpdateCellCount(int npts)
Used in conjunction with InsertNextCell(int npts) and InsertCellPoint() to update the number of point...
GetReturnType< Functor, Args... > Visit(Functor &&functor, Args &&... args) const
vtkNew< vtkIdList > TempCell
vtkIdType GetSize()
Get the size of the allocated connectivity array.
vtkDataArray * GetOffsetsArray()
Return the array used to store cell offsets.
Definition: vtkCellArray.h:624
vtkTypeList::Create< ArrayType32, ArrayType64 > StorageArrayList
List of possible array types used for storage.
Definition: vtkCellArray.h:320
vtkIdType InsertNextCell(const std::initializer_list< vtkIdType > &cell)
Overload that allows InsertNextCell({0, 1, 2}) syntax.
Definition: vtkCellArray.h:772
void Append(vtkCellArray *src, vtkIdType pointOffset=0)
Append cells from src into this.
vtkIdType GetCellSize(const vtkIdType cellId) const
Return the size of the cell at cellId.
bool IsStorage64Bit() const
Definition: vtkCellArray.h:551
void ReplaceCellAtId(vtkIdType cellId, vtkIdType cellSize, const vtkIdType *cellPoints)
Replaces the point ids for the specified cell with the supplied list.
bool CanConvertTo64BitStorage() const
Check if the existing data can safely be converted to use 32- or 64- bit storage.
abstract class to specify cell behavior
Definition: vtkCell.h:150
vtkIdList * GetPointIds()
Return the list of point ids defining the cell.
Definition: vtkCell.h:245
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:165
list of point or cell ids
Definition: vtkIdList.h:143
vtkIdType * GetPointer(const vtkIdType i)
Get a pointer to a particular data index.
Definition: vtkIdList.h:238
vtkIdType GetNumberOfIds() const noexcept
Return the number of id's in the list.
Definition: vtkIdList.h:169
void Reset()
Reset to an empty state but retain previously allocated memory.
Definition: vtkIdList.h:257
void SetNumberOfIds(const vtkIdType number)
Specify the number of ids for this object to hold.
dynamic, self-adjusting array of vtkIdType
a simple class to control print indentation
Definition: vtkIndent.h:119
static vtkMallocingFunction GetCurrentMallocFunction()
static vtkFreeingFunction GetAlternateFreeFunction()
static bool GetUsingMemkind()
A global state flag that controls whether vtkObjects are constructed in the usual way (the default) o...
abstract base class for most VTK objects
Definition: vtkObject.h:82
static vtkSmartPointer< T > New()
Create an instance of a VTK object.
vtkSmartPointer< vtkDataArray > GetData(const Ioss::GroupingEntity *entity, const std::string &fieldname, Ioss::Transform *transform=nullptr, Cache *cache=nullptr, const std::string &cachekey=std::string())
Returns a VTK array for a given field (fieldname) on the chosen block (or set) entity.
vtkSmartPointer< vtkCellArray > GetConnectivity(Ioss::GroupingEntity *group_entity, int &vtk_topology_type, Cache *cache=nullptr)
Read connectivity information from the group_entity.
@ value
Definition: vtkX3D.h:226
@ type
Definition: vtkX3D.h:522
@ data
Definition: vtkX3D.h:321
VisitState< ArrayType32 > & GetArrays32()
const VisitState< ArrayType64 > & GetArrays64() const
const VisitState< ArrayType32 > & GetArrays32() const
VisitState< ArrayType64 > & GetArrays64()
ArrayType * GetOffsets()
Definition: vtkCellArray.h:933
vtkIdType GetNumberOfCells() const
vtkIdType GetEndOffset(vtkIdType cellId) const
vtkSmartPointer< ArrayType > Offsets
Definition: vtkCellArray.h:991
vtkSmartPointer< ArrayType > Connectivity
Definition: vtkCellArray.h:990
static constexpr bool ValueTypeIsSameAsIdType
Definition: vtkCellArray.h:930
ArrayType * GetConnectivity()
Definition: vtkCellArray.h:936
decltype(vtk::DataArrayValueRange< 1 >(std::declval< ArrayType >())) CellRangeType
Definition: vtkCellArray.h:923
CellRangeType GetCellRange(vtkIdType cellId)
const ArrayType * GetOffsets() const
Definition: vtkCellArray.h:934
vtkIdType GetBeginOffset(vtkIdType cellId) const
vtkIdType GetCellSize(vtkIdType cellId) const
const ArrayType * GetConnectivity() const
Definition: vtkCellArray.h:937
typename ArrayType::ValueType ValueType
Definition: vtkCellArray.h:922
void operator()(CellStateT &state, const vtkIdType cellId, vtkIdList *ids)
std::enable_if<!CanShareConnPtr< CellStateT >::value, void >::type operator()(CellStateT &state, const vtkIdType cellId, vtkIdType &cellSize, vtkIdType const *&cellPoints, vtkIdList *temp)
std::enable_if< CanShareConnPtr< CellStateT >::value, void >::type operator()(CellStateT &state, const vtkIdType cellId, vtkIdType &cellSize, vtkIdType const *&cellPoints, vtkIdList *vtkNotUsed(temp))
vtkIdType operator()(CellStateT &state, const vtkIdType cellId)
vtkIdType operator()(CellStateT &state, const vtkIdType npts, const vtkIdType pts[])
vtkIdType operator()(CellStateT &state, const vtkIdType npts)
void operator()(CellStateT &state)
void operator()(CellStateT &state, const vtkIdType npts)
Remove all duplicate types from TypeList TList, storing the new list in Result.
Definition: vtkTypeList.h:126
VisitState< ArrayType32 > * Int32
VisitState< ArrayType64 > * Int64
int vtkTypeBool
Definition: vtkABI.h:69
STL-compatible iterable ranges that provide access to vtkDataArray elements.
int vtkIdType
Definition: vtkType.h:332
#define VTK_SIZEHINT(...)
#define VTK_EXPECTS(x)
#define VTK_NEWINSTANCE