VecGeom/volumes/MultiUnionStruct.h

0001 #ifndef VECGEOM_VOLUMES_MULTIUNIONSTRUCT_H_
0002 #define VECGEOM_VOLUMES_MULTIUNIONSTRUCT_H_
0003
0004 #include "VecGeom/base/Global.h"
0005 #include "VecGeom/management/HybridManager2.h"
0006 #include "VecGeom/navigation/HybridNavigator2.h"
0007 #include "VecGeom/management/ABBoxManager.h"
0008 #ifndef VECCORE_CUDA
0009 #include <atomic>
0010 #endif
0011
0012 namespace vecgeom {
0013
0014 inline namespace VECGEOM_IMPL_NAMESPACE {
0015
0016 /**
0017  @brief Struct containing placed volume components representing a multiple union solid
0018  @author mihaela.gheata@cern.ch
0019 */
0020
0021 struct MultiUnionStruct {
0022   template <typename U>
0023   using vector_t     = vecgeom::Vector<U>;
0024   using BVHStructure = HybridManager2::HybridBoxAccelerationStructure;
0025
0026   vector_t<VPlacedVolume const *> fVolumes; ///< Component placed volumes
0027   BVHStructure *fNavHelper = nullptr;       ///< Navigation helper using bounding boxes
0028
0029   Vector3D<Precision> fMinExtent;      ///< Minimum extent
0030   Vector3D<Precision> fMaxExtent;      ///< Maximum extent
0031   mutable Precision fCapacity    = -1; ///< Capacity of the multiple union
0032   mutable Precision fSurfaceArea = -1; ///< Surface area of the multiple union
0033 #ifndef VECCORE_CUDA
0034   mutable std::atomic<size_t> fLast; ///< Last located component for opportunistic relocation
0035 #endif
0036   size_t **fNeighbours = nullptr; ///< Array of lists of neigbours
0037   size_t *fNneighbours = nullptr; ///< Number of neighbours for each component
0038
0039   size_t *fBuffer = nullptr; ///< Scratch space for storing neighbours
0040
0041   VECCORE_ATT_HOST_DEVICE
0042   MultiUnionStruct()
0043   {
0044     fMinExtent.Set(kInfLength);
0045     fMaxExtent.Set(-kInfLength);
0046 #ifndef VECCORE_CUDA
0047     fLast.store(0);
0048 #endif
0049   }
0050
0051   VECCORE_ATT_HOST_DEVICE
0052   ~MultiUnionStruct()
0053   {
0054     delete[] fNeighbours;
0055     delete[] fNneighbours;
0056     delete[] fBuffer;
0057   }
0058
0059   VECCORE_ATT_HOST_DEVICE
0060   void AddNode(VPlacedVolume const *volume)
0061   {
0062     using vecCore::math::Max;
0063     using vecCore::math::Min;
0064     Vector3D<Precision> amin, amax;
0065     ABBoxManager::ComputeABBox(volume, &amin, &amax);
0066     fMinExtent.Set(Min(fMinExtent.x(), amin.x()), Min(fMinExtent.y(), amin.y()), Min(fMinExtent.z(), amin.z()));
0067     fMaxExtent.Set(Max(fMaxExtent.x(), amax.x()), Max(fMaxExtent.y(), amax.y()), Max(fMaxExtent.z(), amax.z()));
0068     fVolumes.push_back(volume);
0069   }
0070
0071   VECCORE_ATT_HOST_DEVICE
0072   VECGEOM_FORCE_INLINE
0073   bool ABBoxOverlap(Vector3D<Precision> const &amin1, Vector3D<Precision> const &amax1,
0074                     Vector3D<Precision> const &amin2, Vector3D<Precision> const &amax2)
0075   {
0076     // Check if two aligned boxes overlap
0077     if ((amax1 - amin2).Min() < -kTolerance || (amax2 - amin1).Min() < -kTolerance) return false;
0078     return true;
0079   }
0080
0081   VECCORE_ATT_HOST_DEVICE
0082   void Close()
0083   {
0084     // This method prepares the navigation structure
0085     using Boxes_t           = ABBoxManager::ABBoxContainer_t;
0086     using BoxCorner_t       = ABBoxManager::ABBox_s;
0087     size_t nboxes           = fVolumes.size();
0088     BoxCorner_t *boxcorners = new BoxCorner_t[2 * nboxes];
0089     Vector3D<Precision> amin, amax;
0090     for (size_t i = 0; i < nboxes; ++i)
0091       ABBoxManager::ComputeABBox(fVolumes[i], &boxcorners[2 * i], &boxcorners[2 * i + 1]);
0092     Boxes_t boxes = &boxcorners[0];
0093     fNavHelper    = HybridManager2::Instance().BuildStructure(boxes, nboxes);
0094     // Compute the lists of possibly overlapping neighbours
0095     fNeighbours  = new size_t *[nboxes];
0096     fNneighbours = new size_t[nboxes];
0097     memset(fNneighbours, 0, nboxes * sizeof(size_t));
0098     fBuffer = new size_t[nboxes * nboxes];
0099     for (size_t i = 0; i < nboxes; ++i) {
0100       fNeighbours[i] = fBuffer + i * nboxes;
0101     }
0102     size_t newsize = 0;
0103     for (size_t i = 0; i < nboxes - 1; ++i) {
0104       for (size_t j = i + 1; j < nboxes; ++j) {
0105         if (ABBoxOverlap(boxcorners[2 * i], boxcorners[2 * i + 1], boxcorners[2 * j], boxcorners[2 * j + 1])) {
0106           fNeighbours[i][fNneighbours[i]++] = j;
0107           fNeighbours[j][fNneighbours[j]++] = i;
0108           newsize += 2;
0109         }
0110       }
0111     }
0112     // Compacting buffer of neighbours
0113     size_t *buffer  = new size_t[newsize];
0114     size_t *nextloc = buffer;
0115     for (size_t i = 0; i < nboxes; ++i) {
0116       memcpy(nextloc, fNeighbours[i], fNneighbours[i] * sizeof(size_t));
0117       fNeighbours[i] = nextloc;
0118       nextloc += fNneighbours[i];
0119     }
0120     delete[] fBuffer;
0121     fBuffer = buffer;
0122   }
0123
0124 }; // End struct
0125
0126 } // namespace VECGEOM_IMPL_NAMESPACE
0127 } // namespace vecgeom
0128
0129 #endif /* VECGEOM_VOLUMES_MULTIUNIONSTRUCT_H_ */