VTK  9.2.20230207
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 
292 VTK_ABI_NAMESPACE_BEGIN
294 class vtkIdTypeArray;
295 
296 class VTKCOMMONDATAMODEL_EXPORT vtkCellArray : public vtkObject
297 {
298 public:
299  using ArrayType32 = vtkTypeInt32Array;
300  using ArrayType64 = vtkTypeInt64Array;
301 
303 
307  static vtkCellArray* New();
308  vtkTypeMacro(vtkCellArray, vtkObject);
309  void PrintSelf(ostream& os, vtkIndent indent) override;
310  void PrintDebug(ostream& os);
312 
321  using StorageArrayList = vtkTypeList::Create<ArrayType32, ArrayType64>;
322 
334 
343  vtkTypeBool Allocate(vtkIdType sz, vtkIdType vtkNotUsed(ext) = 1000)
344  {
345  return this->AllocateExact(sz, sz) ? 1 : 0;
346  }
347 
357  bool AllocateEstimate(vtkIdType numCells, vtkIdType maxCellSize)
358  {
359  return this->AllocateExact(numCells, numCells * maxCellSize);
360  }
361 
371  bool AllocateExact(vtkIdType numCells, vtkIdType connectivitySize);
372 
383  {
384  return this->AllocateExact(other->GetNumberOfCells(), other->GetNumberOfConnectivityIds());
385  }
386 
396  bool ResizeExact(vtkIdType numCells, vtkIdType connectivitySize);
397 
401  void Initialize();
402 
406  void Reset();
407 
413  void Squeeze();
414 
425  bool IsValid();
426 
431  {
432  if (this->Storage.Is64Bit())
433  {
434  return this->Storage.GetArrays64().Offsets->GetNumberOfValues() - 1;
435  }
436  else
437  {
438  return this->Storage.GetArrays32().Offsets->GetNumberOfValues() - 1;
439  }
440  }
441 
447  {
448  if (this->Storage.Is64Bit())
449  {
450  return this->Storage.GetArrays64().Offsets->GetNumberOfValues();
451  }
452  else
453  {
454  return this->Storage.GetArrays32().Offsets->GetNumberOfValues();
455  }
456  }
457 
462  {
463  if (this->Storage.Is64Bit())
464  {
465  return this->Storage.GetArrays64().Offsets->GetValue(cellId);
466  }
467  else
468  {
469  return this->Storage.GetArrays32().Offsets->GetValue(cellId);
470  }
471  }
472 
480  {
481  if (this->Storage.Is64Bit())
482  {
483  return this->Storage.GetArrays64().Connectivity->GetNumberOfValues();
484  }
485  else
486  {
487  return this->Storage.GetArrays32().Connectivity->GetNumberOfValues();
488  }
489  }
490 
497 
498 #ifndef __VTK_WRAP__ // The wrappers have issues with some of these templates
509  void SetData(vtkTypeInt32Array* offsets, vtkTypeInt32Array* connectivity);
510  void SetData(vtkTypeInt64Array* offsets, vtkTypeInt64Array* connectivity);
511  void SetData(vtkIdTypeArray* offsets, vtkIdTypeArray* connectivity);
514  void SetData(
517 #endif // __VTK_WRAP__
518 
531  bool SetData(vtkDataArray* offsets, vtkDataArray* connectivity);
532 
546  bool SetData(vtkIdType cellSize, vtkDataArray* connectivity);
547 
552  bool IsStorage64Bit() const { return this->Storage.Is64Bit(); }
553 
560  bool IsStorageShareable() const
561  {
562  if (this->Storage.Is64Bit())
563  {
565  }
566  else
567  {
569  }
570  }
571 
626  {
627  if (this->Storage.Is64Bit())
628  {
629  return this->GetOffsetsArray64();
630  }
631  else
632  {
633  return this->GetOffsetsArray32();
634  }
635  }
647  {
648  if (this->Storage.Is64Bit())
649  {
650  return this->GetConnectivityArray64();
651  }
652  else
653  {
654  return this->GetConnectivityArray32();
655  }
656  }
670 
680  void InitTraversal();
681 
696  int GetNextCell(vtkIdType& npts, vtkIdType const*& pts) VTK_SIZEHINT(pts, npts);
697 
708  int GetNextCell(vtkIdList* pts);
709 
720  void GetCellAtId(vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints)
721  VTK_SIZEHINT(cellPoints, cellSize) VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
722 
732  void GetCellAtId(
733  vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints, vtkIdList* ptIds)
734  VTK_SIZEHINT(cellPoints, cellSize) VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
735 
741  void GetCellAtId(vtkIdType cellId, vtkIdList* pts)
742  VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
743 
747  vtkIdType GetCellSize(const vtkIdType cellId) const;
748 
752  vtkIdType InsertNextCell(vtkCell* cell);
753 
758  vtkIdType InsertNextCell(vtkIdType npts, const vtkIdType* pts) VTK_SIZEHINT(pts, npts);
759 
764  vtkIdType InsertNextCell(vtkIdList* pts);
765 
773  vtkIdType InsertNextCell(const std::initializer_list<vtkIdType>& cell)
774  {
775  return this->InsertNextCell(static_cast<vtkIdType>(cell.size()), cell.begin());
776  }
777 
784  vtkIdType InsertNextCell(int npts);
785 
790  void InsertCellPoint(vtkIdType id);
791 
796  void UpdateCellCount(int npts);
797 
812  void ReverseCellAtId(vtkIdType cellId) VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells());
813 
822  void ReplaceCellAtId(vtkIdType cellId, vtkIdList* list);
823  void ReplaceCellAtId(vtkIdType cellId, vtkIdType cellSize, const vtkIdType* cellPoints)
824  VTK_EXPECTS(0 <= cellId && cellId < GetNumberOfCells()) VTK_SIZEHINT(cellPoints, cellSize);
834  void ReplaceCellPointAtId(vtkIdType cellId, vtkIdType cellPointIndex, vtkIdType newPointId);
835 
843  void ReplaceCellAtId(vtkIdType cellId, const std::initializer_list<vtkIdType>& cell)
844  {
845  return this->ReplaceCellAtId(cellId, static_cast<vtkIdType>(cell.size()), cell.begin());
846  }
847 
853 
858 
863 
867  void Append(vtkCellArray* src, vtkIdType pointOffset = 0);
868 
880 
908  void AppendLegacyFormat(const vtkIdType* data, vtkIdType len, vtkIdType ptOffset = 0)
909  VTK_SIZEHINT(data, len);
920  unsigned long GetActualMemorySize() const;
921 
922  // The following code is used to support
923 
924  // The wrappers get understandably confused by some of the template code below
925 #ifndef __VTK_WRAP__
926 
927  // Holds connectivity and offset arrays of the given ArrayType.
928  template <typename ArrayT>
929  struct VisitState
930  {
931  using ArrayType = ArrayT;
932  using ValueType = typename ArrayType::ValueType;
933  using CellRangeType = decltype(vtk::DataArrayValueRange<1>(std::declval<ArrayType>()));
934 
935  // We can't just use is_same here, since binary compatible representations
936  // (e.g. int and long) are distinct types. Instead, ensure that ValueType
937  // is a signed integer the same size as vtkIdType.
938  // If this value is true, ValueType pointers may be safely converted to
939  // vtkIdType pointers via reinterpret cast.
940  static constexpr bool ValueTypeIsSameAsIdType = std::is_integral<ValueType>::value &&
941  std::is_signed<ValueType>::value && (sizeof(ValueType) == sizeof(vtkIdType));
942 
943  ArrayType* GetOffsets() { return this->Offsets; }
944  const ArrayType* GetOffsets() const { return this->Offsets; }
945 
946  ArrayType* GetConnectivity() { return this->Connectivity; }
947  const ArrayType* GetConnectivity() const { return this->Connectivity; }
948 
950 
952 
954 
956 
958 
959  friend class vtkCellArray;
960 
961  protected:
963  {
964  this->Connectivity = vtkSmartPointer<ArrayType>::New();
965  this->Offsets = vtkSmartPointer<ArrayType>::New();
966  this->Offsets->InsertNextValue(0);
968  {
969  this->IsInMemkind = true;
970  }
971  }
972  ~VisitState() = default;
973  void* operator new(size_t nSize)
974  {
975  void* r;
976 #ifdef VTK_USE_MEMKIND
978 #else
979  r = malloc(nSize);
980 #endif
981  return r;
982  }
983  void operator delete(void* p)
984  {
985 #ifdef VTK_USE_MEMKIND
986  VisitState* a = static_cast<VisitState*>(p);
987  if (a->IsInMemkind)
988  {
990  }
991  else
992  {
993  free(p);
994  }
995 #else
996  free(p);
997 #endif
998  }
999 
1002 
1003  private:
1004  VisitState(const VisitState&) = delete;
1005  VisitState& operator=(const VisitState&) = delete;
1006  bool IsInMemkind = false;
1007  };
1008 
1009 private: // Helpers that allow Visit to return a value:
1010  template <typename Functor, typename... Args>
1011  using GetReturnType = decltype(
1012  std::declval<Functor>()(std::declval<VisitState<ArrayType32>&>(), std::declval<Args>()...));
1013 
1014  template <typename Functor, typename... Args>
1015  struct ReturnsVoid : std::is_same<GetReturnType<Functor, Args...>, void>
1016  {
1017  };
1018 
1019 public:
1089  template <typename Functor, typename... Args,
1090  typename = typename std::enable_if<ReturnsVoid<Functor, Args...>::value>::type>
1091  void Visit(Functor&& functor, Args&&... args)
1092  {
1093  if (this->Storage.Is64Bit())
1094  {
1095  // If you get an error on the next line, a call to Visit(functor, Args...)
1096  // is being called with arguments that do not match the functor's call
1097  // signature. See the Visit documentation for details.
1098  functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1099  }
1100  else
1101  {
1102  // If you get an error on the next line, a call to Visit(functor, Args...)
1103  // is being called with arguments that do not match the functor's call
1104  // signature. See the Visit documentation for details.
1105  functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1106  }
1107  }
1108 
1109  template <typename Functor, typename... Args,
1110  typename = typename std::enable_if<ReturnsVoid<Functor, Args...>::value>::type>
1111  void Visit(Functor&& functor, Args&&... args) const
1112  {
1113  if (this->Storage.Is64Bit())
1114  {
1115  // If you get an error on the next line, a call to Visit(functor, Args...)
1116  // is being called with arguments that do not match the functor's call
1117  // signature. See the Visit documentation for details.
1118  functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1119  }
1120  else
1121  {
1122  // If you get an error on the next line, a call to Visit(functor, Args...)
1123  // is being called with arguments that do not match the functor's call
1124  // signature. See the Visit documentation for details.
1125  functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1126  }
1127  }
1128 
1129  template <typename Functor, typename... Args,
1130  typename = typename std::enable_if<!ReturnsVoid<Functor, Args...>::value>::type>
1131  GetReturnType<Functor, Args...> Visit(Functor&& functor, Args&&... args)
1132  {
1133  if (this->Storage.Is64Bit())
1134  {
1135  // If you get an error on the next line, a call to Visit(functor, Args...)
1136  // is being called with arguments that do not match the functor's call
1137  // signature. See the Visit documentation for details.
1138  return functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1139  }
1140  else
1141  {
1142  // If you get an error on the next line, a call to Visit(functor, Args...)
1143  // is being called with arguments that do not match the functor's call
1144  // signature. See the Visit documentation for details.
1145  return functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1146  }
1147  }
1148  template <typename Functor, typename... Args,
1149  typename = typename std::enable_if<!ReturnsVoid<Functor, Args...>::value>::type>
1150  GetReturnType<Functor, Args...> Visit(Functor&& functor, Args&&... args) const
1151  {
1152  if (this->Storage.Is64Bit())
1153  {
1154  // If you get an error on the next line, a call to Visit(functor, Args...)
1155  // is being called with arguments that do not match the functor's call
1156  // signature. See the Visit documentation for details.
1157  return functor(this->Storage.GetArrays64(), std::forward<Args>(args)...);
1158  }
1159  else
1160  {
1161  // If you get an error on the next line, a call to Visit(functor, Args...)
1162  // is being called with arguments that do not match the functor's call
1163  // signature. See the Visit documentation for details.
1164  return functor(this->Storage.GetArrays32(), std::forward<Args>(args)...);
1165  }
1166  }
1167 
1170 #endif // __VTK_WRAP__
1171 
1172  //=================== Begin Legacy Methods ===================================
1173  // These should be deprecated at some point as they are confusing or very slow
1174 
1182 
1194  vtkIdType EstimateSize(vtkIdType numCells, int maxPtsPerCell);
1195 
1205 
1213 
1223  void GetCell(vtkIdType loc, vtkIdType& npts, const vtkIdType*& pts)
1224  VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries()) VTK_SIZEHINT(pts, npts);
1225 
1232  void GetCell(vtkIdType loc, vtkIdList* pts)
1233  VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries());
1234 
1241  vtkIdType GetInsertLocation(int npts);
1242 
1250  vtkIdType GetTraversalLocation();
1251  vtkIdType GetTraversalLocation(vtkIdType npts);
1252  void SetTraversalLocation(vtkIdType loc);
1262  void ReverseCell(vtkIdType loc) VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries());
1263 
1275  void ReplaceCell(vtkIdType loc, int npts, const vtkIdType pts[])
1276  VTK_EXPECTS(0 <= loc && loc < GetNumberOfConnectivityEntries()) VTK_SIZEHINT(pts, npts);
1277 
1292  void SetCells(vtkIdType ncells, vtkIdTypeArray* cells);
1293 
1305 
1306  //=================== End Legacy Methods =====================================
1307 
1308  friend class vtkCellArrayIterator;
1309 
1310 protected:
1312  ~vtkCellArray() override;
1313 
1314  // Encapsulates storage of the internal arrays as a discriminated union
1315  // between 32-bit and 64-bit storage.
1316  struct Storage
1317  {
1318  // Union type that switches 32 and 64 bit array storage
1319  union ArraySwitch {
1320  ArraySwitch() = default; // handled by Storage
1321  ~ArraySwitch() = default; // handle by Storage
1324  };
1325 
1327  {
1328 #ifdef VTK_USE_MEMKIND
1329  this->Arrays =
1331 #else
1332  this->Arrays = new ArraySwitch;
1333 #endif
1334 
1335  // Default to the compile-time setting:
1336 #ifdef VTK_USE_64BIT_IDS
1337 
1338  this->Arrays->Int64 = new VisitState<ArrayType64>;
1339  this->StorageIs64Bit = true;
1340 
1341 #else // VTK_USE_64BIT_IDS
1342 
1343  this->Arrays->Int32 = new VisitState<ArrayType32>;
1344  this->StorageIs64Bit = false;
1345 
1346 #endif // VTK_USE_64BIT_IDS
1347 #ifdef VTK_USE_MEMKIND
1349  {
1350  this->IsInMemkind = true;
1351  }
1352 #else
1353  (void)this->IsInMemkind; // comp warning workaround
1354 #endif
1355  }
1356 
1358  {
1359  if (this->StorageIs64Bit)
1360  {
1361  this->Arrays->Int64->~VisitState();
1362  delete this->Arrays->Int64;
1363  }
1364  else
1365  {
1366  this->Arrays->Int32->~VisitState();
1367  delete this->Arrays->Int32;
1368  }
1369 #ifdef VTK_USE_MEMKIND
1370  if (this->IsInMemkind)
1371  {
1373  }
1374  else
1375  {
1376  free(this->Arrays);
1377  }
1378 #else
1379  delete this->Arrays;
1380 #endif
1381  }
1382 
1383  // Switch the internal arrays to be 32-bit. Any old data is lost. Returns
1384  // true if the storage changes.
1386  {
1387  if (!this->StorageIs64Bit)
1388  {
1389  return false;
1390  }
1391 
1392  this->Arrays->Int64->~VisitState();
1393  delete this->Arrays->Int64;
1394  this->Arrays->Int32 = new VisitState<ArrayType32>;
1395  this->StorageIs64Bit = false;
1396 
1397  return true;
1398  }
1399 
1400  // Switch the internal arrays to be 64-bit. Any old data is lost. Returns
1401  // true if the storage changes.
1403  {
1404  if (this->StorageIs64Bit)
1405  {
1406  return false;
1407  }
1408 
1409  this->Arrays->Int32->~VisitState();
1410  delete this->Arrays->Int32;
1411  this->Arrays->Int64 = new VisitState<ArrayType64>;
1412  this->StorageIs64Bit = true;
1413 
1414  return true;
1415  }
1416 
1417  // Returns true if the storage is currently configured to be 64 bit.
1418  bool Is64Bit() const { return this->StorageIs64Bit; }
1419 
1420  // Get the VisitState for 32-bit arrays
1422  {
1423  assert(!this->StorageIs64Bit);
1424  return *this->Arrays->Int32;
1425  }
1426 
1428  {
1429  assert(!this->StorageIs64Bit);
1430  return *this->Arrays->Int32;
1431  }
1432 
1433  // Get the VisitState for 64-bit arrays
1435  {
1436  assert(this->StorageIs64Bit);
1437  return *this->Arrays->Int64;
1438  }
1439 
1441  {
1442  assert(this->StorageIs64Bit);
1443  return *this->Arrays->Int64;
1444  }
1445 
1446  private:
1447  // Access restricted to ensure proper union construction/destruction thru
1448  // API.
1449  ArraySwitch* Arrays;
1450  bool StorageIs64Bit;
1451  bool IsInMemkind = false;
1452  };
1453 
1456  vtkIdType TraversalCellId{ 0 };
1457 
1459 
1460 private:
1461  vtkCellArray(const vtkCellArray&) = delete;
1462  void operator=(const vtkCellArray&) = delete;
1463 };
1464 
1465 template <typename ArrayT>
1467 {
1468  return this->Offsets->GetNumberOfValues() - 1;
1469 }
1470 
1471 template <typename ArrayT>
1473 {
1474  return static_cast<vtkIdType>(this->Offsets->GetValue(cellId));
1475 }
1476 
1477 template <typename ArrayT>
1479 {
1480  return static_cast<vtkIdType>(this->Offsets->GetValue(cellId + 1));
1481 }
1482 
1483 template <typename ArrayT>
1485 {
1486  return this->GetEndOffset(cellId) - this->GetBeginOffset(cellId);
1487 }
1488 
1489 template <typename ArrayT>
1492 {
1493  return vtk::DataArrayValueRange<1>(
1494  this->GetConnectivity(), this->GetBeginOffset(cellId), this->GetEndOffset(cellId));
1495 }
1496 VTK_ABI_NAMESPACE_END
1497 
1499 {
1500 VTK_ABI_NAMESPACE_BEGIN
1501 
1503 {
1504  // Insert full cell
1505  template <typename CellStateT>
1506  vtkIdType operator()(CellStateT& state, const vtkIdType npts, const vtkIdType pts[])
1507  {
1508  using ValueType = typename CellStateT::ValueType;
1509  auto* conn = state.GetConnectivity();
1510  auto* offsets = state.GetOffsets();
1511 
1512  const vtkIdType cellId = offsets->GetNumberOfValues() - 1;
1513 
1514  offsets->InsertNextValue(static_cast<ValueType>(conn->GetNumberOfValues() + npts));
1515 
1516  for (vtkIdType i = 0; i < npts; ++i)
1517  {
1518  conn->InsertNextValue(static_cast<ValueType>(pts[i]));
1519  }
1520 
1521  return cellId;
1522  }
1523 
1524  // Just update offset table (for incremental API)
1525  template <typename CellStateT>
1526  vtkIdType operator()(CellStateT& state, const vtkIdType npts)
1527  {
1528  using ValueType = typename CellStateT::ValueType;
1529  auto* conn = state.GetConnectivity();
1530  auto* offsets = state.GetOffsets();
1531 
1532  const vtkIdType cellId = offsets->GetNumberOfValues() - 1;
1533 
1534  offsets->InsertNextValue(static_cast<ValueType>(conn->GetNumberOfValues() + npts));
1535 
1536  return cellId;
1537  }
1538 };
1539 
1540 // for incremental API:
1542 {
1543  template <typename CellStateT>
1544  void operator()(CellStateT& state, const vtkIdType npts)
1545  {
1546  using ValueType = typename CellStateT::ValueType;
1547 
1548  auto* offsets = state.GetOffsets();
1549  const ValueType cellBegin = offsets->GetValue(offsets->GetMaxId() - 1);
1550  offsets->SetValue(offsets->GetMaxId(), static_cast<ValueType>(cellBegin + npts));
1551  }
1552 };
1553 
1555 {
1556  template <typename CellStateT>
1557  vtkIdType operator()(CellStateT& state, const vtkIdType cellId)
1558  {
1559  return state.GetCellSize(cellId);
1560  }
1561 };
1562 
1564 {
1565  template <typename CellStateT>
1566  void operator()(CellStateT& state, const vtkIdType cellId, vtkIdList* ids)
1567  {
1568  using ValueType = typename CellStateT::ValueType;
1569 
1570  const vtkIdType beginOffset = state.GetBeginOffset(cellId);
1571  const vtkIdType endOffset = state.GetEndOffset(cellId);
1572  const vtkIdType cellSize = endOffset - beginOffset;
1573  const auto cellConnectivity = state.GetConnectivity()->GetPointer(beginOffset);
1574 
1575  // ValueType differs from vtkIdType, so we have to copy into a temporary buffer:
1576  ids->SetNumberOfIds(cellSize);
1577  vtkIdType* idPtr = ids->GetPointer(0);
1578  for (ValueType i = 0; i < cellSize; ++i)
1579  {
1580  idPtr[i] = static_cast<vtkIdType>(cellConnectivity[i]);
1581  }
1582  }
1583 
1584  // SFINAE helper to check if a VisitState's connectivity array's memory
1585  // can be used as a vtkIdType*.
1586  template <typename CellStateT>
1588  {
1589  private:
1590  using ValueType = typename CellStateT::ValueType;
1591  using ArrayType = typename CellStateT::ArrayType;
1593  static constexpr bool ValueTypeCompat = CellStateT::ValueTypeIsSameAsIdType;
1594  static constexpr bool ArrayTypeCompat = std::is_base_of<AOSArrayType, ArrayType>::value;
1595 
1596  public:
1597  static constexpr bool value = ValueTypeCompat && ArrayTypeCompat;
1598  };
1599 
1600  template <typename CellStateT>
1602  CellStateT& state, const vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints,
1603  vtkIdList* vtkNotUsed(temp))
1604  {
1605  const vtkIdType beginOffset = state.GetBeginOffset(cellId);
1606  const vtkIdType endOffset = state.GetEndOffset(cellId);
1607  cellSize = endOffset - beginOffset;
1608  // This is safe, see CanShareConnPtr helper above.
1609  cellPoints = reinterpret_cast<vtkIdType*>(state.GetConnectivity()->GetPointer(beginOffset));
1610  }
1611 
1612  template <typename CellStateT>
1614  CellStateT& state, const vtkIdType cellId, vtkIdType& cellSize, vtkIdType const*& cellPoints,
1615  vtkIdList* temp)
1616  {
1617  using ValueType = typename CellStateT::ValueType;
1618 
1619  const vtkIdType beginOffset = state.GetBeginOffset(cellId);
1620  const vtkIdType endOffset = state.GetEndOffset(cellId);
1621  cellSize = endOffset - beginOffset;
1622  const ValueType* cellConnectivity = state.GetConnectivity()->GetPointer(beginOffset);
1623 
1624  // ValueType differs from vtkIdType, so we have to copy into a temporary buffer:
1625  temp->SetNumberOfIds(cellSize);
1626  vtkIdType* tempPtr = temp->GetPointer(0);
1627  for (vtkIdType i = 0; i < cellSize; ++i)
1628  {
1629  tempPtr[i] = static_cast<vtkIdType>(cellConnectivity[i]);
1630  }
1631 
1632  cellPoints = temp->GetPointer(0);
1633  }
1634 };
1635 
1637 {
1638  template <typename CellStateT>
1639  void operator()(CellStateT& state)
1640  {
1641  state.GetOffsets()->Reset();
1642  state.GetConnectivity()->Reset();
1643  state.GetOffsets()->InsertNextValue(0);
1644  }
1645 };
1646 
1647 VTK_ABI_NAMESPACE_END
1648 } // end namespace vtkCellArray_detail
1649 
1650 VTK_ABI_NAMESPACE_BEGIN
1651 //----------------------------------------------------------------------------
1653 {
1654  this->TraversalCellId = 0;
1655 }
1656 
1657 //----------------------------------------------------------------------------
1658 inline int vtkCellArray::GetNextCell(vtkIdType& npts, vtkIdType const*& pts) VTK_SIZEHINT(pts, npts)
1659 {
1660  if (this->TraversalCellId < this->GetNumberOfCells())
1661  {
1662  this->GetCellAtId(this->TraversalCellId, npts, pts);
1663  ++this->TraversalCellId;
1664  return 1;
1665  }
1666 
1667  npts = 0;
1668  pts = nullptr;
1669  return 0;
1670 }
1671 
1672 //----------------------------------------------------------------------------
1674 {
1675  if (this->TraversalCellId < this->GetNumberOfCells())
1676  {
1677  this->GetCellAtId(this->TraversalCellId, pts);
1678  ++this->TraversalCellId;
1679  return 1;
1680  }
1681 
1682  pts->Reset();
1683  return 0;
1684 }
1685 //----------------------------------------------------------------------------
1687 {
1688  return this->Visit(vtkCellArray_detail::GetCellSizeImpl{}, cellId);
1689 }
1690 
1691 //----------------------------------------------------------------------------
1692 inline void vtkCellArray::GetCellAtId(vtkIdType cellId, vtkIdType& cellSize,
1693  vtkIdType const*& cellPoints) VTK_SIZEHINT(cellPoints, cellSize)
1694 {
1695  this->Visit(vtkCellArray_detail::GetCellAtIdImpl{}, cellId, cellSize, cellPoints, this->TempCell);
1696 }
1697 
1698 //----------------------------------------------------------------------------
1699 inline void vtkCellArray::GetCellAtId(vtkIdType cellId, vtkIdType& cellSize,
1700  vtkIdType const*& cellPoints, vtkIdList* ptIds) VTK_SIZEHINT(cellPoints, cellSize)
1701 {
1702  this->Visit(vtkCellArray_detail::GetCellAtIdImpl{}, cellId, cellSize, cellPoints, ptIds);
1703 }
1704 
1705 //----------------------------------------------------------------------------
1707 {
1708  this->Visit(vtkCellArray_detail::GetCellAtIdImpl{}, cellId, pts);
1709 }
1710 
1711 //----------------------------------------------------------------------------
1713  VTK_SIZEHINT(pts, npts)
1714 {
1715  return this->Visit(vtkCellArray_detail::InsertNextCellImpl{}, npts, pts);
1716 }
1717 
1718 //----------------------------------------------------------------------------
1720 {
1721  return this->Visit(vtkCellArray_detail::InsertNextCellImpl{}, npts);
1722 }
1723 
1724 //----------------------------------------------------------------------------
1726 {
1727  if (this->Storage.Is64Bit())
1728  {
1729  using ValueType = typename ArrayType64::ValueType;
1730  this->Storage.GetArrays64().Connectivity->InsertNextValue(static_cast<ValueType>(id));
1731  }
1732  else
1733  {
1734  using ValueType = typename ArrayType32::ValueType;
1735  this->Storage.GetArrays32().Connectivity->InsertNextValue(static_cast<ValueType>(id));
1736  }
1737 }
1738 
1739 //----------------------------------------------------------------------------
1740 inline void vtkCellArray::UpdateCellCount(int npts)
1741 {
1743 }
1744 
1745 //----------------------------------------------------------------------------
1747 {
1748  return this->Visit(
1750 }
1751 
1752 //----------------------------------------------------------------------------
1754 {
1755  vtkIdList* pts = cell->GetPointIds();
1756  return this->Visit(
1758 }
1759 
1760 //----------------------------------------------------------------------------
1761 inline void vtkCellArray::Reset()
1762 {
1764 }
1765 
1766 VTK_ABI_NAMESPACE_END
1767 #endif // vtkCellArray.h
Encapsulate traversal logic for vtkCellArray.
object to represent cell connectivity
Definition: vtkCellArray.h:297
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:343
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:461
bool AllocateCopy(vtkCellArray *other)
Pre-allocate memory in internal data structures to match the used size of the input vtkCellArray.
Definition: vtkCellArray.h:382
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:560
ArrayType32 * GetOffsetsArray32()
Return the array used to store cell offsets.
Definition: vtkCellArray.h:636
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:299
void Reset()
Reuse list.
ArrayType64 * GetOffsetsArray64()
Return the array used to store cell offsets.
Definition: vtkCellArray.h:637
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:357
vtkIdType GetNumberOfOffsets() const
Get the number of elements in the offsets array.
Definition: vtkCellArray.h:446
vtkTypeInt64Array ArrayType64
Definition: vtkCellArray.h:300
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:657
ArrayType64 * GetConnectivityArray64()
Return the array used to store the point ids that define the cells' connectivity.
Definition: vtkCellArray.h:658
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:333
vtkDataArray * GetConnectivityArray()
Return the array used to store the point ids that define the cells' connectivity.
Definition: vtkCellArray.h:646
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:430
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:479
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:625
vtkTypeList::Create< ArrayType32, ArrayType64 > StorageArrayList
List of possible array types used for storage.
Definition: vtkCellArray.h:321
vtkIdType InsertNextCell(const std::initializer_list< vtkIdType > &cell)
Overload that allows InsertNextCell({0, 1, 2}) syntax.
Definition: vtkCellArray.h:773
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:552
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:151
vtkIdList * GetPointIds()
Return the list of point ids defining the cell.
Definition: vtkCell.h:246
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:166
list of point or cell ids
Definition: vtkIdList.h:144
vtkIdType * GetPointer(const vtkIdType i)
Get a pointer to a particular data index.
Definition: vtkIdList.h:239
vtkIdType GetNumberOfIds() const noexcept
Return the number of id's in the list.
Definition: vtkIdList.h:170
void Reset()
Reset to an empty state but retain previously allocated memory.
Definition: vtkIdList.h:258
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:120
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:83
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:232
@ type
Definition: vtkX3D.h:528
@ data
Definition: vtkX3D.h:327
VisitState< ArrayType32 > & GetArrays32()
const VisitState< ArrayType64 > & GetArrays64() const
const VisitState< ArrayType32 > & GetArrays32() const
VisitState< ArrayType64 > & GetArrays64()
ArrayType * GetOffsets()
Definition: vtkCellArray.h:943
vtkIdType GetNumberOfCells() const
vtkIdType GetEndOffset(vtkIdType cellId) const
vtkSmartPointer< ArrayType > Offsets
vtkSmartPointer< ArrayType > Connectivity
static constexpr bool ValueTypeIsSameAsIdType
Definition: vtkCellArray.h:940
ArrayType * GetConnectivity()
Definition: vtkCellArray.h:946
decltype(vtk::DataArrayValueRange< 1 >(std::declval< ArrayType >())) CellRangeType
Definition: vtkCellArray.h:933
CellRangeType GetCellRange(vtkIdType cellId)
const ArrayType * GetOffsets() const
Definition: vtkCellArray.h:944
vtkIdType GetBeginOffset(vtkIdType cellId) const
vtkIdType GetCellSize(vtkIdType cellId) const
const ArrayType * GetConnectivity() const
Definition: vtkCellArray.h:947
typename ArrayType::ValueType ValueType
Definition: vtkCellArray.h:932
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:128
VisitState< ArrayType32 > * Int32
VisitState< ArrayType64 > * Int64
int vtkTypeBool
Definition: vtkABI.h:71
STL-compatible iterable ranges that provide access to vtkDataArray elements.
int vtkIdType
Definition: vtkType.h:327
#define VTK_SIZEHINT(...)
#define VTK_EXPECTS(x)
#define VTK_NEWINSTANCE