diff --git a/CMakeLists.txt b/CMakeLists.txt index 2cb65df8e..c6c6f5075 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -25,6 +25,9 @@ option(MANIFOLD_PYBIND "Build python bindings" ON) option(MANIFOLD_CBIND "Build C (FFI) bindings" OFF) option(TRACY_ENABLE "Use tracy profiling" OFF) option(TRACY_MEMORY_USAGE "Track memory allocation with tracy (expensive)" OFF) +option(BUILD_TEST_CGAL "Build CGAL performance comparisons" OFF) +set(MANIFOLD_PAR "NONE" CACHE STRING "Parallel backend, either \"TBB\" or \"NONE\"") +set(MANIFOLD_FLAGS "" CACHE STRING "Manifold compiler flags") if(TRACY_ENABLE) option(CMAKE_BUILD_TYPE "Build type" RelWithDebInfo) @@ -38,7 +41,6 @@ else() option(CMAKE_BUILD_TYPE "Build type" Release) endif() -set(MANIFOLD_PAR "NONE" CACHE STRING "Parallel backend, either \"TBB\" or \"NONE\"") set(MANIFOLD_FLAGS ${MANIFOLD_FLAGS} -O3) if(EMSCRIPTEN) @@ -49,16 +51,15 @@ if(EMSCRIPTEN) endif() if(MANIFOLD_PYBIND) -find_package(Python COMPONENTS Interpreter Development.Module REQUIRED) + set(CMAKE_POSITION_INDEPENDENT_CODE ON) + find_package(Python COMPONENTS Interpreter Development.Module REQUIRED) endif() -option(BUILD_TEST_CGAL "Build CGAL performance comparisons" OFF) -option(BUILD_SHARED_LIBS "Build dynamic/shared libraries" OFF) set(PYBIND11_DIR ${PROJECT_SOURCE_DIR}/bindings/python/third_party/pybind11) set(THRUST_INC_DIR ${PROJECT_SOURCE_DIR}/src/third_party/thrust) -if(CMAKE_EXPORT_COMPILE_COMMANDS AND !EMSCRIPTEN) +if(CMAKE_EXPORT_COMPILE_COMMANDS AND NOT EMSCRIPTEN) # for nixos set(CMAKE_CXX_STANDARD_INCLUDE_DIRECTORIES ${CMAKE_CXX_IMPLICIT_INCLUDE_DIRECTORIES}) @@ -83,6 +84,8 @@ if(CODE_COVERAGE AND NOT MSVC) add_link_options("-coverage") endif() +include(GNUInstallDirs) + add_subdirectory(src) add_subdirectory(bindings) diff --git a/bindings/c/CMakeLists.txt b/bindings/c/CMakeLists.txt index 2899afb28..9fbd9de66 100644 --- a/bindings/c/CMakeLists.txt +++ b/bindings/c/CMakeLists.txt @@ -20,3 +20,14 @@ target_link_libraries( target_include_directories(${PROJECT_NAME} PUBLIC ${PROJECT_SOURCE_DIR}/include) target_compile_options(${PROJECT_NAME} PRIVATE ${MANIFOLD_FLAGS}) target_compile_features(${PROJECT_NAME} PRIVATE cxx_std_17) + +install( + TARGETS ${PROJECT_NAME} + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT bindings +) +install( + FILES include/manifoldc.h include/types.h + DESTINATION ${CMAKE_INSTALL_FULL_INCLUDEDIR}/manifold + COMPONENT devel +) diff --git a/bindings/c/box.cpp b/bindings/c/box.cpp index c6aa5a098..d440ecb1f 100644 --- a/bindings/c/box.cpp +++ b/bindings/c/box.cpp @@ -1,7 +1,20 @@ -#include -#include -#include - +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "conv.h" +#include "manifold.h" +#include "public.h" #include "types.h" ManifoldBox *manifold_box(void *mem, float x1, float y1, float z1, float x2, diff --git a/bindings/c/conv.cpp b/bindings/c/conv.cpp index 84915a694..4f074639e 100644 --- a/bindings/c/conv.cpp +++ b/bindings/c/conv.cpp @@ -1,12 +1,25 @@ -#include -#include -#include +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "conv.h" #include #include "cross_section.h" -#include "include/types.h" +#include "manifold.h" #include "public.h" +#include "sdf.h" #include "types.h" ManifoldManifold *to_c(manifold::Manifold *m) { diff --git a/bindings/c/cross.cpp b/bindings/c/cross.cpp index 152f7778b..2051d44b7 100644 --- a/bindings/c/cross.cpp +++ b/bindings/c/cross.cpp @@ -1,10 +1,23 @@ -#include -#include -#include +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. #include +#include "conv.h" #include "cross_section.h" +#include "manifold.h" +#include "public.h" #include "types.h" ManifoldCrossSection *manifold_cross_section_empty(void *mem) { diff --git a/bindings/c/include/conv.h b/bindings/c/include/conv.h index a769395bc..fee7653a1 100644 --- a/bindings/c/include/conv.h +++ b/bindings/c/include/conv.h @@ -1,13 +1,28 @@ -#pragma once -#include -#include -#include -#include +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. +#pragma once +#include #include #include "cross_section.h" +#include "manifold.h" +#include "public.h" +#include "sdf.h" +#include "types.h" +using namespace manifold; using ManifoldVec = std::vector; using CrossSectionVec = std::vector; diff --git a/bindings/c/include/manifoldc.h b/bindings/c/include/manifoldc.h index 97b691bfa..29fd0b637 100644 --- a/bindings/c/include/manifoldc.h +++ b/bindings/c/include/manifoldc.h @@ -1,6 +1,21 @@ +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + #include #include -#include + +#include "types.h" #ifdef __cplusplus #include diff --git a/bindings/c/include/types.h b/bindings/c/include/types.h index 2c10be5a4..fa9a44e89 100644 --- a/bindings/c/include/types.h +++ b/bindings/c/include/types.h @@ -1,3 +1,17 @@ +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + #pragma once #include diff --git a/bindings/c/manifoldc.cpp b/bindings/c/manifoldc.cpp index f778f0b58..44843b80b 100644 --- a/bindings/c/manifoldc.cpp +++ b/bindings/c/manifoldc.cpp @@ -1,3 +1,17 @@ +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + #include #include #include @@ -29,8 +43,7 @@ ManifoldMeshGL *level_set(void *mem, float (*sdf)(float, float, float), std::function fun = [sdf](glm::vec3 v) { return (sdf(v.x, v.y, v.z)); }; - auto pol = seq ? std::make_optional(ExecutionPolicy::Seq) : std::nullopt; - auto mesh = LevelSet(fun, *from_c(bounds), edge_length, level, pol); + auto mesh = LevelSet(fun, *from_c(bounds), edge_length, level, !seq); return to_c(new (mem) MeshGL(mesh)); } ManifoldMeshGL *level_set_context( @@ -43,8 +56,7 @@ ManifoldMeshGL *level_set_context( std::function fun = [sdf](glm::vec3 v) { return (sdf(v.x, v.y, v.z)); }; - auto pol = seq ? std::make_optional(ExecutionPolicy::Seq) : std::nullopt; - auto mesh = LevelSet(fun, *from_c(bounds), edge_length, level, pol); + auto mesh = LevelSet(fun, *from_c(bounds), edge_length, level, !seq); return to_c(new (mem) MeshGL(mesh)); } } // namespace diff --git a/bindings/c/meshio.cpp b/bindings/c/meshio.cpp index 2bb0c5b39..cb873d166 100644 --- a/bindings/c/meshio.cpp +++ b/bindings/c/meshio.cpp @@ -1,8 +1,21 @@ -#include "meshIO.h" +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. -#include +#include "meshIO.h" -#include "include/types.h" +#include "conv.h" +#include "types.h" // C <-> C++ conversions diff --git a/bindings/c/rect.cpp b/bindings/c/rect.cpp index 1feb02734..4dc717164 100644 --- a/bindings/c/rect.cpp +++ b/bindings/c/rect.cpp @@ -1,8 +1,21 @@ -#include -#include -#include - +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "conv.h" #include "cross_section.h" +#include "manifold.h" +#include "public.h" #include "types.h" ManifoldRect *manifold_rect(void *mem, float x1, float y1, float x2, float y2) { diff --git a/bindings/python/CMakeLists.txt b/bindings/python/CMakeLists.txt index c88e32014..53a565688 100644 --- a/bindings/python/CMakeLists.txt +++ b/bindings/python/CMakeLists.txt @@ -22,3 +22,9 @@ target_compile_features(manifold3d PUBLIC cxx_std_17) target_include_directories(manifold3d PRIVATE ${PYBIND11_DIR}/include ) +set_target_properties(manifold3d PROPERTIES OUTPUT_NAME "manifold3d") +install( + TARGETS manifold3d + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT bindings +) diff --git a/bindings/wasm/bindings.cpp b/bindings/wasm/bindings.cpp index 925e71338..7357bb3ca 100644 --- a/bindings/wasm/bindings.cpp +++ b/bindings/wasm/bindings.cpp @@ -14,14 +14,14 @@ #include #include -#include -#include -#include #include #include "cross_section.h" #include "helpers.cpp" +#include "manifold.h" +#include "polygon.h" +#include "sdf.h" using namespace emscripten; using namespace manifold; diff --git a/bindings/wasm/helpers.cpp b/bindings/wasm/helpers.cpp index 30be83c4e..c8da6b956 100644 --- a/bindings/wasm/helpers.cpp +++ b/bindings/wasm/helpers.cpp @@ -1,12 +1,12 @@ #include #include -#include -#include -#include #include #include "cross_section.h" +#include "manifold.h" +#include "polygon.h" +#include "sdf.h" using namespace emscripten; using namespace manifold; diff --git a/extras/CMakeLists.txt b/extras/CMakeLists.txt index 405bd3ebc..e8914ac4d 100644 --- a/extras/CMakeLists.txt +++ b/extras/CMakeLists.txt @@ -42,7 +42,6 @@ if(BUILD_TEST_CGAL) add_executable(perfTestCGAL perf_test_cgal.cpp) find_package(CGAL REQUIRED COMPONENTS Core) target_compile_definitions(perfTestCGAL PRIVATE CGAL_USE_GMPXX) - # target_compile_definitions(perfTestCGAL PRIVATE CGAL_DEBUG) target_link_libraries(perfTestCGAL manifold CGAL::CGAL CGAL::CGAL_Core) target_compile_options(perfTestCGAL PRIVATE ${MANIFOLD_FLAGS}) target_compile_features(perfTestCGAL PUBLIC cxx_std_17) diff --git a/flake.nix b/flake.nix index 66309c428..30b501c17 100644 --- a/flake.nix +++ b/flake.nix @@ -26,12 +26,15 @@ src = self; nativeBuildInputs = (with pkgs; [ cmake + ninja (python39.withPackages (ps: with ps; [ trimesh ])) gtest ]) ++ build-tools; cmakeFlags = [ "-DMANIFOLD_PYBIND=ON" + "-DMANIFOLD_CBIND=ON" + "-DBUILD_SHARED_LIBS=ON" "-DMANIFOLD_PAR=${pkgs.lib.strings.toUpper parallel-backend}" ]; checkPhase = '' @@ -41,12 +44,6 @@ PYTHONPATH=$PYTHONPATH:$(pwd)/build/bindings/python python3 bindings/python/examples/run_all.py cd build ''; - installPhase = '' - mkdir -p $out - cp src/manifold/libmanifold.a $out/ - cp extras/perfTest $out - cp bindings/python/manifold3d* $out - ''; }; parallelBackends = [ { parallel-backend = "none"; } diff --git a/meshIO/CMakeLists.txt b/meshIO/CMakeLists.txt index c70b0d1fe..ef8ae03a1 100644 --- a/meshIO/CMakeLists.txt +++ b/meshIO/CMakeLists.txt @@ -29,3 +29,14 @@ target_link_libraries(${PROJECT_NAME} target_compile_options(${PROJECT_NAME} PRIVATE ${MANIFOLD_FLAGS}) target_compile_features(${PROJECT_NAME} PUBLIC cxx_std_17) + +install( + TARGETS ${PROJECT_NAME} + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT devel +) +install( + DIRECTORY include/ + DESTINATION ${CMAKE_INSTALL_FULL_INCLUDEDIR}/manifold + COMPONENT devel +) diff --git a/samples/CMakeLists.txt b/samples/CMakeLists.txt index 85f78604d..2cbafc423 100644 --- a/samples/CMakeLists.txt +++ b/samples/CMakeLists.txt @@ -14,32 +14,17 @@ project(samples) -add_library(samples src/bracelet.cpp src/knot.cpp src/menger_sponge.cpp src/rounded_frame.cpp src/scallop.cpp src/tet_puzzle.cpp) +add_library(samples src/bracelet.cpp src/knot.cpp src/menger_sponge.cpp + src/rounded_frame.cpp src/scallop.cpp src/tet_puzzle.cpp + src/gyroid_module.cpp) target_include_directories(samples PUBLIC ${PROJECT_SOURCE_DIR}/include ) target_link_libraries(samples - PUBLIC manifold + PUBLIC manifold sdf ) target_compile_options(samples PRIVATE ${MANIFOLD_FLAGS}) target_compile_features(samples PUBLIC cxx_std_17) - -project(samplesGPU) - -add_library(samplesGPU src/gyroid_module.cpp) - -target_include_directories(samplesGPU - PUBLIC ${PROJECT_SOURCE_DIR}/include -) - -target_link_libraries(samplesGPU - PUBLIC manifold sdf -) - -target_compile_options(samplesGPU PRIVATE ${MANIFOLD_FLAGS}) -target_compile_features(samplesGPU - PUBLIC cxx_std_17 -) diff --git a/samples/src/gyroid_module.cpp b/samples/src/gyroid_module.cpp index cf8b6ed89..0f6bee74c 100644 --- a/samples/src/gyroid_module.cpp +++ b/samples/src/gyroid_module.cpp @@ -12,13 +12,15 @@ // See the License for the specific language governing permissions and // limitations under the License. +#include "manifold.h" #include "samples.h" #include "sdf.h" namespace { +using namespace manifold; struct Gyroid { - __host__ __device__ float operator()(glm::vec3 p) const { + float operator()(glm::vec3 p) const { p -= glm::pi() / 4; return cos(p.x) * sin(p.y) + cos(p.y) * sin(p.z) + cos(p.z) * sin(p.x); } @@ -54,4 +56,4 @@ Manifold GyroidModule(float size, int n) { return result.Rotate(-45, 0, 90).Translate({0, 0, size / glm::sqrt(2.0f)}); } -} // namespace manifold \ No newline at end of file +} // namespace manifold diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index e3af970d6..49ec7a3da 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -12,7 +12,8 @@ # See the License for the specific language governing permissions and # limitations under the License. -add_subdirectory(third_party) +# do not add third_party installations +add_subdirectory(third_party EXCLUDE_FROM_ALL) add_subdirectory(utilities) add_subdirectory(collider) add_subdirectory(cross_section) diff --git a/src/collider/CMakeLists.txt b/src/collider/CMakeLists.txt index e936710ed..f72658fb3 100644 --- a/src/collider/CMakeLists.txt +++ b/src/collider/CMakeLists.txt @@ -23,3 +23,9 @@ target_link_libraries(${PROJECT_NAME} PUBLIC utilities) target_compile_features(${PROJECT_NAME} PUBLIC cxx_std_17 ) + +install( + TARGETS ${PROJECT_NAME} + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT runtime +) diff --git a/src/collider/include/collider.h b/src/collider/include/collider.h index 3decdc6bf..7d81821fe 100644 --- a/src/collider/include/collider.h +++ b/src/collider/include/collider.h @@ -15,7 +15,7 @@ #pragma once #include "public.h" #include "sparse.h" -#include "vec_dh.h" +#include "vec.h" namespace manifold { @@ -23,18 +23,19 @@ namespace manifold { class Collider { public: Collider() {} - Collider(const VecDH& leafBB, const VecDH& leafMorton); + Collider(const VecView& leafBB, + const VecView& leafMorton); bool Transform(glm::mat4x3); - void UpdateBoxes(const VecDH& leafBB); + void UpdateBoxes(const VecView& leafBB); template - SparseIndices Collisions(const VecDH& queriesIn) const; + SparseIndices Collisions(const VecView& queriesIn) const; private: - VecDH nodeBBox_; - VecDH nodeParent_; + Vec nodeBBox_; + Vec nodeParent_; // even nodes are leaves, odd nodes are internal, root is 1 - VecDH> internalChildren_; + Vec> internalChildren_; int NumInternal() const { return internalChildren_.size(); }; int NumLeaves() const { return NumInternal() + 1; }; diff --git a/src/collider/src/collider.cpp b/src/collider/src/collider.cpp index 9b0ce5da1..b90afd3f0 100644 --- a/src/collider/src/collider.cpp +++ b/src/collider/src/collider.cpp @@ -68,9 +68,9 @@ int Node2Leaf(int node) { return node / 2; } int Leaf2Node(int leaf) { return leaf * 2; } struct CreateRadixTree { - int* nodeParent_; - thrust::pair* internalChildren_; - const VecDc leafMorton_; + VecView nodeParent_; + VecView> internalChildren_; + const VecView leafMorton_; int PrefixLength(uint32_t a, uint32_t b) const { // count-leading-zeros is used to find the number of identical highest-order @@ -154,9 +154,9 @@ struct CreateRadixTree { template struct FindCollisions { - const Box* nodeBBox_; - const thrust::pair* internalChildren_; - const Recorder& recorder; + VecView nodeBBox_; + VecView> internalChildren_; + Recorder recorder; int RecordCollision(int node, thrust::tuple& query) { const T& queryObj = thrust::get<0>(query); @@ -205,11 +205,11 @@ struct FindCollisions { }; struct CountCollisions { - int* counts; - char* empty; - void record(int queryIdx, int _leafIdx) const { counts[queryIdx]++; } - bool earlyexit(int _queryIdx) const { return false; } - void end(int queryIdx) const { + VecView counts; + VecView empty; + void record(int queryIdx, int _leafIdx) { counts[queryIdx]++; } + bool earlyexit(int _queryIdx) { return false; } + void end(int queryIdx) { if (counts[queryIdx] == 0) empty[queryIdx] = 1; } }; @@ -230,9 +230,9 @@ struct SeqCollisionRecorder { template struct ParCollisionRecorder { SparseIndices& queryTri; - int* counts; - char* empty; - void record(int queryIdx, int leafIdx) const { + VecView counts; + VecView empty; + void record(int queryIdx, int leafIdx) { int pos = counts[queryIdx]++; if (inverted) queryTri.Set(pos, leafIdx, queryIdx); @@ -244,10 +244,10 @@ struct ParCollisionRecorder { }; struct BuildInternalBoxes { - Box* nodeBBox_; - int* counter_; - const int* nodeParent_; - const thrust::pair* internalChildren_; + VecView nodeBBox_; + VecView counter_; + const VecView nodeParent_; + const VecView> internalChildren_; void operator()(int leaf) { int node = Leaf2Node(leaf); @@ -274,8 +274,8 @@ namespace manifold { * bounding boxes and corresponding Morton codes. It is assumed these vectors * are already sorted by increasing Morton code. */ -Collider::Collider(const VecDH& leafBB, - const VecDH& leafMorton) { +Collider::Collider(const VecView& leafBB, + const VecView& leafMorton) { ASSERT(leafBB.size() == leafMorton.size(), userErr, "vectors must be the same length"); int num_nodes = 2 * leafBB.size() - 1; @@ -285,8 +285,7 @@ Collider::Collider(const VecDH& leafBB, internalChildren_.resize(leafBB.size() - 1, thrust::make_pair(-1, -1)); // organize tree for_each_n(autoPolicy(NumInternal()), countAt(0), NumInternal(), - CreateRadixTree( - {nodeParent_.ptrD(), internalChildren_.ptrD(), leafMorton})); + CreateRadixTree({nodeParent_, internalChildren_, leafMorton})); UpdateBoxes(leafBB); } @@ -299,7 +298,7 @@ Collider::Collider(const VecDH& leafBB, * then not report any collisions between an index and itself. */ template -SparseIndices Collider::Collisions(const VecDH& queriesIn) const { +SparseIndices Collider::Collisions(const VecView& queriesIn) const { // note that the length is 1 larger than the number of queries so the last // element can store the sum when using exclusive scan if (queriesIn.size() < kSequentialThreshold) { @@ -307,19 +306,17 @@ SparseIndices Collider::Collisions(const VecDH& queriesIn) const { for_each_n(ExecutionPolicy::Seq, zip(queriesIn.cbegin(), countAt(0)), queriesIn.size(), FindCollisions>{ - nodeBBox_.ptrD(), internalChildren_.ptrD(), {queryTri}}); + nodeBBox_, internalChildren_, {queryTri}}); return queryTri; } else { // compute the number of collisions to determine the size for allocation and // offset, this avoids the need for atomic - VecDH counts(queriesIn.size() + 1, 0); - VecDH empty(queriesIn.size(), 0); + Vec counts(queriesIn.size() + 1, 0); + Vec empty(queriesIn.size(), 0); for_each_n(ExecutionPolicy::Par, zip(queriesIn.cbegin(), countAt(0)), queriesIn.size(), FindCollisions{ - nodeBBox_.ptrD(), - internalChildren_.ptrD(), - {counts.ptrD(), empty.ptrD()}}); + nodeBBox_, internalChildren_, {counts, empty}}); // compute start index for each query and total count exclusive_scan(ExecutionPolicy::Par, counts.begin(), counts.end(), counts.begin(), 0, std::plus()); @@ -329,9 +326,7 @@ SparseIndices Collider::Collisions(const VecDH& queriesIn) const { for_each_n(ExecutionPolicy::Par, zip(queriesIn.cbegin(), countAt(0)), queriesIn.size(), FindCollisions>{ - nodeBBox_.ptrD(), - internalChildren_.ptrD(), - {queryTri, counts.ptrD(), empty.ptrD()}}); + nodeBBox_, internalChildren_, {queryTri, counts, empty}}); return queryTri; } } @@ -340,20 +335,19 @@ SparseIndices Collider::Collisions(const VecDH& queriesIn) const { * Recalculate the collider's internal bounding boxes without changing the * hierarchy. */ -void Collider::UpdateBoxes(const VecDH& leafBB) { +void Collider::UpdateBoxes(const VecView& leafBB) { ASSERT(leafBB.size() == NumLeaves(), userErr, "must have the same number of updated boxes as original"); // copy in leaf node Boxes - strided_range::Iter> leaves(nodeBBox_.begin(), nodeBBox_.end(), 2); + strided_range::Iter> leaves(nodeBBox_.begin(), nodeBBox_.end(), 2); auto policy = autoPolicy(NumInternal()); copy(policy, leafBB.cbegin(), leafBB.cend(), leaves.begin()); // create global counters - VecDH counter(NumInternal(), 0); + Vec counter(NumInternal(), 0); // kernel over leaves to save internal Boxes for_each_n( policy, countAt(0), NumLeaves(), - BuildInternalBoxes({nodeBBox_.ptrD(), counter.ptrD(), nodeParent_.cptrD(), - internalChildren_.cptrD()})); + BuildInternalBoxes({nodeBBox_, counter, nodeParent_, internalChildren_})); } /** @@ -377,21 +371,21 @@ bool Collider::Transform(glm::mat4x3 transform) { } template SparseIndices Collider::Collisions( - const VecDH&) const; + const VecView&) const; template SparseIndices Collider::Collisions( - const VecDH&) const; + const VecView&) const; template SparseIndices Collider::Collisions( - const VecDH&) const; + const VecView&) const; template SparseIndices Collider::Collisions( - const VecDH&) const; + const VecView&) const; template SparseIndices Collider::Collisions( - const VecDH&) const; + const VecView&) const; template SparseIndices Collider::Collisions( - const VecDH&) const; + const VecView&) const; } // namespace manifold diff --git a/src/cross_section/CMakeLists.txt b/src/cross_section/CMakeLists.txt index 4f9c28d00..f78194ba6 100644 --- a/src/cross_section/CMakeLists.txt +++ b/src/cross_section/CMakeLists.txt @@ -26,3 +26,15 @@ target_link_libraries( ${PROJECT_NAME} target_compile_options(${PROJECT_NAME} PRIVATE ${MANIFOLD_FLAGS}) target_compile_features(${PROJECT_NAME} PUBLIC cxx_std_17) + +install( + TARGETS ${PROJECT_NAME} + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT runtime +) +install( + DIRECTORY include/ + DESTINATION ${CMAKE_INSTALL_FULL_INCLUDEDIR}/manifold + COMPONENT devel +) + diff --git a/src/cross_section/include/cross_section.h b/src/cross_section/include/cross_section.h index 4cfdba95a..b8d241555 100644 --- a/src/cross_section/include/cross_section.h +++ b/src/cross_section/include/cross_section.h @@ -14,19 +14,14 @@ #pragma once -#include - +#include #include #include -#include "clipper2/clipper.core.h" -#include "clipper2/clipper.offset.h" #include "glm/ext/matrix_float3x2.hpp" #include "glm/ext/vector_float2.hpp" #include "public.h" -namespace C2 = Clipper2Lib; - namespace manifold { class Rect; @@ -35,6 +30,8 @@ class Rect; * @{ */ +struct PathImpl; + /** * Two-dimensional cross sections guaranteed to be without self-intersections, * or overlaps between polygons (from construction onwards). This class makes @@ -169,10 +166,10 @@ class CrossSection { ///@} private: - mutable std::shared_ptr paths_; + mutable std::shared_ptr paths_; mutable glm::mat3x2 transform_ = glm::mat3x2(1.0f); - CrossSection(C2::PathsD paths); - C2::PathsD GetPaths() const; + CrossSection(std::shared_ptr paths); + std::shared_ptr GetPaths() const; }; /** @} */ diff --git a/src/cross_section/src/cross_section.cpp b/src/cross_section/src/cross_section.cpp index 02880b4f7..8d41db74c 100644 --- a/src/cross_section/src/cross_section.cpp +++ b/src/cross_section/src/cross_section.cpp @@ -14,8 +14,22 @@ #include "cross_section.h" +#include "clipper2/clipper.core.h" +#include "clipper2/clipper.h" +#include "clipper2/clipper.offset.h" + +namespace C2 = Clipper2Lib; + using namespace manifold; +namespace manifold { +struct PathImpl { + PathImpl(const C2::PathsD paths_) : paths_(paths_) {} + operator const C2::PathsD&() const { return paths_; } + const C2::PathsD paths_; +}; +} // namespace manifold + namespace { const int precision_ = 8; @@ -96,8 +110,8 @@ C2::PathsD transform(const C2::PathsD ps, const glm::mat3x2 m) { return transformed; } -std::shared_ptr shared_paths(const C2::PathsD& ps) { - return std::make_shared(ps); +std::shared_ptr shared_paths(const C2::PathsD& ps) { + return std::make_shared(ps); } // forward declaration for mutual recursion @@ -191,7 +205,9 @@ namespace manifold { /** * The default constructor is an empty cross-section (containing no contours). */ -CrossSection::CrossSection() { paths_ = shared_paths(C2::PathsD()); } +CrossSection::CrossSection() { + paths_ = std::make_shared(C2::PathsD()); +} CrossSection::~CrossSection() = default; CrossSection::CrossSection(CrossSection&&) noexcept = default; @@ -218,7 +234,7 @@ CrossSection& CrossSection::operator=(const CrossSection& other) { }; // Private, skips unioning. -CrossSection::CrossSection(C2::PathsD ps) { paths_ = shared_paths(ps); } +CrossSection::CrossSection(std::shared_ptr ps) { paths_ = ps; } /** * Create a 2d cross-section from a single contour. A boolean union operation @@ -271,13 +287,13 @@ CrossSection::CrossSection(const Rect& rect) { // Private // All access to paths_ should be done through the GetPaths() method, which // applies the accumulated transform_ -C2::PathsD CrossSection::GetPaths() const { +std::shared_ptr CrossSection::GetPaths() const { if (transform_ == glm::mat3x2(1.0f)) { - return *paths_; + return paths_; } - paths_ = shared_paths(transform(*paths_, transform_)); + paths_ = shared_paths(transform(paths_->paths_, transform_)); transform_ = glm::mat3x2(1.0f); - return *paths_; + return paths_; } /** @@ -309,7 +325,7 @@ CrossSection CrossSection::Square(const glm::vec2 size, bool center) { p[2] = C2::PointD(x, y); p[3] = C2::PointD(0.0, y); } - return CrossSection(C2::PathsD{p}); + return CrossSection(shared_paths(C2::PathsD{p})); } /** @@ -330,7 +346,7 @@ CrossSection CrossSection::Circle(float radius, int circularSegments) { for (int i = 0; i < n; ++i) { circle[i] = C2::PointD(radius * cosd(dPhi * i), radius * sind(dPhi * i)); } - return CrossSection(C2::PathsD{circle}); + return CrossSection(shared_paths(C2::PathsD{circle})); } /** @@ -339,9 +355,9 @@ CrossSection CrossSection::Circle(float radius, int circularSegments) { CrossSection CrossSection::Boolean(const CrossSection& second, OpType op) const { auto ct = cliptype_of_op(op); - auto res = C2::BooleanOp(ct, C2::FillRule::Positive, GetPaths(), - second.GetPaths(), precision_); - return CrossSection(res); + auto res = C2::BooleanOp(ct, C2::FillRule::Positive, GetPaths()->paths_, + second.GetPaths()->paths_, precision_); + return CrossSection(shared_paths(res)); } /** @@ -358,19 +374,19 @@ CrossSection CrossSection::BatchBoolean( auto subjs = crossSections[0].GetPaths(); int n_clips = 0; for (int i = 1; i < crossSections.size(); ++i) { - n_clips += crossSections[i].GetPaths().size(); + n_clips += crossSections[i].GetPaths()->paths_.size(); } auto clips = C2::PathsD(); clips.reserve(n_clips); for (int i = 1; i < crossSections.size(); ++i) { auto ps = crossSections[i].GetPaths(); - clips.insert(clips.end(), ps.begin(), ps.end()); + clips.insert(clips.end(), ps->paths_.begin(), ps->paths_.end()); } auto ct = cliptype_of_op(op); - auto res = - C2::BooleanOp(ct, C2::FillRule::Positive, subjs, clips, precision_); - return CrossSection(res); + auto res = C2::BooleanOp(ct, C2::FillRule::Positive, subjs->paths_, clips, + precision_); + return CrossSection(shared_paths(res)); } /** @@ -441,7 +457,7 @@ std::vector CrossSection::Decompose() const { } C2::PolyTreeD tree; - C2::BooleanOp(C2::ClipType::Union, C2::FillRule::Positive, GetPaths(), + C2::BooleanOp(C2::ClipType::Union, C2::FillRule::Positive, GetPaths()->paths_, C2::PathsD(), tree, precision_); auto polys = std::vector(); @@ -451,7 +467,7 @@ std::vector CrossSection::Decompose() const { auto comps = std::vector(n_polys); // reverse the stack while wrapping for (int i = 0; i < n_polys; ++i) { - comps[n_polys - i - 1] = CrossSection(polys[i]); + comps[n_polys - i - 1] = CrossSection(shared_paths(polys[i])); } return comps; @@ -541,8 +557,8 @@ CrossSection CrossSection::Warp( std::function warpFunc) const { auto paths = GetPaths(); auto warped = C2::PathsD(); - warped.reserve(paths.size()); - for (auto path : paths) { + warped.reserve(paths->paths_.size()); + for (auto path : paths->paths_) { auto sz = path.size(); auto s = C2::PathD(sz); for (int i = 0; i < sz; ++i) { @@ -552,7 +568,8 @@ CrossSection CrossSection::Warp( } warped.push_back(s); } - return CrossSection(C2::Union(warped, C2::FillRule::Positive, precision_)); + return CrossSection( + shared_paths(C2::Union(warped, C2::FillRule::Positive, precision_))); } /** @@ -568,8 +585,8 @@ CrossSection CrossSection::Warp( * offseting operations are to be performed, which would compound the issue. */ CrossSection CrossSection::Simplify(double epsilon) const { - auto ps = SimplifyPaths(GetPaths(), epsilon, false); - return CrossSection(ps); + auto ps = SimplifyPaths(GetPaths()->paths_, epsilon, false); + return CrossSection(shared_paths(ps)); } /** @@ -607,9 +624,9 @@ CrossSection CrossSection::Offset(double delta, JoinType jointype, arc_tol = (std::cos(Clipper2Lib::PI / n) - 1) * -scaled_delta; } auto ps = - C2::InflatePaths(GetPaths(), delta, jt(jointype), C2::EndType::Polygon, - miter_limit, precision_, arc_tol); - return CrossSection(ps); + C2::InflatePaths(GetPaths()->paths_, delta, jt(jointype), + C2::EndType::Polygon, miter_limit, precision_, arc_tol); + return CrossSection(shared_paths(ps)); } /** @@ -625,14 +642,14 @@ CrossSection CrossSection::Hull( SimplePolygon pts; pts.reserve(n); for (auto cs : crossSections) { - auto paths = cs.GetPaths(); + auto paths = cs.GetPaths()->paths_; for (auto path : paths) { for (auto p : path) { pts.push_back(v2_of_pd(p)); } } } - return CrossSection(C2::PathsD{HullImpl(pts)}); + return CrossSection(shared_paths(C2::PathsD{HullImpl(pts)})); } /** @@ -650,7 +667,7 @@ CrossSection CrossSection::Hull() const { * hull. */ CrossSection CrossSection::Hull(SimplePolygon pts) { - return CrossSection(C2::PathsD{HullImpl(pts)}); + return CrossSection(shared_paths(C2::PathsD{HullImpl(pts)})); } /** @@ -674,14 +691,14 @@ CrossSection CrossSection::Hull(const Polygons polys) { * Return the total area covered by complex polygons making up the * CrossSection. */ -double CrossSection::Area() const { return C2::Area(GetPaths()); } +double CrossSection::Area() const { return C2::Area(GetPaths()->paths_); } /** * Return the number of vertices in the CrossSection. */ int CrossSection::NumVert() const { int n = 0; - auto paths = GetPaths(); + auto paths = GetPaths()->paths_; for (auto p : paths) { n += p.size(); } @@ -692,19 +709,19 @@ int CrossSection::NumVert() const { * Return the number of contours (both outer and inner paths) in the * CrossSection. */ -int CrossSection::NumContour() const { return GetPaths().size(); } +int CrossSection::NumContour() const { return GetPaths()->paths_.size(); } /** * Does the CrossSection contain any contours? */ -bool CrossSection::IsEmpty() const { return GetPaths().empty(); } +bool CrossSection::IsEmpty() const { return GetPaths()->paths_.empty(); } /** * Returns the axis-aligned bounding rectangle of all the CrossSections' * vertices. */ Rect CrossSection::Bounds() const { - auto r = C2::GetBounds(GetPaths()); + auto r = C2::GetBounds(GetPaths()->paths_); return Rect({r.left, r.bottom}, {r.right, r.top}); } @@ -713,7 +730,7 @@ Rect CrossSection::Bounds() const { */ Polygons CrossSection::ToPolygons() const { auto polys = Polygons(); - auto paths = GetPaths(); + auto paths = GetPaths()->paths_; polys.reserve(paths.size()); for (auto p : paths) { auto sp = SimplePolygon(); diff --git a/src/manifold/CMakeLists.txt b/src/manifold/CMakeLists.txt index 687cfb554..518d990dd 100644 --- a/src/manifold/CMakeLists.txt +++ b/src/manifold/CMakeLists.txt @@ -26,3 +26,15 @@ target_link_libraries(${PROJECT_NAME} target_compile_features(${PROJECT_NAME} PUBLIC cxx_std_17 ) + +install( + TARGETS ${PROJECT_NAME} + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT runtime +) +install( + DIRECTORY include/ + DESTINATION ${CMAKE_INSTALL_FULL_INCLUDEDIR}/manifold + COMPONENT devel +) + diff --git a/src/manifold/src/boolean3.cpp b/src/manifold/src/boolean3.cpp index a7cb8014d..44c70ec3d 100644 --- a/src/manifold/src/boolean3.cpp +++ b/src/manifold/src/boolean3.cpp @@ -71,7 +71,7 @@ struct CopyFaceEdges { // const int *p1q1; // x can be either vert or edge (0 or 1). SparseIndices &pXq1; - const Halfedge *halfedgesQ; + VecView halfedgesQ; void operator()(thrust::tuple in) { int idx = 3 * thrust::get<0>(in); @@ -95,9 +95,9 @@ SparseIndices Filter11(const Manifold::Impl &inP, const Manifold::Impl &inQ, ExecutionPolicy policy) { SparseIndices p1q1(3 * p1q2.size() + 3 * p2q1.size()); for_each_n(policy, zip(countAt(0), countAt(0)), p1q2.size(), - CopyFaceEdges({p1q2, p1q1, inQ.halfedge_.cptrD()})); + CopyFaceEdges({p1q2, p1q1, inQ.halfedge_})); for_each_n(policy, zip(countAt(p1q2.size()), countAt(0)), p2q1.size(), - CopyFaceEdges({p2q1, p1q1, inP.halfedge_.cptrD()})); + CopyFaceEdges({p2q1, p1q1, inP.halfedge_})); p1q1.Unique(policy); return p1q1; } @@ -116,9 +116,9 @@ bool Shadows(float p, float q, float dir) { return p == q ? dir < 0 : p < q; } * by the shared policy_ member. */ thrust::pair Shadow01( - const int p0, const int q1, const glm::vec3 *vertPosP, - const glm::vec3 *vertPosQ, const Halfedge *halfedgeQ, const float expandP, - const glm::vec3 *normalP, const bool reverse) { + const int p0, const int q1, VecView vertPosP, + VecView vertPosQ, VecView halfedgeQ, + const float expandP, VecView normalP, const bool reverse) { const int q1s = halfedgeQ[q1].startVert; const int q1e = halfedgeQ[q1].endVert; const float p0x = vertPosP[p0].x; @@ -148,13 +148,12 @@ thrust::pair Shadow01( // https://github.com/scandum/binary_search/blob/master/README.md // much faster than standard binary search on large arrays -int monobound_quaternary_search(const int64_t *array, unsigned int array_size, - int64_t key) { - if (array_size == 0) { +int monobound_quaternary_search(VecView array, int64_t key) { + if (array.size() == 0) { return -1; } unsigned int bot = 0; - unsigned int top = array_size; + unsigned int top = array.size(); while (top >= 65536) { unsigned int mid = top / 4; top -= mid * 3; @@ -187,12 +186,12 @@ int monobound_quaternary_search(const int64_t *array, unsigned int array_size, } struct Kernel11 { - const glm::vec3 *vertPosP; - const glm::vec3 *vertPosQ; - const Halfedge *halfedgeP; - const Halfedge *halfedgeQ; + VecView vertPosP; + VecView vertPosQ; + VecView halfedgeP; + VecView halfedgeQ; float expandP; - const glm::vec3 *normalP; + VecView normalP; const SparseIndices &p1q1; void operator()(thrust::tuple inout) { @@ -269,18 +268,17 @@ struct Kernel11 { } }; -std::tuple, VecDH> Shadow11(SparseIndices &p1q1, - const Manifold::Impl &inP, - const Manifold::Impl &inQ, - float expandP, - ExecutionPolicy policy) { - VecDH s11(p1q1.size()); - VecDH xyzz11(p1q1.size()); +std::tuple, Vec> Shadow11(SparseIndices &p1q1, + const Manifold::Impl &inP, + const Manifold::Impl &inQ, + float expandP, + ExecutionPolicy policy) { + Vec s11(p1q1.size()); + Vec xyzz11(p1q1.size()); for_each_n(policy, zip(countAt(0), xyzz11.begin(), s11.begin()), p1q1.size(), - Kernel11({inP.vertPos_.cptrD(), inQ.vertPos_.cptrD(), - inP.halfedge_.cptrD(), inQ.halfedge_.cptrD(), expandP, - inP.vertNormal_.cptrD(), p1q1})); + Kernel11({inP.vertPos_, inQ.vertPos_, inP.halfedge_, inQ.halfedge_, + expandP, inP.vertNormal_, p1q1})); p1q1.KeepFinite(xyzz11, s11); @@ -288,12 +286,12 @@ std::tuple, VecDH> Shadow11(SparseIndices &p1q1, }; struct Kernel02 { - const glm::vec3 *vertPosP; - const Halfedge *halfedgeQ; - const glm::vec3 *vertPosQ; + VecView vertPosP; + VecView halfedgeQ; + VecView vertPosQ; const bool forward; const float expandP; - const glm::vec3 *vertNormalP; + VecView vertNormalP; const SparseIndices &p0q2; void operator()(thrust::tuple inout) { @@ -364,20 +362,18 @@ struct Kernel02 { } }; -std::tuple, VecDH> Shadow02(const Manifold::Impl &inP, - const Manifold::Impl &inQ, - SparseIndices &p0q2, bool forward, - float expandP, - ExecutionPolicy policy) { - VecDH s02(p0q2.size()); - VecDH z02(p0q2.size()); +std::tuple, Vec> Shadow02(const Manifold::Impl &inP, + const Manifold::Impl &inQ, + SparseIndices &p0q2, bool forward, + float expandP, + ExecutionPolicy policy) { + Vec s02(p0q2.size()); + Vec z02(p0q2.size()); - auto vertNormalP = - forward ? inP.vertNormal_.cptrD() : inQ.vertNormal_.cptrD(); - for_each_n( - policy, zip(countAt(0), s02.begin(), z02.begin()), p0q2.size(), - Kernel02({inP.vertPos_.cptrD(), inQ.halfedge_.cptrD(), - inQ.vertPos_.cptrD(), forward, expandP, vertNormalP, p0q2})); + auto vertNormalP = forward ? inP.vertNormal_ : inQ.vertNormal_; + for_each_n(policy, zip(countAt(0), s02.begin(), z02.begin()), p0q2.size(), + Kernel02({inP.vertPos_, inQ.halfedge_, inQ.vertPos_, forward, + expandP, vertNormalP, p0q2})); p0q2.KeepFinite(z02, s02); @@ -385,17 +381,15 @@ std::tuple, VecDH> Shadow02(const Manifold::Impl &inP, }; struct Kernel12 { - const VecDH &p0q2; - const int *s02; - const float *z02; - const int size02; - const VecDH &p1q1; - const int *s11; - const glm::vec4 *xyzz11; - const int size11; - const Halfedge *halfedgesP; - const Halfedge *halfedgesQ; - const glm::vec3 *vertPosP; + const Vec &p0q2; + VecView s02; + VecView z02; + const Vec &p1q1; + VecView s11; + VecView xyzz11; + VecView halfedgesP; + VecView halfedgesQ; + VecView vertPosP; const bool forward; const SparseIndices &p1q2; @@ -419,7 +413,7 @@ struct Kernel12 { for (int vert : {edge.startVert, edge.endVert}) { const int64_t key = forward ? SparseIndices::EncodePQ(vert, q2) : SparseIndices::EncodePQ(q2, vert); - const int idx = monobound_quaternary_search(p0q2.ptrD(), size02, key); + const int idx = monobound_quaternary_search(p0q2, key); if (idx != -1) { const int s = s02[idx]; x12 += s * ((vert == edge.startVert) == forward ? 1 : -1); @@ -440,7 +434,7 @@ struct Kernel12 { const int q1F = edge.IsForward() ? q1 : edge.pairedHalfedge; const int64_t key = forward ? SparseIndices::EncodePQ(p1, q1F) : SparseIndices::EncodePQ(q1F, p1); - const int idx = monobound_quaternary_search(p1q1.ptrD(), size11, key); + const int idx = monobound_quaternary_search(p1q1, key); if (idx != -1) { // s is implicitly zero for anything not found const int s = s11[idx]; x12 -= s * (edge.IsForward() ? 1 : -1); @@ -476,33 +470,32 @@ struct Kernel12 { } }; -std::tuple, VecDH> Intersect12( - const Manifold::Impl &inP, const Manifold::Impl &inQ, const VecDH &s02, - const SparseIndices &p0q2, const VecDH &s11, const SparseIndices &p1q1, - const VecDH &z02, const VecDH &xyzz11, - SparseIndices &p1q2, bool forward, ExecutionPolicy policy) { - VecDH x12(p1q2.size()); - VecDH v12(p1q2.size()); +std::tuple, Vec> Intersect12( + const Manifold::Impl &inP, const Manifold::Impl &inQ, const Vec &s02, + const SparseIndices &p0q2, const Vec &s11, const SparseIndices &p1q1, + const Vec &z02, const Vec &xyzz11, SparseIndices &p1q2, + bool forward, ExecutionPolicy policy) { + Vec x12(p1q2.size()); + Vec v12(p1q2.size()); - for_each_n(policy, zip(countAt(0), x12.begin(), v12.begin()), p1q2.size(), - Kernel12({p0q2.AsVec64(), s02.ptrD(), z02.cptrD(), p0q2.size(), - p1q1.AsVec64(), s11.ptrD(), xyzz11.cptrD(), p1q1.size(), - inP.halfedge_.cptrD(), inQ.halfedge_.cptrD(), - inP.vertPos_.cptrD(), forward, p1q2})); + for_each_n( + policy, zip(countAt(0), x12.begin(), v12.begin()), p1q2.size(), + Kernel12({p0q2.AsVec64(), s02, z02, p1q1.AsVec64(), s11, xyzz11, + inP.halfedge_, inQ.halfedge_, inP.vertPos_, forward, p1q2})); p1q2.KeepFinite(v12, x12); return std::make_tuple(x12, v12); }; -VecDH Winding03(const Manifold::Impl &inP, VecDH &vertices, - VecDH &s02, bool reverse, ExecutionPolicy policy) { +Vec Winding03(const Manifold::Impl &inP, Vec &vertices, Vec &s02, + bool reverse, ExecutionPolicy policy) { // verts that are not shadowed (not in p0q2) have winding number zero. - VecDH w03(inP.NumVert(), 0); + Vec w03(inP.NumVert(), 0); // checking is slow, so just sort and reduce stable_sort_by_key(policy, vertices.begin(), vertices.end(), s02.begin()); - VecDH w03val(w03.size()); - VecDH w03vert(w03.size()); + Vec w03val(w03.size()); + Vec w03vert(w03.size()); // sum known s02 values into w03 (winding number) auto endPair = reduce_by_key< thrust::pair>( @@ -575,20 +568,20 @@ Boolean3::Boolean3(const Manifold::Impl &inP, const Manifold::Impl &inQ, // Level 2 // Build up XY-projection intersection of two edges, including the z-value for // each edge, keeping only those whose intersection exists. - VecDH s11; - VecDH xyzz11; + Vec s11; + Vec xyzz11; std::tie(s11, xyzz11) = Shadow11(p1q1, inP, inQ, expandP_, policy_); PRINT("s11 size = " << s11.size()); // Build up Z-projection of vertices onto triangles, keeping only those that // fall inside the triangle. - VecDH s02; - VecDH z02; + Vec s02; + Vec z02; std::tie(s02, z02) = Shadow02(inP, inQ, p0q2, true, expandP_, policy_); PRINT("s02 size = " << s02.size()); - VecDH s20; - VecDH z20; + Vec s20; + Vec z20; std::tie(s20, z20) = Shadow02(inQ, inP, p2q0, false, expandP_, policy_); PRINT("s20 size = " << s20.size()); @@ -604,9 +597,9 @@ Boolean3::Boolean3(const Manifold::Impl &inP, const Manifold::Impl &inQ, xyzz11, p2q1_, false, policy_); PRINT("x21 size = " << x21_.size()); - VecDH p0 = p0q2.Copy(false); + Vec p0 = p0q2.Copy(false); p0q2.Resize(0); - VecDH q0 = p2q0.Copy(true); + Vec q0 = p2q0.Copy(true); p2q0.Resize(0); // Sum up the winding numbers of all vertices. w03_ = Winding03(inP, p0, s02, false, policy_); diff --git a/src/manifold/src/boolean3.h b/src/manifold/src/boolean3.h index 84fd2839e..d9d222551 100644 --- a/src/manifold/src/boolean3.h +++ b/src/manifold/src/boolean3.h @@ -54,8 +54,8 @@ class Boolean3 { const Manifold::Impl &inP_, &inQ_; const float expandP_; SparseIndices p1q2_, p2q1_; - VecDH x12_, x21_, w03_, w30_; - VecDH v12_, v21_; + Vec x12_, x21_, w03_, w30_; + Vec v12_, v21_; ExecutionPolicy policy_; }; } // namespace manifold diff --git a/src/manifold/src/boolean_result.cpp b/src/manifold/src/boolean_result.cpp index 179d029c0..b6d4d3a1f 100644 --- a/src/manifold/src/boolean_result.cpp +++ b/src/manifold/src/boolean_result.cpp @@ -49,7 +49,7 @@ struct AbsSum : public thrust::binary_function { }; struct DuplicateVerts { - glm::vec3 *vertPosR; + VecView vertPosR; void operator()(thrust::tuple in) { int inclusion = abs(thrust::get<0>(in)); @@ -63,8 +63,8 @@ struct DuplicateVerts { }; struct CountVerts { - int *count; - const int *inclusion; + VecView count; + VecView inclusion; void operator()(const Halfedge &edge) { AtomicAdd(count[edge.face], glm::abs(inclusion[edge.startVert])); @@ -73,10 +73,10 @@ struct CountVerts { template struct CountNewVerts { - int *countP; - int *countQ; + VecView countP; + VecView countQ; const SparseIndices &pq; - const Halfedge *halfedges; + VecView halfedges; void operator()(thrust::tuple in) { int edgeP = pq.Get(thrust::get<0>(in), inverted); @@ -94,27 +94,27 @@ struct NotZero : public thrust::unary_function { int operator()(int x) const { return x > 0 ? 1 : 0; } }; -std::tuple, VecDH> SizeOutput( +std::tuple, Vec> SizeOutput( Manifold::Impl &outR, const Manifold::Impl &inP, const Manifold::Impl &inQ, - const VecDH &i03, const VecDH &i30, const VecDH &i12, - const VecDH &i21, const SparseIndices &p1q2, const SparseIndices &p2q1, + const Vec &i03, const Vec &i30, const Vec &i12, + const Vec &i21, const SparseIndices &p1q2, const SparseIndices &p2q1, bool invertQ, ExecutionPolicy policy) { - VecDH sidesPerFacePQ(inP.NumTri() + inQ.NumTri(), 0); - auto sidesPerFaceP = sidesPerFacePQ.ptrD(); - auto sidesPerFaceQ = sidesPerFacePQ.ptrD() + inP.NumTri(); + Vec sidesPerFacePQ(inP.NumTri() + inQ.NumTri(), 0); + auto sidesPerFaceP = sidesPerFacePQ.view(0, inP.NumTri()); + auto sidesPerFaceQ = sidesPerFacePQ.view(inP.NumTri(), inQ.NumTri()); for_each(policy, inP.halfedge_.begin(), inP.halfedge_.end(), - CountVerts({sidesPerFaceP, i03.cptrD()})); + CountVerts({sidesPerFaceP, i03})); for_each(policy, inQ.halfedge_.begin(), inQ.halfedge_.end(), - CountVerts({sidesPerFaceQ, i30.cptrD()})); + CountVerts({sidesPerFaceQ, i30})); for_each_n(policy, zip(countAt(0), i12.begin()), i12.size(), CountNewVerts( - {sidesPerFaceP, sidesPerFaceQ, p1q2, inP.halfedge_.cptrD()})); - for_each_n(policy, zip(countAt(0), i21.begin()), i21.size(), - CountNewVerts( - {sidesPerFaceQ, sidesPerFaceP, p2q1, inQ.halfedge_.cptrD()})); + {sidesPerFaceP, sidesPerFaceQ, p1q2, inP.halfedge_})); + for_each_n( + policy, zip(countAt(0), i21.begin()), i21.size(), + CountNewVerts({sidesPerFaceQ, sidesPerFaceP, p2q1, inQ.halfedge_})); - VecDH facePQ2R(inP.NumTri() + inQ.NumTri() + 1, 0); + Vec facePQ2R(inP.NumTri() + inQ.NumTri() + 1, 0); auto keepFace = thrust::make_transform_iterator(sidesPerFacePQ.begin(), NotZero()); inclusive_scan(policy, keepFace, keepFace + sidesPerFacePQ.size(), @@ -142,7 +142,7 @@ std::tuple, VecDH> SizeOutput( auto newEnd = remove( policy, sidesPerFacePQ.begin(), sidesPerFacePQ.end(), 0); - VecDH faceEdge(newEnd - sidesPerFacePQ.begin() + 1, 0); + Vec faceEdge(newEnd - sidesPerFacePQ.begin() + 1, 0); inclusive_scan(policy, sidesPerFacePQ.begin(), newEnd, faceEdge.begin() + 1); outR.halfedge_.resize(faceEdge.back()); @@ -159,15 +159,15 @@ void AddNewEdgeVerts( // we need concurrent_map because we will be adding things concurrently concurrent_map> &edgesP, concurrent_map, std::vector> &edgesNew, - const SparseIndices &p1q2, const VecDH &i12, const VecDH &v12R, - const VecDH &halfedgeP, bool forward) { + const SparseIndices &p1q2, const Vec &i12, const Vec &v12R, + const Vec &halfedgeP, bool forward) { // For each edge of P that intersects a face of Q (p1q2), add this vertex to // P's corresponding edge vector and to the two new edges, which are // intersections between the face of Q and the two faces of P attached to the // edge. The direction and duplicity are given by i12, while v12R remaps to // the output vert index. When forward is false, all is reversed. - const VecDH &p1 = p1q2.Copy(!forward); - const VecDH &q2 = p1q2.Copy(forward); + const Vec &p1 = p1q2.Copy(!forward); + const Vec &q2 = p1q2.Copy(forward); auto process = [&](std::function lock, std::function unlock, int i) { const int edgeP = p1[i]; @@ -246,20 +246,20 @@ std::vector PairUp(std::vector &edgePos) { return edges; } -void AppendPartialEdges(Manifold::Impl &outR, VecDH &wholeHalfedgeP, - VecDH &facePtrR, +void AppendPartialEdges(Manifold::Impl &outR, Vec &wholeHalfedgeP, + Vec &facePtrR, concurrent_map> &edgesP, - VecDH &halfedgeRef, const Manifold::Impl &inP, - const VecDH &i03, const VecDH &vP2R, - const VecDH::IterC faceP2R, bool forward) { + Vec &halfedgeRef, const Manifold::Impl &inP, + const Vec &i03, const Vec &vP2R, + const Vec::IterC faceP2R, bool forward) { // Each edge in the map is partially retained; for each of these, look up // their original verts and include them based on their winding number (i03), // while remapping them to the output using vP2R. Use the verts position // projected along the edge vector to pair them up, then distribute these // edges to their faces. - VecDH &halfedgeR = outR.halfedge_; - const VecDH &vertPosP = inP.vertPos_; - const VecDH &halfedgeP = inP.halfedge_; + Vec &halfedgeR = outR.halfedge_; + const Vec &vertPosP = inP.vertPos_; + const Vec &halfedgeP = inP.halfedge_; for (auto &value : edgesP) { const int edgeP = value.first; @@ -330,13 +330,12 @@ void AppendPartialEdges(Manifold::Impl &outR, VecDH &wholeHalfedgeP, } void AppendNewEdges( - Manifold::Impl &outR, VecDH &facePtrR, + Manifold::Impl &outR, Vec &facePtrR, concurrent_map, std::vector> &edgesNew, - VecDH &halfedgeRef, const VecDH &facePQ2R, - const int numFaceP) { + Vec &halfedgeRef, const Vec &facePQ2R, const int numFaceP) { // Pair up each edge's verts and distribute to faces based on indices in key. - VecDH &halfedgeR = outR.halfedge_; - VecDH &vertPosR = outR.vertPos_; + Vec &halfedgeR = outR.halfedge_; + Vec &vertPosR = outR.vertPos_; for (auto &value : edgesNew) { const int faceP = value.first.first; @@ -383,13 +382,13 @@ void AppendNewEdges( } struct DuplicateHalfedges { - Halfedge *halfedgesR; - TriRef *halfedgeRef; - int *facePtr; - const Halfedge *halfedgesP; - const int *i03; - const int *vP2R; - const int *faceP2R; + VecView halfedgesR; + VecView halfedgeRef; + VecView facePtr; + VecView halfedgesP; + VecView i03; + VecView vP2R; + VecView faceP2R; const bool forward; void operator()(thrust::tuple in) { @@ -434,22 +433,21 @@ struct DuplicateHalfedges { } }; -void AppendWholeEdges(Manifold::Impl &outR, VecDH &facePtrR, - VecDH &halfedgeRef, const Manifold::Impl &inP, - const VecDH wholeHalfedgeP, const VecDH &i03, - const VecDH &vP2R, const int *faceP2R, bool forward, - ExecutionPolicy policy) { +void AppendWholeEdges(Manifold::Impl &outR, Vec &facePtrR, + Vec &halfedgeRef, const Manifold::Impl &inP, + const Vec wholeHalfedgeP, const Vec &i03, + const Vec &vP2R, VecView faceP2R, + bool forward, ExecutionPolicy policy) { for_each_n(policy, zip(wholeHalfedgeP.begin(), inP.halfedge_.begin(), countAt(0)), inP.halfedge_.size(), - DuplicateHalfedges({outR.halfedge_.ptrD(), halfedgeRef.ptrD(), - facePtrR.ptrD(), inP.halfedge_.cptrD(), - i03.cptrD(), vP2R.cptrD(), faceP2R, forward})); + DuplicateHalfedges({outR.halfedge_, halfedgeRef, facePtrR, + inP.halfedge_, i03, vP2R, faceP2R, forward})); } struct MapTriRef { - const TriRef *triRefP; - const TriRef *triRefQ; + VecView triRefP; + VecView triRefQ; const int offsetQ; void operator()(TriRef &triRef) { @@ -460,14 +458,14 @@ struct MapTriRef { } }; -VecDH UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP, - const Manifold::Impl &inQ, bool invertQ, - ExecutionPolicy policy) { - VecDH refPQ = outR.meshRelation_.triRef; +Vec UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP, + const Manifold::Impl &inQ, bool invertQ, + ExecutionPolicy policy) { + Vec refPQ = outR.meshRelation_.triRef; const int offsetQ = Manifold::Impl::meshIDCounter_; - for_each_n(policy, outR.meshRelation_.triRef.begin(), outR.NumTri(), - MapTriRef({inP.meshRelation_.triRef.cptrD(), - inQ.meshRelation_.triRef.cptrD(), offsetQ})); + for_each_n( + policy, outR.meshRelation_.triRef.begin(), outR.NumTri(), + MapTriRef({inP.meshRelation_.triRef, inQ.meshRelation_.triRef, offsetQ})); for (const auto &pair : inP.meshRelation_.meshIDtransform) { outR.meshRelation_.meshIDtransform[pair.first] = pair.second; @@ -481,13 +479,13 @@ VecDH UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP, } struct Barycentric { - glm::vec3 *uvw; - const glm::vec3 *vertPosP; - const glm::vec3 *vertPosQ; - const glm::vec3 *vertPosR; - const Halfedge *halfedgeP; - const Halfedge *halfedgeQ; - const Halfedge *halfedgeR; + VecView uvw; + VecView vertPosP; + VecView vertPosQ; + VecView vertPosR; + VecView halfedgeP; + VecView halfedgeQ; + VecView halfedgeR; const float precision; void operator()(thrust::tuple in) { @@ -511,7 +509,7 @@ struct Barycentric { } }; -void CreateProperties(Manifold::Impl &outR, const VecDH &refPQ, +void CreateProperties(Manifold::Impl &outR, const Vec &refPQ, const Manifold::Impl &inP, const Manifold::Impl &inQ, ExecutionPolicy policy) { const int numPropP = inP.NumProp(); @@ -523,12 +521,11 @@ void CreateProperties(Manifold::Impl &outR, const VecDH &refPQ, const int numTri = outR.NumTri(); outR.meshRelation_.triProperties.resize(numTri); - VecDH bary(outR.halfedge_.size()); + Vec bary(outR.halfedge_.size()); for_each_n(policy, zip(countAt(0), refPQ.cbegin()), numTri, - Barycentric({bary.ptrD(), inP.vertPos_.cptrD(), - inQ.vertPos_.cptrD(), outR.vertPos_.cptrD(), - inP.halfedge_.cptrD(), inQ.halfedge_.cptrD(), - outR.halfedge_.cptrD(), outR.precision_})); + Barycentric({bary, inP.vertPos_, inQ.vertPos_, outR.vertPos_, + inP.halfedge_, inQ.halfedge_, outR.halfedge_, + outR.precision_})); using Entry = std::pair; int idMissProp = outR.NumVert(); @@ -632,28 +629,28 @@ Manifold::Impl Boolean3::Result(OpType op) const { const bool invertQ = op == OpType::Subtract; // Convert winding numbers to inclusion values based on operation type. - VecDH i12(x12_.size()); - VecDH i21(x21_.size()); - VecDH i03(w03_.size()); - VecDH i30(w30_.size()); + Vec i12(x12_.size()); + Vec i21(x21_.size()); + Vec i03(w03_.size()); + Vec i30(w30_.size()); transform(policy_, x12_.begin(), x12_.end(), i12.begin(), c3 * _1); transform(policy_, x21_.begin(), x21_.end(), i21.begin(), c3 * _1); transform(policy_, w03_.begin(), w03_.end(), i03.begin(), c1 + c3 * _1); transform(policy_, w30_.begin(), w30_.end(), i30.begin(), c2 + c3 * _1); - VecDH vP2R(inP_.NumVert()); + Vec vP2R(inP_.NumVert()); exclusive_scan(policy_, i03.begin(), i03.end(), vP2R.begin(), 0, AbsSum()); int numVertR = AbsSum()(vP2R.back(), i03.back()); const int nPv = numVertR; - VecDH vQ2R(inQ_.NumVert()); + Vec vQ2R(inQ_.NumVert()); exclusive_scan(policy_, i30.begin(), i30.end(), vQ2R.begin(), numVertR, AbsSum()); numVertR = AbsSum()(vQ2R.back(), i30.back()); const int nQv = numVertR - nPv; - VecDH v12R(v12_.size()); + Vec v12R(v12_.size()); if (v12_.size() > 0) { exclusive_scan(policy_, i12.begin(), i12.end(), v12R.begin(), numVertR, AbsSum()); @@ -661,7 +658,7 @@ Manifold::Impl Boolean3::Result(OpType op) const { } const int n12 = numVertR - nPv - nQv; - VecDH v21R(v21_.size()); + Vec v21R(v21_.size()); if (v21_.size() > 0) { exclusive_scan(policy_, i21.begin(), i21.end(), v21R.begin(), numVertR, AbsSum()); @@ -680,14 +677,14 @@ Manifold::Impl Boolean3::Result(OpType op) const { // Add vertices, duplicating for inclusion numbers not in [-1, 1]. // Retained vertices from P and Q: for_each_n(policy_, zip(i03.begin(), vP2R.begin(), inP_.vertPos_.begin()), - inP_.NumVert(), DuplicateVerts({outR.vertPos_.ptrD()})); + inP_.NumVert(), DuplicateVerts({outR.vertPos_})); for_each_n(policy_, zip(i30.begin(), vQ2R.begin(), inQ_.vertPos_.begin()), - inQ_.NumVert(), DuplicateVerts({outR.vertPos_.ptrD()})); + inQ_.NumVert(), DuplicateVerts({outR.vertPos_})); // New vertices created from intersections: for_each_n(policy_, zip(i12.begin(), v12R.begin(), v12_.begin()), i12.size(), - DuplicateVerts({outR.vertPos_.ptrD()})); + DuplicateVerts({outR.vertPos_})); for_each_n(policy_, zip(i21.begin(), v21R.begin(), v21_.begin()), i21.size(), - DuplicateVerts({outR.vertPos_.ptrD()})); + DuplicateVerts({outR.vertPos_})); PRINT(nPv << " verts from inP"); PRINT(nQv << " verts from inQ"); @@ -709,21 +706,21 @@ Manifold::Impl Boolean3::Result(OpType op) const { AddNewEdgeVerts(edgesQ, edgesNew, p2q1_, i21, v21R, inQ_.halfedge_, false); // Level 4 - VecDH faceEdge; - VecDH facePQ2R; + Vec faceEdge; + Vec facePQ2R; std::tie(faceEdge, facePQ2R) = SizeOutput( outR, inP_, inQ_, i03, i30, i12, i21, p1q2_, p2q1_, invertQ, policy_); const int numFaceR = faceEdge.size() - 1; // This gets incremented for each halfedge that's added to a face so that the // next one knows where to slot in. - VecDH facePtrR = faceEdge; + Vec facePtrR = faceEdge; // Intersected halfedges are marked false. - VecDH wholeHalfedgeP(inP_.halfedge_.size(), true); - VecDH wholeHalfedgeQ(inQ_.halfedge_.size(), true); + Vec wholeHalfedgeP(inP_.halfedge_.size(), true); + Vec wholeHalfedgeQ(inQ_.halfedge_.size(), true); // The halfedgeRef contains the data that will become triRef once the faces // are triangulated. - VecDH halfedgeRef(2 * outR.NumEdge()); + Vec halfedgeRef(2 * outR.NumEdge()); AppendPartialEdges(outR, wholeHalfedgeP, facePtrR, edgesP, halfedgeRef, inP_, i03, vP2R, facePQ2R.begin(), true); @@ -734,9 +731,10 @@ Manifold::Impl Boolean3::Result(OpType op) const { inP_.NumTri()); AppendWholeEdges(outR, facePtrR, halfedgeRef, inP_, wholeHalfedgeP, i03, vP2R, - facePQ2R.cptrD(), true, policy_); + facePQ2R.cview(0, inP_.NumTri()), true, policy_); AppendWholeEdges(outR, facePtrR, halfedgeRef, inQ_, wholeHalfedgeQ, i30, vQ2R, - facePQ2R.cptrD() + inP_.NumTri(), false, policy_); + facePQ2R.cview(inP_.NumTri(), inQ_.NumTri()), false, + policy_); #ifdef MANIFOLD_DEBUG assemble.Stop(); @@ -760,7 +758,7 @@ Manifold::Impl Boolean3::Result(OpType op) const { if (ManifoldParams().intermediateChecks) ASSERT(outR.IsManifold(), logicErr, "triangulated mesh is not manifold!"); - VecDH refPQ = UpdateReference(outR, inP_, inQ_, invertQ, policy_); + Vec refPQ = UpdateReference(outR, inP_, inQ_, invertQ, policy_); outR.SimplifyTopology(); diff --git a/src/manifold/src/constructors.cpp b/src/manifold/src/constructors.cpp index d1c9d32d3..186478a1b 100644 --- a/src/manifold/src/constructors.cpp +++ b/src/manifold/src/constructors.cpp @@ -40,8 +40,8 @@ struct Equals { }; struct RemoveFace { - const Halfedge* halfedge; - const int* vertLabel; + VecView halfedge; + VecView vertLabel; const int keepLabel; bool operator()(int face) { @@ -247,7 +247,7 @@ Manifold Manifold::Extrude(const CrossSection& crossSection, float height, auto pImpl_ = std::make_shared(); ++nDivisions; auto& vertPos = pImpl_->vertPos_; - VecDH triVertsDH; + Vec triVertsDH; auto& triVerts = triVertsDH; int nCrossSection = 0; bool isCone = scaleTop.x == 0.0 && scaleTop.y == 0.0; @@ -359,7 +359,7 @@ Manifold Manifold::Revolve(const CrossSection& crossSection, auto pImpl_ = std::make_shared(); auto& vertPos = pImpl_->vertPos_; - VecDH triVertsDH; + Vec triVertsDH; auto& triVerts = triVertsDH; std::vector startPoses; @@ -482,7 +482,7 @@ std::vector Manifold::Decompose() const { meshes[0] = *this; return meshes; } - VecDH vertLabel(componentIndices); + Vec vertLabel(componentIndices); std::vector meshes; for (int i = 0; i < numComponents; ++i) { @@ -490,7 +490,7 @@ std::vector Manifold::Decompose() const { // inherit original object's precision impl->precision_ = pImpl_->precision_; impl->vertPos_.resize(NumVert()); - VecDH vertNew2Old(NumVert()); + Vec vertNew2Old(NumVert()); auto policy = autoPolicy(NumVert()); auto start = zip(impl->vertPos_.begin(), vertNew2Old.begin()); int nVert = @@ -501,14 +501,13 @@ std::vector Manifold::Decompose() const { start; impl->vertPos_.resize(nVert); - VecDH faceNew2Old(NumTri()); + Vec faceNew2Old(NumTri()); sequence(policy, faceNew2Old.begin(), faceNew2Old.end()); - int nFace = - remove_if( - policy, faceNew2Old.begin(), faceNew2Old.end(), - RemoveFace({pImpl_->halfedge_.cptrD(), vertLabel.cptrD(), i})) - - faceNew2Old.begin(); + int nFace = remove_if( + policy, faceNew2Old.begin(), faceNew2Old.end(), + RemoveFace({pImpl_->halfedge_, vertLabel, i})) - + faceNew2Old.begin(); faceNew2Old.resize(nFace); impl->GatherFaces(*pImpl_, faceNew2Old); diff --git a/src/manifold/src/csg_tree.cpp b/src/manifold/src/csg_tree.cpp index 383275b0e..c89c7d761 100644 --- a/src/manifold/src/csg_tree.cpp +++ b/src/manifold/src/csg_tree.cpp @@ -72,7 +72,7 @@ struct UpdateMeshIDs { }; struct CheckOverlap { - const Box *boxes; + VecView boxes; const size_t i; bool operator()(int j) { return boxes[i].DoesOverlap(boxes[j]); } }; @@ -233,9 +233,9 @@ Manifold::Impl CsgLeafNode::Compose( auto &oldProp = node->pImpl_->meshRelation_.properties; auto &newProp = combined.meshRelation_.properties; for (int p = 0; p < numProp; ++p) { - strided_range::IterC> oldRange( - oldProp.begin() + p, oldProp.end(), numProp); - strided_range::Iter> newRange( + strided_range::IterC> oldRange(oldProp.begin() + p, + oldProp.end(), numProp); + strided_range::Iter> newRange( newProp.begin() + numPropOut * propVertIndices[i] + p, newProp.end(), numPropOut); copy(policy, oldRange.begin(), oldRange.end(), newRange.begin()); @@ -275,14 +275,13 @@ Manifold::Impl CsgLeafNode::Compose( countAt(0)), node->pImpl_->halfedgeTangent_.size(), TransformTangents{glm::mat3(node->transform_), invert, - node->pImpl_->halfedgeTangent_.cptrD(), - node->pImpl_->halfedge_.cptrD()}); + node->pImpl_->halfedgeTangent_, + node->pImpl_->halfedge_}); if (invert) for_each_n(policy, zip(combined.meshRelation_.triRef.begin(), countAt(triIndices[i])), - node->pImpl_->NumTri(), - FlipTris({combined.halfedge_.ptrD()})); + node->pImpl_->NumTri(), FlipTris({combined.halfedge_})); } // Since the nodes may be copies containing the same meshIDs, it is // important to add an offset so that each node instance gets @@ -532,22 +531,21 @@ void CsgOpNode::BatchUnion() const { const int start = (children_.size() > kMaxUnionSize) ? (children_.size() - kMaxUnionSize) : 0; - VecDH boxes; + Vec boxes; boxes.reserve(children_.size() - start); for (int i = start; i < children_.size(); i++) { boxes.push_back(std::dynamic_pointer_cast(children_[i]) ->GetImpl() ->bBox_); } - const Box *boxesD = boxes.cptrD(); // partition the children into a set of disjoint sets // each set contains a set of children that are pairwise disjoint - std::vector> disjointSets; + std::vector> disjointSets; for (size_t i = 0; i < boxes.size(); i++) { - auto lambda = [boxesD, i](const VecDH &set) { + auto lambda = [&boxes, i](const Vec &set) { return find_if( autoPolicy(set.size()), set.begin(), set.end(), - CheckOverlap({boxesD, i})) == set.end(); + CheckOverlap({boxes, i})) == set.end(); }; auto it = std::find_if(disjointSets.begin(), disjointSets.end(), lambda); if (it == disjointSets.end()) { diff --git a/src/manifold/src/edge_op.cpp b/src/manifold/src/edge_op.cpp index 3a023722c..4e3db366b 100644 --- a/src/manifold/src/edge_op.cpp +++ b/src/manifold/src/edge_op.cpp @@ -45,8 +45,8 @@ struct DuplicateEdge { }; struct ShortEdge { - const Halfedge* halfedge; - const glm::vec3* vertPos; + VecView halfedge; + VecView vertPos; const float precision; bool operator()(int edge) const { @@ -59,8 +59,8 @@ struct ShortEdge { }; struct FlagEdge { - const Halfedge* halfedge; - const TriRef* triRef; + VecView halfedge; + VecView triRef; bool operator()(int edge) const { if (halfedge[edge].pairedHalfedge < 0) return false; @@ -82,9 +82,9 @@ struct FlagEdge { }; struct SwappableEdge { - const Halfedge* halfedge; - const glm::vec3* vertPos; - const glm::vec3* triNormal; + VecView halfedge; + VecView vertPos; + VecView triNormal; const float precision; bool operator()(int edge) const { @@ -147,17 +147,15 @@ void Manifold::Impl::SimplifyTopology() { int nbEdges = halfedge_.size(); int numFlagged = 0; - VecDH bflags(nbEdges); - uint8_t* bflagsPtr = bflags.ptrD(); + Vec bflags(nbEdges); - VecDH entries(nbEdges); - SortEntry* entriesPtr = entries.ptrD(); + Vec entries(nbEdges); auto policy = autoPolicy(halfedge_.size()); - for_each_n(policy, countAt(0), nbEdges, [=] __host__ __device__(int i) { - entriesPtr[i].start = halfedge_[i].startVert; - entriesPtr[i].end = halfedge_[i].endVert; - entriesPtr[i].index = i; + for_each_n(policy, countAt(0), nbEdges, [&](int i) { + entries[i].start = halfedge_[i].startVert; + entries[i].end = halfedge_[i].endVert; + entries[i].index = i; }); sort(policy, entries.begin(), entries.end()); @@ -180,9 +178,8 @@ void Manifold::Impl::SimplifyTopology() { scratchBuffer.reserve(10); { numFlagged = 0; - ShortEdge se{halfedge_.cptrD(), vertPos_.cptrD(), precision_}; - for_each_n(policy, countAt(0), nbEdges, - [=] __host__ __device__(int i) { bflagsPtr[i] = se(i); }); + ShortEdge se{halfedge_, vertPos_, precision_}; + for_each_n(policy, countAt(0), nbEdges, [&](int i) { bflags[i] = se(i); }); for (int i = 0; i < nbEdges; ++i) { if (bflags[i]) { CollapseEdge(i, scratchBuffer); @@ -201,9 +198,8 @@ void Manifold::Impl::SimplifyTopology() { { numFlagged = 0; - FlagEdge se{halfedge_.cptrD(), meshRelation_.triRef.cptrD()}; - for_each_n(policy, countAt(0), nbEdges, - [=] __host__ __device__(int i) { bflagsPtr[i] = se(i); }); + FlagEdge se{halfedge_, meshRelation_.triRef}; + for_each_n(policy, countAt(0), nbEdges, [&](int i) { bflags[i] = se(i); }); for (int i = 0; i < nbEdges; ++i) { if (bflags[i]) { CollapseEdge(i, scratchBuffer); @@ -222,10 +218,8 @@ void Manifold::Impl::SimplifyTopology() { { numFlagged = 0; - SwappableEdge se{halfedge_.cptrD(), vertPos_.cptrD(), faceNormal_.cptrD(), - precision_}; - for_each_n(policy, countAt(0), nbEdges, - [=] __host__ __device__(int i) { bflagsPtr[i] = se(i); }); + SwappableEdge se{halfedge_, vertPos_, faceNormal_, precision_}; + for_each_n(policy, countAt(0), nbEdges, [&](int i) { bflags[i] = se(i); }); std::vector edgeSwapStack; std::vector visited(halfedge_.size(), -1); int tag = 0; @@ -386,8 +380,8 @@ void Manifold::Impl::RemoveIfFolded(int edge) { } void Manifold::Impl::CollapseEdge(const int edge, std::vector& edges) { - VecDH& triRef = meshRelation_.triRef; - VecDH& triProp = meshRelation_.triProperties; + Vec& triRef = meshRelation_.triRef; + Vec& triProp = meshRelation_.triProperties; const Halfedge toRemove = halfedge_[edge]; if (toRemove.pairedHalfedge < 0) return; @@ -499,7 +493,7 @@ void Manifold::Impl::RecursiveEdgeSwap(const int edge, int& tag, std::vector& visited, std::vector& edgeSwapStack, std::vector& edges) { - VecDH& triRef = meshRelation_.triRef; + Vec& triRef = meshRelation_.triRef; if (edge < 0) return; const int pair = halfedge_[edge].pairedHalfedge; @@ -548,8 +542,8 @@ void Manifold::Impl::RecursiveEdgeSwap(const int edge, int& tag, const float a = glm::max(0.0f, glm::min(1.0f, l02 / l01)); // Update properties if applicable if (meshRelation_.properties.size() > 0) { - VecDH& triProp = meshRelation_.triProperties; - VecDH& prop = meshRelation_.properties; + Vec& triProp = meshRelation_.triProperties; + Vec& prop = meshRelation_.properties; triProp[tri0] = triProp[tri1]; triProp[tri0][perm0[1]] = triProp[tri1][perm1[0]]; triProp[tri0][perm0[0]] = triProp[tri1][perm1[2]]; diff --git a/src/manifold/src/face_op.cpp b/src/manifold/src/face_op.cpp index 13e571857..6d94a3576 100644 --- a/src/manifold/src/face_op.cpp +++ b/src/manifold/src/face_op.cpp @@ -37,11 +37,11 @@ using AddTriangle = std::function; * represents the mesh as a set of triangles as usual. In this process the * faceNormal_ values are retained, repeated as necessary. */ -void Manifold::Impl::Face2Tri(const VecDH& faceEdge, - const VecDH& halfedgeRef) { - VecDH triVerts; - VecDH triNormal; - VecDH& triRef = meshRelation_.triRef; +void Manifold::Impl::Face2Tri(const Vec& faceEdge, + const Vec& halfedgeRef) { + Vec triVerts; + Vec triNormal; + Vec& triRef = meshRelation_.triRef; triRef.resize(0); auto processFace = [&](GeneralTriangulation general, AddTriangle addTri, int face) { @@ -133,7 +133,7 @@ void Manifold::Impl::Face2Tri(const VecDH& faceEdge, tbb::task_group group; // map from face to triangle tbb::concurrent_unordered_map> results; - VecDH triCount(faceEdge.size()); + Vec triCount(faceEdge.size()); triCount.back() = 0; // precompute number of triangles per face, and launch async tasks to // triangulate complex faces @@ -196,7 +196,7 @@ void Manifold::Impl::Face2Tri(const VecDH& faceEdge, * projection of the vertices. */ PolygonsIdx Manifold::Impl::Face2Polygons(int face, glm::mat3x2 projection, - const VecDH& faceEdge) const { + const Vec& faceEdge) const { const int firstEdge = faceEdge[face]; const int lastEdge = faceEdge[face + 1]; diff --git a/src/manifold/src/impl.cpp b/src/manifold/src/impl.cpp index c52d97479..be03c2601 100644 --- a/src/manifold/src/impl.cpp +++ b/src/manifold/src/impl.cpp @@ -54,9 +54,9 @@ struct Transform4x3 { }; struct AssignNormals { - glm::vec3* vertNormal; - const glm::vec3* vertPos; - const Halfedge* halfedges; + VecView vertNormal; + VecView vertPos; + VecView halfedges; const float precision; const bool calculateTriNormal; @@ -94,8 +94,8 @@ struct AssignNormals { }; struct Tri2Halfedges { - Halfedge* halfedges; - glm::uint64_t* edges; + VecView halfedges; + VecView edges; void operator()(thrust::tuple in) { const int tri = thrust::get<0>(in); @@ -114,8 +114,8 @@ struct Tri2Halfedges { }; struct LinkHalfedges { - Halfedge* halfedges; - const int* ids; + VecView halfedges; + VecView ids; const int numEdge; void operator()(int i) { @@ -127,18 +127,18 @@ struct LinkHalfedges { }; struct MarkVerts { - int* vert; + VecView vert; void operator()(glm::ivec3 triVerts) { for (int i : {0, 1, 2}) { - reinterpret_cast*>(vert)[triVerts[i]].store( - 1, std::memory_order_relaxed); + reinterpret_cast*>(&vert[triVerts[i]]) + ->store(1, std::memory_order_relaxed); } } }; struct ReindexTriVerts { - const int* old2new; + VecView old2new; void operator()(glm::ivec3& triVerts) { for (int i : {0, 1, 2}) { @@ -149,7 +149,7 @@ struct ReindexTriVerts { struct InitializeTriRef { const int meshID; - const Halfedge* halfedge; + VecView halfedge; void operator()(thrust::tuple inOut) { TriRef& baryRef = thrust::get<0>(inOut); @@ -168,13 +168,13 @@ struct UpdateMeshID { }; struct CoplanarEdge { - float* triArea; - const Halfedge* halfedge; - const glm::vec3* vertPos; - const TriRef* triRef; - const glm::ivec3* triProp; - const float* prop; - const float* propTol; + VecView triArea; + VecView halfedge; + VecView vertPos; + VecView triRef; + VecView triProp; + VecView prop; + VecView propTol; const int numProp; const float precision; @@ -271,10 +271,10 @@ struct CoplanarEdge { }; struct CheckCoplanarity { - int* comp2tri; - const Halfedge* halfedge; - const glm::vec3* vertPos; - const int* components; + VecView comp2tri; + VecView halfedge; + VecView vertPos; + VecView components; const float precision; void operator()(int tri) { @@ -293,8 +293,8 @@ struct CheckCoplanarity { // triangle, unmark the entire component so that none of its triangles are // marked coplanar. if (glm::abs(glm::dot(normal, vert - origin)) > precision) { - reinterpret_cast*>(comp2tri)[component].store( - -1, std::memory_order_relaxed); + reinterpret_cast*>(&comp2tri[component]) + ->store(-1, std::memory_order_relaxed); break; } } @@ -302,7 +302,7 @@ struct CheckCoplanarity { }; struct EdgeBox { - const glm::vec3* vertPos; + VecView vertPos; void operator()(thrust::tuple inout) { const TmpEdge& edge = thrust::get<1>(inout); @@ -311,7 +311,7 @@ struct EdgeBox { }; int GetLabels(std::vector& components, - const VecDH>& edges, int numNodes) { + const Vec>& edges, int numNodes) { Graph graph; for (int i = 0; i < numNodes; ++i) { graph.add_nodes(i); @@ -325,8 +325,8 @@ int GetLabels(std::vector& components, return ConnectedComponents(components, graph); } -void DedupePropVerts(manifold::VecDH& triProp, - const VecDH>& vert2vert) { +void DedupePropVerts(manifold::Vec& triProp, + const Vec>& vert2vert) { std::vector vertLabels; const int numLabels = GetLabels(vertLabels, vert2vert, vert2vert.size()); @@ -488,7 +488,7 @@ Manifold::Impl::Impl(const Mesh& mesh, const MeshRelationD& relation, meshRelation_ = {relation.originalID, relation.numProp, relation.properties, relation.meshIDtransform}; - VecDH triVerts; + Vec triVerts; for (int i = 0; i < mesh.triVerts.size(); ++i) { const glm::ivec3 tri = mesh.triVerts[i]; // Remove topological degenerates @@ -587,13 +587,13 @@ Manifold::Impl::Impl(Shape shape) { CreateFaces(); } -void Manifold::Impl::RemoveUnreferencedVerts(VecDH& triVerts) { - VecDH vertOld2New(NumVert() + 1, 0); +void Manifold::Impl::RemoveUnreferencedVerts(Vec& triVerts) { + Vec vertOld2New(NumVert() + 1, 0); auto policy = autoPolicy(NumVert()); for_each(policy, triVerts.cbegin(), triVerts.cend(), - MarkVerts({vertOld2New.ptrD() + 1})); + MarkVerts({vertOld2New.view(1)})); - const VecDH oldVertPos = vertPos_; + const Vec oldVertPos = vertPos_; vertPos_.resize(copy_if( policy, oldVertPos.cbegin(), oldVertPos.cend(), vertOld2New.cbegin() + 1, vertPos_.begin(), @@ -604,7 +604,7 @@ void Manifold::Impl::RemoveUnreferencedVerts(VecDH& triVerts) { vertOld2New.begin() + 1); for_each(policy, triVerts.begin(), triVerts.end(), - ReindexTriVerts({vertOld2New.cptrD()})); + ReindexTriVerts({vertOld2New})); } void Manifold::Impl::InitializeOriginal() { @@ -614,28 +614,25 @@ void Manifold::Impl::InitializeOriginal() { meshRelation_.triRef.resize(NumTri()); for_each_n(autoPolicy(NumTri()), zip(meshRelation_.triRef.begin(), countAt(0)), NumTri(), - InitializeTriRef({meshID, halfedge_.cptrD()})); + InitializeTriRef({meshID, halfedge_})); meshRelation_.meshIDtransform.clear(); meshRelation_.meshIDtransform[meshID] = {meshID}; } void Manifold::Impl::CreateFaces(const std::vector& propertyTolerance) { - VecDH propertyToleranceD = - propertyTolerance.empty() - ? VecDH(meshRelation_.numProp, kTolerance) - : propertyTolerance; - - VecDH> face2face(halfedge_.size(), {-1, -1}); - VecDH> vert2vert(halfedge_.size(), {-1, -1}); - VecDH triArea(NumTri()); + Vec propertyToleranceD = + propertyTolerance.empty() ? Vec(meshRelation_.numProp, kTolerance) + : propertyTolerance; + + Vec> face2face(halfedge_.size(), {-1, -1}); + Vec> vert2vert(halfedge_.size(), {-1, -1}); + Vec triArea(NumTri()); for_each_n( autoPolicy(halfedge_.size()), zip(face2face.begin(), vert2vert.begin(), countAt(0)), halfedge_.size(), - CoplanarEdge( - {triArea.ptrD(), halfedge_.cptrD(), vertPos_.cptrD(), - meshRelation_.triRef.cptrD(), meshRelation_.triProperties.cptrD(), - meshRelation_.properties.cptrD(), propertyToleranceD.cptrD(), - meshRelation_.numProp, precision_})); + CoplanarEdge({triArea, halfedge_, vertPos_, meshRelation_.triRef, + meshRelation_.triProperties, meshRelation_.properties, + propertyToleranceD, meshRelation_.numProp, precision_})); if (meshRelation_.triProperties.size() > 0) { DedupePropVerts(meshRelation_.triProperties, vert2vert); @@ -654,14 +651,13 @@ void Manifold::Impl::CreateFaces(const std::vector& propertyTolerance) { } } - VecDH componentsD(components); - VecDH comp2triD(comp2tri); - for_each_n( - autoPolicy(halfedge_.size()), countAt(0), NumTri(), - CheckCoplanarity({comp2triD.ptrD(), halfedge_.cptrD(), vertPos_.cptrD(), - componentsD.cptrD(), precision_})); + Vec componentsD(components); + Vec comp2triD(comp2tri); + for_each_n(autoPolicy(halfedge_.size()), countAt(0), NumTri(), + CheckCoplanarity( + {comp2triD, halfedge_, vertPos_, componentsD, precision_})); - VecDH& triRef = meshRelation_.triRef; + Vec& triRef = meshRelation_.triRef; for (int tri = 0; tri < NumTri(); ++tri) { const int referenceTri = comp2triD[components[tri]]; if (referenceTri >= 0) { @@ -673,18 +669,18 @@ void Manifold::Impl::CreateFaces(const std::vector& propertyTolerance) { /** * Create the halfedge_ data structure from an input triVerts array like Mesh. */ -void Manifold::Impl::CreateHalfedges(const VecDH& triVerts) { +void Manifold::Impl::CreateHalfedges(const Vec& triVerts) { const int numTri = triVerts.size(); const int numHalfedge = 3 * numTri; // drop the old value first to avoid copy halfedge_.resize(0); halfedge_.resize(numHalfedge); - VecDH edge(numHalfedge); - VecDH ids(numHalfedge); + Vec edge(numHalfedge); + Vec ids(numHalfedge); auto policy = autoPolicy(numTri); sequence(policy, ids.begin(), ids.end()); for_each_n(policy, zip(countAt(0), triVerts.begin()), numTri, - Tri2Halfedges({halfedge_.ptrD(), edge.ptrD()})); + Tri2Halfedges({halfedge_, edge})); // Stable sort is required here so that halfedges from the same face are // paired together (the triangles were created in face order). In some // degenerate situations the triangulator can add the same internal edge in @@ -695,7 +691,7 @@ void Manifold::Impl::CreateHalfedges(const VecDH& triVerts) { // correspond to their backward pair at the same offset in the second half // of the range. for_each_n(policy, countAt(0), numHalfedge / 2, - LinkHalfedges({halfedge_.ptrD(), ids.ptrD(), numHalfedge / 2})); + LinkHalfedges({halfedge_, ids, numHalfedge / 2})); } /** @@ -704,8 +700,8 @@ void Manifold::Impl::CreateHalfedges(const VecDH& triVerts) { */ void Manifold::Impl::Update() { CalculateBBox(); - VecDH faceBox; - VecDH faceMorton; + Vec faceBox; + Vec faceMorton; GetFaceBoxMorton(faceBox, faceMorton); collider_.UpdateBoxes(faceBox); } @@ -768,16 +764,15 @@ Manifold::Impl Manifold::Impl::Transform(const glm::mat4x3& transform_) const { const bool invert = glm::determinant(glm::mat3(transform_)) < 0; if (halfedgeTangent_.size() > 0) { - for_each_n( - policy, zip(result.halfedgeTangent_.begin(), countAt(0)), - halfedgeTangent_.size(), - TransformTangents({glm::mat3(transform_), invert, - halfedgeTangent_.cptrD(), halfedge_.cptrD()})); + for_each_n(policy, zip(result.halfedgeTangent_.begin(), countAt(0)), + halfedgeTangent_.size(), + TransformTangents({glm::mat3(transform_), invert, + halfedgeTangent_, halfedge_})); } if (invert) { for_each_n(policy, zip(result.meshRelation_.triRef.begin(), countAt(0)), - result.NumTri(), FlipTris({result.halfedge_.ptrD()})); + result.NumTri(), FlipTris({result.halfedge_})); } // This optimization does a cheap collider update if the transform is @@ -824,10 +819,9 @@ void Manifold::Impl::CalculateNormals() { faceNormal_.resize(NumTri()); calculateTriNormal = true; } - for_each_n( - policy, zip(faceNormal_.begin(), countAt(0)), NumTri(), - AssignNormals({vertNormal_.ptrD(), vertPos_.cptrD(), halfedge_.cptrD(), - precision_, calculateTriNormal})); + for_each_n(policy, zip(faceNormal_.begin(), countAt(0)), NumTri(), + AssignNormals({vertNormal_, vertPos_, halfedge_, precision_, + calculateTriNormal})); for_each(policy, vertNormal_.begin(), vertNormal_.end(), Normalize()); } @@ -859,25 +853,25 @@ void Manifold::Impl::IncrementMeshIDs() { */ SparseIndices Manifold::Impl::EdgeCollisions(const Impl& Q, bool inverted) const { - VecDH edges = CreateTmpEdges(Q.halfedge_); + Vec edges = CreateTmpEdges(Q.halfedge_); const int numEdge = edges.size(); - VecDH QedgeBB(numEdge); + Vec QedgeBB(numEdge); auto policy = autoPolicy(numEdge); for_each_n(policy, zip(QedgeBB.begin(), edges.cbegin()), numEdge, - EdgeBox({Q.vertPos_.cptrD()})); + EdgeBox({Q.vertPos_})); SparseIndices q1p2(0); if (inverted) - q1p2 = collider_.Collisions(QedgeBB); + q1p2 = collider_.Collisions(QedgeBB.cview()); else - q1p2 = collider_.Collisions(QedgeBB); + q1p2 = collider_.Collisions(QedgeBB.cview()); if (inverted) for_each(policy, countAt(0), countAt(q1p2.size()), - ReindexEdge({edges.cptrD(), q1p2})); + ReindexEdge({edges, q1p2})); else for_each(policy, countAt(0), countAt(q1p2.size()), - ReindexEdge({edges.cptrD(), q1p2})); + ReindexEdge({edges, q1p2})); return q1p2; } @@ -885,8 +879,8 @@ SparseIndices Manifold::Impl::EdgeCollisions(const Impl& Q, * Returns a sparse array of the input vertices that project inside the XY * bounding boxes of the faces of this manifold. */ -SparseIndices Manifold::Impl::VertexCollisionsZ(const VecDH& vertsIn, - bool inverted) const { +SparseIndices Manifold::Impl::VertexCollisionsZ( + VecView vertsIn, bool inverted) const { if (inverted) return collider_.Collisions(vertsIn); else diff --git a/src/manifold/src/impl.h b/src/manifold/src/impl.h index 375bd8a45..28fbfb961 100644 --- a/src/manifold/src/impl.h +++ b/src/manifold/src/impl.h @@ -22,7 +22,7 @@ #include "shared.h" #include "sparse.h" #include "utils.h" -#include "vec_dh.h" +#include "vec.h" namespace manifold { @@ -37,20 +37,20 @@ struct Manifold::Impl { /// The originalID of this Manifold if it is an original; -1 otherwise. int originalID = -1; int numProp = 0; - VecDH properties; + Vec properties; std::map meshIDtransform; - VecDH triRef; - VecDH triProperties; + Vec triRef; + Vec triProperties; }; Box bBox_; float precision_ = -1; Error status_ = Error::NoError; - VecDH vertPos_; - VecDH halfedge_; - VecDH vertNormal_; - VecDH faceNormal_; - VecDH halfedgeTangent_; + Vec vertPos_; + Vec halfedge_; + Vec vertNormal_; + Vec faceNormal_; + Vec halfedgeTangent_; MeshRelationD meshRelation_; Collider collider_; @@ -67,9 +67,9 @@ struct Manifold::Impl { bool hasFaceIDs = false); void CreateFaces(const std::vector& propertyTolerance = {}); - void RemoveUnreferencedVerts(VecDH& triVerts); + void RemoveUnreferencedVerts(Vec& triVerts); void InitializeOriginal(); - void CreateHalfedges(const VecDH& triVerts); + void CreateHalfedges(const Vec& triVerts); void CalculateNormals(); void IncrementMeshIDs(); @@ -78,7 +78,7 @@ struct Manifold::Impl { void Warp(std::function warpFunc); Impl Transform(const glm::mat4x3& transform) const; SparseIndices EdgeCollisions(const Impl& B, bool inverted = false) const; - SparseIndices VertexCollisionsZ(const VecDH& vertsIn, + SparseIndices VertexCollisionsZ(VecView vertsIn, bool inverted = false) const; bool IsEmpty() const { return NumVert() == 0; } @@ -96,7 +96,7 @@ struct Manifold::Impl { void CalculateCurvature(int gaussianIdx, int meanIdx); void CalculateBBox(); bool IsFinite() const; - bool IsIndexInBounds(const VecDH& triVerts) const; + bool IsIndexInBounds(VecView triVerts) const; void SetPrecision(float minPrecision = -1); bool IsManifold() const; bool Is2Manifold() const; @@ -106,17 +106,17 @@ struct Manifold::Impl { // sort.cu void Finish(); void SortVerts(); - void ReindexVerts(const VecDH& vertNew2Old, int numOldVert); + void ReindexVerts(const Vec& vertNew2Old, int numOldVert); void CompactProps(); - void GetFaceBoxMorton(VecDH& faceBox, VecDH& faceMorton) const; - void SortFaces(VecDH& faceBox, VecDH& faceMorton); - void GatherFaces(const VecDH& faceNew2Old); - void GatherFaces(const Impl& old, const VecDH& faceNew2Old); + void GetFaceBoxMorton(Vec& faceBox, Vec& faceMorton) const; + void SortFaces(Vec& faceBox, Vec& faceMorton); + void GatherFaces(const Vec& faceNew2Old); + void GatherFaces(const Impl& old, const Vec& faceNew2Old); // face_op.cu - void Face2Tri(const VecDH& faceEdge, const VecDH& halfedgeRef); + void Face2Tri(const Vec& faceEdge, const Vec& halfedgeRef); PolygonsIdx Face2Polygons(int face, glm::mat3x2 projection, - const VecDH& faceEdge) const; + const Vec& faceEdge) const; // edge_op.cu void SimplifyTopology(); @@ -134,7 +134,7 @@ struct Manifold::Impl { // smoothing.cu void CreateTangents(const std::vector&); - VecDH Subdivide(int n); + Vec Subdivide(int n); void Refine(int n); }; } // namespace manifold diff --git a/src/manifold/src/manifold.cpp b/src/manifold/src/manifold.cpp index 35d2ac22a..2753242b8 100644 --- a/src/manifold/src/manifold.cpp +++ b/src/manifold/src/manifold.cpp @@ -12,11 +12,11 @@ // See the License for the specific language governing permissions and // limitations under the License. -#include #include #include #include +#include "QuickHull.hpp" #include "boolean3.h" #include "csg_tree.h" #include "impl.h" @@ -29,7 +29,7 @@ using namespace thrust::placeholders; ExecutionParams manifoldParams; struct MakeTri { - const Halfedge* halfedges; + VecView halfedges; void operator()(thrust::tuple inOut) { glm::ivec3& tri = thrust::get<0>(inOut); @@ -170,7 +170,7 @@ Mesh Manifold::GetMesh() const { result.triVerts.resize(NumTri()); // note that `triVerts` is `std::vector`, so we cannot use thrust::device thrust::for_each_n(thrust::host, zip(result.triVerts.begin(), countAt(0)), - NumTri(), MakeTri({impl.halfedge_.cptrH()})); + NumTri(), MakeTri({impl.halfedge_})); return result; } @@ -218,7 +218,7 @@ MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const { out.faceID.resize(numTri); std::vector triNew2Old(numTri); std::iota(triNew2Old.begin(), triNew2Old.end(), 0); - const TriRef* triRef = impl.meshRelation_.triRef.cptrD(); + VecView triRef = impl.meshRelation_.triRef; // Don't sort originals - keep them in order if (!isOriginal) { std::sort(triNew2Old.begin(), triNew2Old.end(), [triRef](int a, int b) { @@ -583,7 +583,7 @@ Manifold Manifold::SetProperties( propFunc) const { auto pImpl = std::make_shared(*GetCsgLeafNode().GetImpl()); const int oldNumProp = NumProp(); - const VecDH oldProperties = pImpl->meshRelation_.properties; + const Vec oldProperties = pImpl->meshRelation_.properties; auto& triProperties = pImpl->meshRelation_.triProperties; if (numProp == 0) { @@ -599,17 +599,16 @@ Manifold Manifold::SetProperties( triProperties[i][j] = idx++; } } - pImpl->meshRelation_.properties = VecDH(numProp * idx, 0); + pImpl->meshRelation_.properties = Vec(numProp * idx, 0); } else { - pImpl->meshRelation_.properties = - VecDH(numProp * NumPropVert(), 0); + pImpl->meshRelation_.properties = Vec(numProp * NumPropVert(), 0); } thrust::for_each_n( thrust::host, countAt(0), NumTri(), - UpdateProperties({pImpl->meshRelation_.properties.ptrH(), numProp, - oldProperties.ptrH(), oldNumProp, - pImpl->vertPos_.ptrH(), triProperties.ptrH(), - pImpl->halfedge_.ptrH(), propFunc})); + UpdateProperties({pImpl->meshRelation_.properties.data(), numProp, + oldProperties.data(), oldNumProp, + pImpl->vertPos_.data(), triProperties.data(), + pImpl->halfedge_.data(), propFunc})); } pImpl->meshRelation_.numProp = numProp; diff --git a/src/manifold/src/mesh_fixes.h b/src/manifold/src/mesh_fixes.h index af0e7a387..d8cb4de89 100644 --- a/src/manifold/src/mesh_fixes.h +++ b/src/manifold/src/mesh_fixes.h @@ -22,8 +22,8 @@ struct TransformNormals { struct TransformTangents { const glm::mat3 transform; const bool invert; - const glm::vec4* oldTangents; - const Halfedge* halfedge; + VecView oldTangents; + VecView halfedge; void operator()(thrust::tuple inOut) { glm::vec4& tangent = thrust::get<0>(inOut); @@ -38,7 +38,7 @@ struct TransformTangents { }; struct FlipTris { - Halfedge* halfedge; + VecView halfedge; void operator()(thrust::tuple inOut) { TriRef& bary = thrust::get<0>(inOut); diff --git a/src/manifold/src/properties.cpp b/src/manifold/src/properties.cpp index 2196a766c..42dd6d37a 100644 --- a/src/manifold/src/properties.cpp +++ b/src/manifold/src/properties.cpp @@ -25,8 +25,8 @@ namespace { using namespace manifold; struct FaceAreaVolume { - const Halfedge* halfedges; - const glm::vec3* vertPos; + VecView halfedges; + VecView vertPos; const float precision; thrust::pair operator()(int face) { @@ -97,13 +97,13 @@ struct SumPair : public thrust::binary_function, }; struct CurvatureAngles { - float* meanCurvature; - float* gaussianCurvature; - float* area; - float* degree; - const Halfedge* halfedge; - const glm::vec3* vertPos; - const glm::vec3* triNormal; + VecView meanCurvature; + VecView gaussianCurvature; + VecView area; + VecView degree; + VecView halfedge; + VecView vertPos; + VecView triNormal; void operator()(int tri) { glm::vec3 edge[3]; @@ -152,12 +152,12 @@ struct NormalizeCurvature { }; struct UpdateProperties { - float* properties; + VecView properties; - const float* oldProperties; - const Halfedge* halfedge; - const float* meanCurvature; - const float* gaussianCurvature; + VecView oldProperties; + VecView halfedge; + VecView meanCurvature; + VecView gaussianCurvature; const int oldNumProp; const int numProp; const int gaussianIdx; @@ -191,7 +191,7 @@ struct UpdateProperties { }; struct CheckHalfedges { - const Halfedge* halfedges; + VecView halfedges; bool operator()(int edge) { const Halfedge halfedge = halfedges[edge]; @@ -209,7 +209,7 @@ struct CheckHalfedges { }; struct NoDuplicates { - const Halfedge* halfedges; + VecView halfedges; bool operator()(int edge) { const Halfedge halfedge = halfedges[edge]; @@ -222,9 +222,9 @@ struct NoDuplicates { }; struct CheckCCW { - const Halfedge* halfedges; - const glm::vec3* vertPos; - const glm::vec3* triNormal; + VecView halfedges; + VecView vertPos; + VecView triNormal; const float tol; bool operator()(int face) { @@ -276,7 +276,7 @@ bool Manifold::Impl::IsManifold() const { auto policy = autoPolicy(halfedge_.size()); return all_of(policy, countAt(0), countAt(halfedge_.size()), - CheckHalfedges({halfedge_.cptrD()})); + CheckHalfedges({halfedge_})); } /** @@ -289,11 +289,11 @@ bool Manifold::Impl::Is2Manifold() const { if (!IsManifold()) return false; - VecDH halfedge(halfedge_); + Vec halfedge(halfedge_); sort(policy, halfedge.begin(), halfedge.end()); return all_of(policy, countAt(0), countAt(2 * NumEdge() - 1), - NoDuplicates({halfedge.cptrD()})); + NoDuplicates({halfedge})); } /** @@ -302,8 +302,7 @@ bool Manifold::Impl::Is2Manifold() const { bool Manifold::Impl::MatchesTriNormals() const { if (halfedge_.size() == 0 || faceNormal_.size() != NumTri()) return true; return all_of(autoPolicy(NumTri()), countAt(0), countAt(NumTri()), - CheckCCW({halfedge_.cptrD(), vertPos_.cptrD(), - faceNormal_.cptrD(), 2 * precision_})); + CheckCCW({halfedge_, vertPos_, faceNormal_, 2 * precision_})); } /** @@ -311,16 +310,16 @@ bool Manifold::Impl::MatchesTriNormals() const { */ int Manifold::Impl::NumDegenerateTris() const { if (halfedge_.size() == 0 || faceNormal_.size() != NumTri()) return true; - return count_if(autoPolicy(NumTri()), countAt(0), countAt(NumTri()), - CheckCCW({halfedge_.cptrD(), vertPos_.cptrD(), - faceNormal_.cptrD(), -1 * precision_ / 2})); + return count_if( + autoPolicy(NumTri()), countAt(0), countAt(NumTri()), + CheckCCW({halfedge_, vertPos_, faceNormal_, -1 * precision_ / 2})); } Properties Manifold::Impl::GetProperties() const { if (IsEmpty()) return {0, 0}; auto areaVolume = transform_reduce>( autoPolicy(NumTri()), countAt(0), countAt(NumTri()), - FaceAreaVolume({halfedge_.cptrD(), vertPos_.cptrD(), precision_}), + FaceAreaVolume({halfedge_, vertPos_, precision_}), thrust::make_pair(0.0f, 0.0f), SumPair()); return {areaVolume.first, areaVolume.second}; } @@ -328,16 +327,14 @@ Properties Manifold::Impl::GetProperties() const { void Manifold::Impl::CalculateCurvature(int gaussianIdx, int meanIdx) { if (IsEmpty()) return; if (gaussianIdx < 0 && meanIdx < 0) return; - VecDH vertMeanCurvature(NumVert(), 0); - VecDH vertGaussianCurvature(NumVert(), glm::two_pi()); - VecDH vertArea(NumVert(), 0); - VecDH degree(NumVert(), 0); + Vec vertMeanCurvature(NumVert(), 0); + Vec vertGaussianCurvature(NumVert(), glm::two_pi()); + Vec vertArea(NumVert(), 0); + Vec degree(NumVert(), 0); auto policy = autoPolicy(NumTri()); - for_each( - policy, countAt(0), countAt(NumTri()), - CurvatureAngles({vertMeanCurvature.ptrD(), vertGaussianCurvature.ptrD(), - vertArea.ptrD(), degree.ptrD(), halfedge_.cptrD(), - vertPos_.cptrD(), faceNormal_.cptrD()})); + for_each(policy, countAt(0), countAt(NumTri()), + CurvatureAngles({vertMeanCurvature, vertGaussianCurvature, vertArea, + degree, halfedge_, vertPos_, faceNormal_})); for_each_n(policy, zip(vertMeanCurvature.begin(), vertGaussianCurvature.begin(), vertArea.begin(), degree.begin()), @@ -345,8 +342,8 @@ void Manifold::Impl::CalculateCurvature(int gaussianIdx, int meanIdx) { const int oldNumProp = NumProp(); const int numProp = glm::max(oldNumProp, glm::max(gaussianIdx, meanIdx) + 1); - const VecDH oldProperties = meshRelation_.properties; - meshRelation_.properties = VecDH(numProp * NumPropVert(), 0); + const Vec oldProperties = meshRelation_.properties; + meshRelation_.properties = Vec(numProp * NumPropVert(), 0); meshRelation_.numProp = numProp; if (meshRelation_.triProperties.size() == 0) { meshRelation_.triProperties.resize(NumTri()); @@ -354,10 +351,9 @@ void Manifold::Impl::CalculateCurvature(int gaussianIdx, int meanIdx) { for_each_n( policy, zip(meshRelation_.triProperties.begin(), countAt(0)), NumTri(), - UpdateProperties({meshRelation_.properties.ptrD(), oldProperties.cptrD(), - halfedge_.cptrD(), vertMeanCurvature.cptrD(), - vertGaussianCurvature.cptrD(), oldNumProp, numProp, - gaussianIdx, meanIdx})); + UpdateProperties({meshRelation_.properties, oldProperties, halfedge_, + vertMeanCurvature, vertGaussianCurvature, oldNumProp, + numProp, gaussianIdx, meanIdx})); CreateFaces(); SimplifyTopology(); @@ -394,7 +390,7 @@ bool Manifold::Impl::IsFinite() const { * Checks that the input triVerts array has all indices inside bounds of the * vertPos_ array. */ -bool Manifold::Impl::IsIndexInBounds(const VecDH& triVerts) const { +bool Manifold::Impl::IsIndexInBounds(VecView triVerts) const { auto policy = autoPolicy(triVerts.size()); glm::ivec2 minmax = transform_reduce( policy, triVerts.begin(), triVerts.end(), MakeMinMax(), diff --git a/src/manifold/src/shared.h b/src/manifold/src/shared.h index d1304dfac..8bcf19b25 100644 --- a/src/manifold/src/shared.h +++ b/src/manifold/src/shared.h @@ -17,7 +17,7 @@ #include "par.h" #include "sparse.h" #include "utils.h" -#include "vec_dh.h" +#include "vec.h" namespace manifold { @@ -193,8 +193,8 @@ struct TmpInvalid { bool operator()(const TmpEdge& edge) { return edge.halfedgeIdx < 0; } }; -VecDH inline CreateTmpEdges(const VecDH& halfedge) { - VecDH edges(halfedge.size()); +Vec inline CreateTmpEdges(const Vec& halfedge) { + Vec edges(halfedge.size()); for_each_n(autoPolicy(edges.size()), zip(edges.begin(), halfedge.begin(), countAt(0)), edges.size(), Halfedge2Tmp()); @@ -209,7 +209,7 @@ VecDH inline CreateTmpEdges(const VecDH& halfedge) { template struct ReindexEdge { - const TmpEdge* edges; + VecView edges; SparseIndices& indices; void operator()(int i) { diff --git a/src/manifold/src/smoothing.cpp b/src/manifold/src/smoothing.cpp index 170114a04..8272eb0d3 100644 --- a/src/manifold/src/smoothing.cpp +++ b/src/manifold/src/smoothing.cpp @@ -39,8 +39,8 @@ int VertsPerTri(int edgeVerts) { * CPU for simplicity for now. Using AtomicCAS on .tri should work for a GPU * version if desired. */ -void FillRetainedVerts(VecDH& vertBary, - const VecDH& halfedge_) { +void FillRetainedVerts(Vec& vertBary, + const Vec& halfedge_) { const int numTri = halfedge_.size() / 3; for (int tri = 0; tri < numTri; ++tri) { for (const int i : {0, 1, 2}) { @@ -52,7 +52,7 @@ void FillRetainedVerts(VecDH& vertBary, } struct ReindexHalfedge { - int* half2Edge; + VecView half2Edge; void operator()(thrust::tuple in) { const int edge = thrust::get<0>(in); @@ -63,8 +63,8 @@ struct ReindexHalfedge { }; struct EdgeVerts { - glm::vec3* vertPos; - Barycentric* vertBary; + VecView vertPos; + VecView vertBary; const int startIdx; const int n; @@ -91,11 +91,11 @@ struct EdgeVerts { }; struct InteriorVerts { - glm::vec3* vertPos; - Barycentric* vertBary; + VecView vertPos; + VecView vertBary; const int startIdx; const int n; - const Halfedge* halfedge; + VecView halfedge; void operator()(int tri) { const float invTotal = 1.0f / n; @@ -119,9 +119,9 @@ struct InteriorVerts { }; struct SplitTris { - glm::ivec3* triVerts; - const Halfedge* halfedge; - const int* half2Edge; + VecView triVerts; + VecView halfedge; + VecView half2Edge; const int edgeIdx; const int triIdx; const int n; @@ -182,10 +182,10 @@ struct SplitTris { }; struct SmoothBezier { - const glm::vec3* vertPos; - const glm::vec3* triNormal; - const glm::vec3* vertNormal; - const Halfedge* halfedge; + VecView vertPos; + VecView triNormal; + VecView vertNormal; + VecView halfedge; void operator()(thrust::tuple inOut) { glm::vec4& tangent = thrust::get<0>(inOut); @@ -212,9 +212,9 @@ struct SmoothBezier { }; struct InterpTri { - const Halfedge* halfedge; - const glm::vec4* halfedgeTangent; - const glm::vec3* vertPos; + VecView halfedge; + VecView halfedgeTangent; + VecView vertPos; glm::vec4 Homogeneous(glm::vec4 v) const { v.x *= v.w; @@ -325,11 +325,10 @@ void Manifold::Impl::CreateTangents( for_each_n(autoPolicy(numHalfedge), zip(halfedgeTangent_.begin(), halfedge_.cbegin()), numHalfedge, - SmoothBezier({vertPos_.cptrD(), faceNormal_.cptrD(), - vertNormal_.cptrD(), halfedge_.cptrD()})); + SmoothBezier({vertPos_, faceNormal_, vertNormal_, halfedge_})); if (!sharpenedEdges.empty()) { - const VecDH& triRef = meshRelation_.triRef; + const Vec& triRef = meshRelation_.triRef; // sharpenedEdges are referenced to the input Mesh, but the triangles have // been sorted in creating the Manifold, so the indices are converted using @@ -362,7 +361,7 @@ void Manifold::Impl::CreateTangents( {edge.second, edge.first}); } - VecDH& tangent = halfedgeTangent_; + Vec& tangent = halfedgeTangent_; for (const auto& value : vertTangents) { const std::vector& vert = value.second; // Sharp edges that end are smooth at their terminal vert. @@ -423,8 +422,8 @@ void Manifold::Impl::CreateTangents( * run after the new vertices have moved, which is a likely scenario after * refinement (smoothing). */ -VecDH Manifold::Impl::Subdivide(int n) { - if (n < 2) return VecDH(); +Vec Manifold::Impl::Subdivide(int n) { + if (n < 2) return Vec(); faceNormal_.resize(0); vertNormal_.resize(0); const int numVert = NumVert(); @@ -434,28 +433,27 @@ VecDH Manifold::Impl::Subdivide(int n) { const int vertsPerEdge = n - 1; const int triVertStart = numVert + numEdge * vertsPerEdge; vertPos_.resize(triVertStart + numTri * VertsPerTri(n - 2)); - VecDH vertBary(vertPos_.size()); + Vec vertBary(vertPos_.size()); FillRetainedVerts(vertBary, halfedge_); MeshRelationD oldMeshRelation = std::move(meshRelation_); meshRelation_.triRef.resize(n * n * numTri); meshRelation_.originalID = oldMeshRelation.originalID; - VecDH edges = CreateTmpEdges(halfedge_); - VecDH half2Edge(2 * numEdge); + Vec edges = CreateTmpEdges(halfedge_); + Vec half2Edge(2 * numEdge); auto policy = autoPolicy(numEdge); for_each_n(policy, zip(countAt(0), edges.begin()), numEdge, - ReindexHalfedge({half2Edge.ptrD()})); + ReindexHalfedge({half2Edge})); for_each_n(policy, zip(countAt(0), edges.begin()), numEdge, - EdgeVerts({vertPos_.ptrD(), vertBary.ptrD(), numVert, n})); + EdgeVerts({vertPos_, vertBary, numVert, n})); for_each_n(policy, countAt(0), numTri, - InteriorVerts({vertPos_.ptrD(), vertBary.ptrD(), triVertStart, n, - halfedge_.cptrD()})); + InteriorVerts({vertPos_, vertBary, triVertStart, n, halfedge_})); // Create sub-triangles - VecDH triVerts(n * n * numTri); - for_each_n(policy, countAt(0), numTri, - SplitTris({triVerts.ptrD(), halfedge_.cptrD(), half2Edge.cptrD(), - numVert, triVertStart, n})); + Vec triVerts(n * n * numTri); + for_each_n( + policy, countAt(0), numTri, + SplitTris({triVerts, halfedge_, half2Edge, numVert, triVertStart, n})); CreateHalfedges(triVerts); // Make original since the subdivided faces are intended to be warped into // being non-coplanar, and hence not being related to the original faces. @@ -475,14 +473,13 @@ VecDH Manifold::Impl::Subdivide(int n) { void Manifold::Impl::Refine(int n) { Manifold::Impl old = *this; - VecDH vertBary = Subdivide(n); + Vec vertBary = Subdivide(n); if (vertBary.size() == 0) return; if (old.halfedgeTangent_.size() == old.halfedge_.size()) { for_each_n(autoPolicy(NumTri()), zip(vertPos_.begin(), vertBary.begin()), NumVert(), - InterpTri({old.halfedge_.cptrD(), old.halfedgeTangent_.cptrD(), - old.vertPos_.cptrD()})); + InterpTri({old.halfedge_, old.halfedgeTangent_, old.vertPos_})); } halfedgeTangent_.resize(0); diff --git a/src/manifold/src/sort.cpp b/src/manifold/src/sort.cpp index 2ac1988e6..13deab440 100644 --- a/src/manifold/src/sort.cpp +++ b/src/manifold/src/sort.cpp @@ -82,8 +82,8 @@ struct Morton { }; struct FaceMortonBox { - const Halfedge* halfedge; - const glm::vec3* vertPos; + VecView halfedge; + VecView vertPos; const Box bBox; void operator()(thrust::tuple inout) { @@ -113,7 +113,7 @@ struct FaceMortonBox { }; struct Reindex { - const int* indexInv; + VecView indexInv; void operator()(Halfedge& edge) { if (edge.startVert < 0) return; @@ -123,19 +123,19 @@ struct Reindex { }; struct MarkProp { - int* keep; + VecView keep; void operator()(glm::ivec3 triProp) { for (const int i : {0, 1, 2}) { - reinterpret_cast*>(keep)[triProp[i]].store( - 1, std::memory_order_relaxed); + reinterpret_cast*>(&keep[triProp[i]]) + ->store(1, std::memory_order_relaxed); } } }; struct GatherProps { - float* properties; - const float* oldProperties; + VecView properties; + VecView oldProperties; const int numProp; void operator()(thrust::tuple in) { @@ -150,7 +150,7 @@ struct GatherProps { }; struct ReindexProps { - const int* old2new; + VecView old2new; void operator()(glm::ivec3& triProp) { for (const int i : {0, 1, 2}) { @@ -160,23 +160,23 @@ struct ReindexProps { }; template -void Permute(VecDH& inOut, const VecDH& new2Old) { - VecDH tmp(std::move(inOut)); +void Permute(Vec& inOut, const Vec& new2Old) { + Vec tmp(std::move(inOut)); inOut.resize(new2Old.size()); gather(autoPolicy(new2Old.size()), new2Old.begin(), new2Old.end(), tmp.begin(), inOut.begin()); } -template void Permute(VecDH&, const VecDH&); -template void Permute(VecDH&, const VecDH&); +template void Permute(Vec&, const Vec&); +template void Permute(Vec&, const Vec&); struct ReindexFace { - Halfedge* halfedge; - glm::vec4* halfedgeTangent; - const Halfedge* oldHalfedge; - const glm::vec4* oldHalfedgeTangent; - const int* faceNew2Old; - const int* faceOld2New; + VecView halfedge; + VecView halfedgeTangent; + VecView oldHalfedge; + VecView oldHalfedgeTangent; + VecView faceNew2Old; + VecView faceOld2New; void operator()(int newFace) { const int oldFace = faceNew2Old[newFace]; @@ -189,7 +189,7 @@ struct ReindexFace { edge.pairedHalfedge = 3 * faceOld2New[pairedFace] + offset; const int newEdge = 3 * newFace + i; halfedge[newEdge] = edge; - if (oldHalfedgeTangent != nullptr) { + if (!oldHalfedgeTangent.empty()) { halfedgeTangent[newEdge] = oldHalfedgeTangent[oldEdge]; } } @@ -197,7 +197,7 @@ struct ReindexFace { }; struct VertMortonBox { - const float* vertProperties; + VecView vertProperties; const uint32_t numProp; const float tol; const Box bBox; @@ -254,8 +254,8 @@ void Manifold::Impl::Finish() { } SortVerts(); - VecDH faceBox; - VecDH faceMorton; + Vec faceBox; + Vec faceMorton; GetFaceBoxMorton(faceBox, faceMorton); SortFaces(faceBox, faceMorton); if (halfedge_.size() == 0) return; @@ -304,12 +304,12 @@ void Manifold::Impl::Finish() { */ void Manifold::Impl::SortVerts() { const int numVert = NumVert(); - VecDH vertMorton(numVert); + Vec vertMorton(numVert); auto policy = autoPolicy(numVert); for_each_n(policy, zip(vertMorton.begin(), vertPos_.cbegin()), numVert, Morton({bBox_})); - VecDH vertNew2Old(numVert); + Vec vertNew2Old(numVert); sequence(policy, vertNew2Old.begin(), vertNew2Old.end()); sort_by_key(policy, vertMorton.begin(), vertMorton.end(), @@ -337,13 +337,12 @@ void Manifold::Impl::SortVerts() { * vertNew2Old. This may be a subset, so the total number of original verts is * also given. */ -void Manifold::Impl::ReindexVerts(const VecDH& vertNew2Old, - int oldNumVert) { - VecDH vertOld2New(oldNumVert); +void Manifold::Impl::ReindexVerts(const Vec& vertNew2Old, int oldNumVert) { + Vec vertOld2New(oldNumVert); scatter(autoPolicy(oldNumVert), countAt(0), countAt(NumVert()), vertNew2Old.begin(), vertOld2New.begin()); for_each(autoPolicy(oldNumVert), halfedge_.begin(), halfedge_.end(), - Reindex({vertOld2New.cptrD()})); + Reindex({vertOld2New})); } /** @@ -353,23 +352,22 @@ void Manifold::Impl::CompactProps() { if (meshRelation_.numProp == 0) return; const int numVerts = meshRelation_.properties.size() / meshRelation_.numProp; - VecDH keep(numVerts, 0); + Vec keep(numVerts, 0); auto policy = autoPolicy(numVerts); for_each(policy, meshRelation_.triProperties.cbegin(), - meshRelation_.triProperties.cend(), MarkProp({keep.ptrD()})); - VecDH propOld2New(numVerts + 1, 0); + meshRelation_.triProperties.cend(), MarkProp({keep})); + Vec propOld2New(numVerts + 1, 0); inclusive_scan(policy, keep.begin(), keep.end(), propOld2New.begin() + 1); - VecDH oldProp = meshRelation_.properties; + Vec oldProp = meshRelation_.properties; const int numVertsNew = propOld2New[numVerts]; meshRelation_.properties.resize(meshRelation_.numProp * numVertsNew); - for_each_n(policy, zip(countAt(0), propOld2New.cbegin(), keep.cbegin()), - numVerts, - GatherProps({meshRelation_.properties.ptrD(), oldProp.cptrD(), - meshRelation_.numProp})); + for_each_n( + policy, zip(countAt(0), propOld2New.cbegin(), keep.cbegin()), numVerts, + GatherProps({meshRelation_.properties, oldProp, meshRelation_.numProp})); for_each_n(policy, meshRelation_.triProperties.begin(), NumTri(), - ReindexProps({propOld2New.cptrD()})); + ReindexProps({propOld2New})); } /** @@ -377,22 +375,21 @@ void Manifold::Impl::CompactProps() { * codes of the faces, respectively. The Morton code is based on the center of * the bounding box. */ -void Manifold::Impl::GetFaceBoxMorton(VecDH& faceBox, - VecDH& faceMorton) const { +void Manifold::Impl::GetFaceBoxMorton(Vec& faceBox, + Vec& faceMorton) const { faceBox.resize(NumTri()); faceMorton.resize(NumTri()); for_each_n(autoPolicy(NumTri()), zip(faceMorton.begin(), faceBox.begin(), countAt(0)), NumTri(), - FaceMortonBox({halfedge_.cptrD(), vertPos_.cptrD(), bBox_})); + FaceMortonBox({halfedge_, vertPos_, bBox_})); } /** * Sorts the faces of this manifold according to their input Morton code. The * bounding box and Morton code arrays are also sorted accordingly. */ -void Manifold::Impl::SortFaces(VecDH& faceBox, - VecDH& faceMorton) { - VecDH faceNew2Old(NumTri()); +void Manifold::Impl::SortFaces(Vec& faceBox, Vec& faceMorton) { + Vec faceNew2Old(NumTri()); auto policy = autoPolicy(faceNew2Old.size()); sequence(policy, faceNew2Old.begin(), faceNew2Old.end()); @@ -417,7 +414,7 @@ void Manifold::Impl::SortFaces(VecDH& faceBox, * another manifold, given by oldHalfedge. Input faceNew2Old defines the old * faces to gather into this. */ -void Manifold::Impl::GatherFaces(const VecDH& faceNew2Old) { +void Manifold::Impl::GatherFaces(const Vec& faceNew2Old) { const int numTri = faceNew2Old.size(); if (meshRelation_.triRef.size() == NumTri()) Permute(meshRelation_.triRef, faceNew2Old); @@ -425,9 +422,9 @@ void Manifold::Impl::GatherFaces(const VecDH& faceNew2Old) { Permute(meshRelation_.triProperties, faceNew2Old); if (faceNormal_.size() == NumTri()) Permute(faceNormal_, faceNew2Old); - VecDH oldHalfedge(std::move(halfedge_)); - VecDH oldHalfedgeTangent(std::move(halfedgeTangent_)); - VecDH faceOld2New(oldHalfedge.size() / 3); + Vec oldHalfedge(std::move(halfedge_)); + Vec oldHalfedgeTangent(std::move(halfedgeTangent_)); + Vec faceOld2New(oldHalfedge.size() / 3); auto policy = autoPolicy(numTri); scatter(policy, countAt(0), countAt(numTri), faceNew2Old.begin(), faceOld2New.begin()); @@ -435,13 +432,11 @@ void Manifold::Impl::GatherFaces(const VecDH& faceNew2Old) { halfedge_.resize(3 * numTri); if (oldHalfedgeTangent.size() != 0) halfedgeTangent_.resize(3 * numTri); for_each_n(policy, countAt(0), numTri, - ReindexFace({halfedge_.ptrD(), halfedgeTangent_.ptrD(), - oldHalfedge.cptrD(), oldHalfedgeTangent.cptrD(), - faceNew2Old.cptrD(), faceOld2New.cptrD()})); + ReindexFace({halfedge_, halfedgeTangent_, oldHalfedge, + oldHalfedgeTangent, faceNew2Old, faceOld2New})); } -void Manifold::Impl::GatherFaces(const Impl& old, - const VecDH& faceNew2Old) { +void Manifold::Impl::GatherFaces(const Impl& old, const Vec& faceNew2Old) { const int numTri = faceNew2Old.size(); auto policy = autoPolicy(numTri); @@ -467,16 +462,15 @@ void Manifold::Impl::GatherFaces(const Impl& old, old.faceNormal_.begin(), faceNormal_.begin()); } - VecDH faceOld2New(old.NumTri()); + Vec faceOld2New(old.NumTri()); scatter(policy, countAt(0), countAt(numTri), faceNew2Old.begin(), faceOld2New.begin()); halfedge_.resize(3 * numTri); if (old.halfedgeTangent_.size() != 0) halfedgeTangent_.resize(3 * numTri); for_each_n(policy, countAt(0), numTri, - ReindexFace({halfedge_.ptrD(), halfedgeTangent_.ptrD(), - old.halfedge_.cptrD(), old.halfedgeTangent_.cptrD(), - faceNew2Old.cptrD(), faceOld2New.cptrD()})); + ReindexFace({halfedge_, halfedgeTangent_, old.halfedge_, + old.halfedgeTangent_, faceNew2Old, faceOld2New})); } /// Constructs a position-only MeshGL from the input Mesh. @@ -542,18 +536,18 @@ bool MeshGL::Merge() { } const int numOpenVert = openEdges.size(); - VecDH openVerts(numOpenVert); + Vec openVerts(numOpenVert); int i = 0; for (const auto& edge : openEdges) { const int vert = edge.first; openVerts[i++] = vert; } - VecDH vertPropD(vertProperties); + Vec vertPropD(vertProperties); Box bBox; for (const int i : {0, 1, 2}) { - strided_range::Iter> iPos(vertPropD.begin() + i, - vertPropD.end(), numProp); + strided_range::Iter> iPos(vertPropD.begin() + i, vertPropD.end(), + numProp); auto minMax = transform_reduce>( autoPolicy(numVert), iPos.begin(), iPos.end(), Duplicate(), thrust::make_pair(std::numeric_limits::infinity(), @@ -566,19 +560,18 @@ bool MeshGL::Merge() { if (precision < 0) return false; auto policy = autoPolicy(numOpenVert); - VecDH vertBox(numOpenVert); - VecDH vertMorton(numOpenVert); + Vec vertBox(numOpenVert); + Vec vertMorton(numOpenVert); for_each_n(policy, zip(vertMorton.begin(), vertBox.begin(), openVerts.cbegin()), - numOpenVert, - VertMortonBox({vertPropD.cptrD(), numProp, precision, bBox})); + numOpenVert, VertMortonBox({vertPropD, numProp, precision, bBox})); sort_by_key(policy, vertMorton.begin(), vertMorton.end(), zip(vertBox.begin(), openVerts.begin())); Collider collider(vertBox, vertMorton); - SparseIndices toMerge = collider.Collisions(vertBox); + SparseIndices toMerge = collider.Collisions(vertBox.cview()); Graph graph; for (int i = 0; i < numVert; ++i) { diff --git a/src/sdf/CMakeLists.txt b/src/sdf/CMakeLists.txt index 4e0ad3a99..73dacc267 100644 --- a/src/sdf/CMakeLists.txt +++ b/src/sdf/CMakeLists.txt @@ -14,11 +14,23 @@ project(sdf) -add_library(${PROJECT_NAME} INTERFACE) +add_library(${PROJECT_NAME} src/sdf.cpp) target_include_directories(${PROJECT_NAME} - INTERFACE + PUBLIC ${PROJECT_SOURCE_DIR}/include ) -target_link_libraries(${PROJECT_NAME} INTERFACE utilities) +target_link_libraries(${PROJECT_NAME} PUBLIC utilities) + +install( + TARGETS ${PROJECT_NAME} + LIBRARY DESTINATION ${CMAKE_INSTALL_FULL_LIBDIR} + COMPONENT runtime +) +install( + DIRECTORY include/ + DESTINATION ${CMAKE_INSTALL_FULL_INCLUDEDIR}/manifold + COMPONENT devel +) + diff --git a/src/sdf/include/sdf.h b/src/sdf/include/sdf.h index 7d71d47b6..d65afe803 100644 --- a/src/sdf/include/sdf.h +++ b/src/sdf/include/sdf.h @@ -1,4 +1,4 @@ -// Copyright 2022 The Manifold Authors. +// Copyright 2023 The Manifold Authors. // // Licensed under the Apache License, Version 2.0 (the "License"); // you may not use this file except in compliance with the License. @@ -14,361 +14,11 @@ #pragma once -#include "hashtable.h" -#include "optional" -#include "par.h" -#include "public.h" -#include "utils.h" -#include "vec_dh.h" - -namespace { -using namespace manifold; -Uint64 identity(Uint64 x) { return x; } - -glm::ivec3 TetTri0(int i) { - constexpr glm::ivec3 tetTri0[16] = {{-1, -1, -1}, // - {0, 3, 4}, // - {0, 1, 5}, // - {1, 5, 3}, // - {1, 4, 2}, // - {1, 0, 3}, // - {2, 5, 0}, // - {5, 3, 2}, // - {2, 3, 5}, // - {0, 5, 2}, // - {3, 0, 1}, // - {2, 4, 1}, // - {3, 5, 1}, // - {5, 1, 0}, // - {4, 3, 0}, // - {-1, -1, -1}}; - return tetTri0[i]; -} - -glm::ivec3 TetTri1(int i) { - constexpr glm::ivec3 tetTri1[16] = {{-1, -1, -1}, // - {-1, -1, -1}, // - {-1, -1, -1}, // - {3, 4, 1}, // - {-1, -1, -1}, // - {3, 2, 1}, // - {0, 4, 2}, // - {-1, -1, -1}, // - {-1, -1, -1}, // - {2, 4, 0}, // - {1, 2, 3}, // - {-1, -1, -1}, // - {1, 4, 3}, // - {-1, -1, -1}, // - {-1, -1, -1}, // - {-1, -1, -1}}; - return tetTri1[i]; -} - -glm::ivec4 Neighbors(int i) { - constexpr glm::ivec4 neighbors[7] = {{0, 0, 0, 1}, // - {1, 0, 0, 0}, // - {0, 1, 0, 0}, // - {0, 0, 1, 0}, // - {-1, 0, 0, 1}, // - {0, -1, 0, 1}, // - {0, 0, -1, 1}}; - return neighbors[i]; -} - -Uint64 SpreadBits3(Uint64 v) { - v = v & 0x1fffff; - v = (v | v << 32) & 0x1f00000000ffff; - v = (v | v << 16) & 0x1f0000ff0000ff; - v = (v | v << 8) & 0x100f00f00f00f00f; - v = (v | v << 4) & 0x10c30c30c30c30c3; - v = (v | v << 2) & 0x1249249249249249; - return v; -} - -Uint64 SqueezeBits3(Uint64 v) { - v = v & 0x1249249249249249; - v = (v ^ v >> 2) & 0x10c30c30c30c30c3; - v = (v ^ v >> 4) & 0x100f00f00f00f00f; - v = (v ^ v >> 8) & 0x1f0000ff0000ff; - v = (v ^ v >> 16) & 0x1f00000000ffff; - v = (v ^ v >> 32) & 0x1fffff; - return v; -} - -// This is a modified 3D MortonCode, where the xyz code is shifted by one bit -// and the w bit is added as the least significant. This allows 21 bits per x, -// y, and z channel and 1 for w, filling the 64 bit total. -Uint64 MortonCode(const glm::ivec4& index) { - return static_cast(index.w) | (SpreadBits3(index.x) << 1) | - (SpreadBits3(index.y) << 2) | (SpreadBits3(index.z) << 3); -} - -glm::ivec4 DecodeMorton(Uint64 code) { - glm::ivec4 index; - index.x = SqueezeBits3(code >> 1); - index.y = SqueezeBits3(code >> 2); - index.z = SqueezeBits3(code >> 3); - index.w = code & 0x1u; - return index; -} - -struct GridVert { - float distance = NAN; - int edgeVerts[7] = {-1, -1, -1, -1, -1, -1, -1}; - - int Inside() const { return distance > 0 ? 1 : -1; } - - int NeighborInside(int i) const { - return Inside() * (edgeVerts[i] < 0 ? 1 : -1); - } -}; - -template -struct ComputeVerts { - glm::vec3* vertPos; - int* vertIndex; - HashTableD gridVerts; - const Func sdf; - const glm::vec3 origin; - const glm::ivec3 gridSize; - const glm::vec3 spacing; - const float level; - - inline glm::vec3 Position(glm::ivec4 gridIndex) const { - return origin + - spacing * (glm::vec3(gridIndex) + (gridIndex.w == 1 ? 0.0f : -0.5f)); - } - - inline float BoundedSDF(glm::ivec4 gridIndex) const { - const float d = sdf(Position(gridIndex)) - level; - - const glm::ivec3 xyz(gridIndex); - const bool onLowerBound = glm::any(glm::lessThanEqual(xyz, glm::ivec3(0))); - const bool onUpperBound = glm::any(glm::greaterThanEqual(xyz, gridSize)); - const bool onHalfBound = - gridIndex.w == 1 && glm::any(glm::greaterThanEqual(xyz, gridSize - 1)); - if (onLowerBound || onUpperBound || onHalfBound) return glm::min(d, 0.0f); - - return d; - } - - inline void operator()(Uint64 mortonCode) { - if (gridVerts.Full()) return; - - const glm::ivec4 gridIndex = DecodeMorton(mortonCode); - - if (glm::any(glm::greaterThan(glm::ivec3(gridIndex), gridSize))) return; - - const glm::vec3 position = Position(gridIndex); - - GridVert gridVert; - gridVert.distance = BoundedSDF(gridIndex); - - bool keep = false; - // These seven edges are uniquely owned by this gridVert; any of them - // which intersect the surface create a vert. - for (int i = 0; i < 7; ++i) { - glm::ivec4 neighborIndex = gridIndex + Neighbors(i); - if (neighborIndex.w == 2) { - neighborIndex += 1; - neighborIndex.w = 0; - } - const float val = BoundedSDF(neighborIndex); - if ((val > 0) == (gridVert.distance > 0)) continue; - keep = true; +#include - const int idx = AtomicAdd(*vertIndex, 1); - vertPos[idx] = - (val * position - gridVert.distance * Position(neighborIndex)) / - (val - gridVert.distance); - gridVert.edgeVerts[i] = idx; - } - - if (keep) gridVerts.Insert(mortonCode, gridVert); - } -}; - -struct BuildTris { - glm::ivec3* triVerts; - int* triIndex; - const HashTableD gridVerts; - - void CreateTri(const glm::ivec3& tri, const int edges[6]) { - if (tri[0] < 0) return; - int idx = AtomicAdd(*triIndex, 1); - triVerts[idx] = {edges[tri[0]], edges[tri[1]], edges[tri[2]]}; - } - - void CreateTris(const glm::ivec4& tet, const int edges[6]) { - const int i = (tet[0] > 0 ? 1 : 0) + (tet[1] > 0 ? 2 : 0) + - (tet[2] > 0 ? 4 : 0) + (tet[3] > 0 ? 8 : 0); - CreateTri(TetTri0(i), edges); - CreateTri(TetTri1(i), edges); - } - - void operator()(int idx) { - Uint64 basekey = gridVerts.KeyAt(idx); - if (basekey == kOpen) return; - - const GridVert& base = gridVerts.At(idx); - const glm::ivec4 baseIndex = DecodeMorton(basekey); - - glm::ivec4 leadIndex = baseIndex; - if (leadIndex.w == 0) - leadIndex.w = 1; - else { - leadIndex += 1; - leadIndex.w = 0; - } - - // This GridVert is in charge of the 6 tetrahedra surrounding its edge in - // the (1,1,1) direction (edge 0). - glm::ivec4 tet(base.NeighborInside(0), base.Inside(), -2, -2); - glm::ivec4 thisIndex = baseIndex; - thisIndex.x += 1; - - GridVert thisVert = gridVerts[MortonCode(thisIndex)]; - - tet[2] = base.NeighborInside(1); - for (const int i : {0, 1, 2}) { - thisIndex = leadIndex; - --thisIndex[Prev3(i)]; - // MortonCodes take unsigned input, so check for negatives, given the - // decrement. - GridVert nextVert = thisIndex[Prev3(i)] < 0 - ? GridVert() - : gridVerts[MortonCode(thisIndex)]; - tet[3] = base.NeighborInside(Prev3(i) + 4); - - const int edges1[6] = {base.edgeVerts[0], - base.edgeVerts[i + 1], - nextVert.edgeVerts[Next3(i) + 4], - nextVert.edgeVerts[Prev3(i) + 1], - thisVert.edgeVerts[i + 4], - base.edgeVerts[Prev3(i) + 4]}; - thisVert = nextVert; - CreateTris(tet, edges1); - - thisIndex = baseIndex; - ++thisIndex[Next3(i)]; - nextVert = gridVerts[MortonCode(thisIndex)]; - tet[2] = tet[3]; - tet[3] = base.NeighborInside(Next3(i) + 1); - - const int edges2[6] = {base.edgeVerts[0], - edges1[5], - thisVert.edgeVerts[i + 4], - nextVert.edgeVerts[Next3(i) + 4], - edges1[3], - base.edgeVerts[Next3(i) + 1]}; - thisVert = nextVert; - CreateTris(tet, edges2); - - tet[2] = tet[3]; - } - } -}; - -template -void for_each_n_wrapper(ExecutionPolicy policy, int morton, - ComputeVerts cv) { - for_each_n(policy, countAt(0), morton, cv); -} -template <> -void for_each_n_wrapper>( - ExecutionPolicy policy, int morton, - ComputeVerts> cv) { - for_each_n(policy, countAt(0), morton, cv); -} - -} // namespace +#include "public.h" namespace manifold { - -/** @addtogroup Core - * @{ - */ - -/** - * Constructs a level-set Mesh from the input Signed-Distance Function (SDF). - * This uses a form of Marching Tetrahedra (akin to Marching Cubes, but better - * for manifoldness). Instead of using a cubic grid, it uses a body-centered - * cubic grid (two shifted cubic grids). This means if your function's interior - * exceeds the given bounds, you will see a kind of egg-crate shape closing off - * the manifold, which is due to the underlying grid. - * - * @param sdf The signed-distance functor, containing this function signature: - * `float operator()(glm::vec3 point)`, which returns the - * signed distance of a given point in R^3. Positive values are inside, - * negative outside. - * @param bounds An axis-aligned box that defines the extent of the grid. - * @param edgeLength Approximate maximum edge length of the triangles in the - * final result. This affects grid spacing, and hence has a strong effect on - * performance. - * @param level You can inset your Mesh by using a positive value, or outset - * it with a negative value. - * @return Mesh This class does not depend on Manifold, so it just returns a - * Mesh, but it is guaranteed to be manifold and so can always be used as - * input to the Manifold constructor for further operations. - */ -template -inline Mesh LevelSet(Func sdf, Box bounds, float edgeLength, float level = 0, - std::optional policy = std::nullopt) { - Mesh out; - - const glm::vec3 dim = bounds.Size(); - const float maxDim = std::max(dim[0], std::max(dim[1], dim[2])); - const glm::ivec3 gridSize(dim / edgeLength); - const glm::vec3 spacing = dim / (glm::vec3(gridSize)); - - const Uint64 maxMorton = MortonCode(glm::ivec4(gridSize + 1, 1)); - - // Parallel policies violate will crash language runtimes with runtime locks - // that expect to not be called back by unregistered threads. This allows - // bindings use LevelSet despite being compiled with MANIFOLD_PAR - // active. - const auto pol = !policy.has_value() ? autoPolicy(maxMorton) : policy.value(); - - int tableSize = glm::min( - 2 * maxMorton, static_cast(10 * glm::pow(maxMorton, 0.667))); - HashTable gridVerts(tableSize); - VecDH vertPos(gridVerts.Size() * 7); - - while (1) { - VecDH index(1, 0); - for_each_n_wrapper( - pol, maxMorton + 1, - ComputeVerts({vertPos.ptrD(), index.ptrD(), gridVerts.D(), sdf, - bounds.min, gridSize + 1, spacing, level})); - - if (gridVerts.Full()) { // Resize HashTable - const glm::vec3 lastVert = vertPos[index[0] - 1]; - const Uint64 lastMorton = - MortonCode(glm::ivec4((lastVert - bounds.min) / spacing, 1)); - const float ratio = static_cast(maxMorton) / lastMorton; - if (ratio > 1000) // do not trust the ratio if it is too large - tableSize *= 2; - else - tableSize *= ratio; - gridVerts = HashTable(tableSize); - vertPos = VecDH(gridVerts.Size() * 7); - } else { // Success - vertPos.resize(index[0]); - break; - } - } - - VecDH triVerts(gridVerts.Entries() * 12); // worst case - - VecDH index(1, 0); - for_each_n(pol, countAt(0), gridVerts.Size(), - BuildTris({triVerts.ptrD(), index.ptrD(), gridVerts.D()})); - triVerts.resize(index[0]); - - out.vertPos.insert(out.vertPos.end(), vertPos.begin(), vertPos.end()); - out.triVerts.insert(out.triVerts.end(), triVerts.begin(), triVerts.end()); - return out; +Mesh LevelSet(std::function sdf, Box bounds, float edgeLength, + float level = 0, bool canParallel = true); } -/** @} */ -} // namespace manifold diff --git a/src/sdf/src/sdf.cpp b/src/sdf/src/sdf.cpp new file mode 100644 index 000000000..4acb0774e --- /dev/null +++ b/src/sdf/src/sdf.cpp @@ -0,0 +1,360 @@ +// Copyright 2023 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "sdf.h" + +#include "hashtable.h" +#include "par.h" +#include "utils.h" +#include "vec.h" + +namespace { +using namespace manifold; +Uint64 identity(Uint64 x) { return x; } + +glm::ivec3 TetTri0(int i) { + constexpr glm::ivec3 tetTri0[16] = {{-1, -1, -1}, // + {0, 3, 4}, // + {0, 1, 5}, // + {1, 5, 3}, // + {1, 4, 2}, // + {1, 0, 3}, // + {2, 5, 0}, // + {5, 3, 2}, // + {2, 3, 5}, // + {0, 5, 2}, // + {3, 0, 1}, // + {2, 4, 1}, // + {3, 5, 1}, // + {5, 1, 0}, // + {4, 3, 0}, // + {-1, -1, -1}}; + return tetTri0[i]; +} + +glm::ivec3 TetTri1(int i) { + constexpr glm::ivec3 tetTri1[16] = {{-1, -1, -1}, // + {-1, -1, -1}, // + {-1, -1, -1}, // + {3, 4, 1}, // + {-1, -1, -1}, // + {3, 2, 1}, // + {0, 4, 2}, // + {-1, -1, -1}, // + {-1, -1, -1}, // + {2, 4, 0}, // + {1, 2, 3}, // + {-1, -1, -1}, // + {1, 4, 3}, // + {-1, -1, -1}, // + {-1, -1, -1}, // + {-1, -1, -1}}; + return tetTri1[i]; +} + +glm::ivec4 Neighbors(int i) { + constexpr glm::ivec4 neighbors[7] = {{0, 0, 0, 1}, // + {1, 0, 0, 0}, // + {0, 1, 0, 0}, // + {0, 0, 1, 0}, // + {-1, 0, 0, 1}, // + {0, -1, 0, 1}, // + {0, 0, -1, 1}}; + return neighbors[i]; +} + +Uint64 SpreadBits3(Uint64 v) { + v = v & 0x1fffff; + v = (v | v << 32) & 0x1f00000000ffff; + v = (v | v << 16) & 0x1f0000ff0000ff; + v = (v | v << 8) & 0x100f00f00f00f00f; + v = (v | v << 4) & 0x10c30c30c30c30c3; + v = (v | v << 2) & 0x1249249249249249; + return v; +} + +Uint64 SqueezeBits3(Uint64 v) { + v = v & 0x1249249249249249; + v = (v ^ v >> 2) & 0x10c30c30c30c30c3; + v = (v ^ v >> 4) & 0x100f00f00f00f00f; + v = (v ^ v >> 8) & 0x1f0000ff0000ff; + v = (v ^ v >> 16) & 0x1f00000000ffff; + v = (v ^ v >> 32) & 0x1fffff; + return v; +} + +// This is a modified 3D MortonCode, where the xyz code is shifted by one bit +// and the w bit is added as the least significant. This allows 21 bits per x, +// y, and z channel and 1 for w, filling the 64 bit total. +Uint64 MortonCode(const glm::ivec4& index) { + return static_cast(index.w) | (SpreadBits3(index.x) << 1) | + (SpreadBits3(index.y) << 2) | (SpreadBits3(index.z) << 3); +} + +glm::ivec4 DecodeMorton(Uint64 code) { + glm::ivec4 index; + index.x = SqueezeBits3(code >> 1); + index.y = SqueezeBits3(code >> 2); + index.z = SqueezeBits3(code >> 3); + index.w = code & 0x1u; + return index; +} + +struct GridVert { + float distance = NAN; + int edgeVerts[7] = {-1, -1, -1, -1, -1, -1, -1}; + + int Inside() const { return distance > 0 ? 1 : -1; } + + int NeighborInside(int i) const { + return Inside() * (edgeVerts[i] < 0 ? 1 : -1); + } +}; + +struct ComputeVerts { + VecView vertPos; + VecView vertIndex; + HashTableD gridVerts; + const std::function sdf; + const glm::vec3 origin; + const glm::ivec3 gridSize; + const glm::vec3 spacing; + const float level; + + inline glm::vec3 Position(glm::ivec4 gridIndex) const { + return origin + + spacing * (glm::vec3(gridIndex) + (gridIndex.w == 1 ? 0.0f : -0.5f)); + } + + inline float BoundedSDF(glm::ivec4 gridIndex) const { + const float d = sdf(Position(gridIndex)) - level; + + const glm::ivec3 xyz(gridIndex); + const bool onLowerBound = glm::any(glm::lessThanEqual(xyz, glm::ivec3(0))); + const bool onUpperBound = glm::any(glm::greaterThanEqual(xyz, gridSize)); + const bool onHalfBound = + gridIndex.w == 1 && glm::any(glm::greaterThanEqual(xyz, gridSize - 1)); + if (onLowerBound || onUpperBound || onHalfBound) return glm::min(d, 0.0f); + + return d; + } + + inline void operator()(Uint64 mortonCode) { + if (gridVerts.Full()) return; + + const glm::ivec4 gridIndex = DecodeMorton(mortonCode); + + if (glm::any(glm::greaterThan(glm::ivec3(gridIndex), gridSize))) return; + + const glm::vec3 position = Position(gridIndex); + + GridVert gridVert; + gridVert.distance = BoundedSDF(gridIndex); + + bool keep = false; + // These seven edges are uniquely owned by this gridVert; any of them + // which intersect the surface create a vert. + for (int i = 0; i < 7; ++i) { + glm::ivec4 neighborIndex = gridIndex + Neighbors(i); + if (neighborIndex.w == 2) { + neighborIndex += 1; + neighborIndex.w = 0; + } + const float val = BoundedSDF(neighborIndex); + if ((val > 0) == (gridVert.distance > 0)) continue; + keep = true; + + const int idx = AtomicAdd(vertIndex[0], 1); + vertPos[idx] = + (val * position - gridVert.distance * Position(neighborIndex)) / + (val - gridVert.distance); + gridVert.edgeVerts[i] = idx; + } + + if (keep) gridVerts.Insert(mortonCode, gridVert); + } +}; + +struct BuildTris { + VecView triVerts; + VecView triIndex; + const HashTableD gridVerts; + + void CreateTri(const glm::ivec3& tri, const int edges[6]) { + if (tri[0] < 0) return; + int idx = AtomicAdd(triIndex[0], 1); + triVerts[idx] = {edges[tri[0]], edges[tri[1]], edges[tri[2]]}; + } + + void CreateTris(const glm::ivec4& tet, const int edges[6]) { + const int i = (tet[0] > 0 ? 1 : 0) + (tet[1] > 0 ? 2 : 0) + + (tet[2] > 0 ? 4 : 0) + (tet[3] > 0 ? 8 : 0); + CreateTri(TetTri0(i), edges); + CreateTri(TetTri1(i), edges); + } + + void operator()(int idx) { + Uint64 basekey = gridVerts.KeyAt(idx); + if (basekey == kOpen) return; + + const GridVert& base = gridVerts.At(idx); + const glm::ivec4 baseIndex = DecodeMorton(basekey); + + glm::ivec4 leadIndex = baseIndex; + if (leadIndex.w == 0) + leadIndex.w = 1; + else { + leadIndex += 1; + leadIndex.w = 0; + } + + // This GridVert is in charge of the 6 tetrahedra surrounding its edge in + // the (1,1,1) direction (edge 0). + glm::ivec4 tet(base.NeighborInside(0), base.Inside(), -2, -2); + glm::ivec4 thisIndex = baseIndex; + thisIndex.x += 1; + + GridVert thisVert = gridVerts[MortonCode(thisIndex)]; + + tet[2] = base.NeighborInside(1); + for (const int i : {0, 1, 2}) { + thisIndex = leadIndex; + --thisIndex[Prev3(i)]; + // MortonCodes take unsigned input, so check for negatives, given the + // decrement. + GridVert nextVert = thisIndex[Prev3(i)] < 0 + ? GridVert() + : gridVerts[MortonCode(thisIndex)]; + tet[3] = base.NeighborInside(Prev3(i) + 4); + + const int edges1[6] = {base.edgeVerts[0], + base.edgeVerts[i + 1], + nextVert.edgeVerts[Next3(i) + 4], + nextVert.edgeVerts[Prev3(i) + 1], + thisVert.edgeVerts[i + 4], + base.edgeVerts[Prev3(i) + 4]}; + thisVert = nextVert; + CreateTris(tet, edges1); + + thisIndex = baseIndex; + ++thisIndex[Next3(i)]; + nextVert = gridVerts[MortonCode(thisIndex)]; + tet[2] = tet[3]; + tet[3] = base.NeighborInside(Next3(i) + 1); + + const int edges2[6] = {base.edgeVerts[0], + edges1[5], + thisVert.edgeVerts[i + 4], + nextVert.edgeVerts[Next3(i) + 4], + edges1[3], + base.edgeVerts[Next3(i) + 1]}; + thisVert = nextVert; + CreateTris(tet, edges2); + + tet[2] = tet[3]; + } + } +}; +} // namespace + +namespace manifold { + +/** @addtogroup Core + * @{ + */ + +/** + * Constructs a level-set Mesh from the input Signed-Distance Function (SDF). + * This uses a form of Marching Tetrahedra (akin to Marching Cubes, but better + * for manifoldness). Instead of using a cubic grid, it uses a body-centered + * cubic grid (two shifted cubic grids). This means if your function's interior + * exceeds the given bounds, you will see a kind of egg-crate shape closing off + * the manifold, which is due to the underlying grid. + * + * @param sdf The signed-distance functor, containing this function signature: + * `float operator()(glm::vec3 point)`, which returns the + * signed distance of a given point in R^3. Positive values are inside, + * negative outside. + * @param bounds An axis-aligned box that defines the extent of the grid. + * @param edgeLength Approximate maximum edge length of the triangles in the + * final result. This affects grid spacing, and hence has a strong effect on + * performance. + * @param level You can inset your Mesh by using a positive value, or outset + * it with a negative value. + * @param canParallel Parallel policies violate will crash language runtimes + * with runtime locks that expect to not be called back by unregistered threads. + * This allows bindings use LevelSet despite being compiled with MANIFOLD_PAR + * active. + * @return Mesh This class does not depend on Manifold, so it just returns a + * Mesh, but it is guaranteed to be manifold and so can always be used as + * input to the Manifold constructor for further operations. + */ +Mesh LevelSet(std::function sdf, Box bounds, float edgeLength, + float level, bool canParallel) { + Mesh out; + + const glm::vec3 dim = bounds.Size(); + const float maxDim = std::max(dim[0], std::max(dim[1], dim[2])); + const glm::ivec3 gridSize(dim / edgeLength); + const glm::vec3 spacing = dim / (glm::vec3(gridSize)); + + const Uint64 maxMorton = MortonCode(glm::ivec4(gridSize + 1, 1)); + + // Parallel policies violate will crash language runtimes with runtime locks + // that expect to not be called back by unregistered threads. This allows + // bindings use LevelSet despite being compiled with MANIFOLD_PAR + // active. + const auto pol = canParallel ? autoPolicy(maxMorton) : ExecutionPolicy::Seq; + + int tableSize = glm::min( + 2 * maxMorton, static_cast(10 * glm::pow(maxMorton, 0.667))); + HashTable gridVerts(tableSize); + Vec vertPos(gridVerts.Size() * 7); + + while (1) { + Vec index(1, 0); + for_each_n(pol, countAt(0), maxMorton + 1, + ComputeVerts({vertPos, index, gridVerts.D(), sdf, bounds.min, + gridSize + 1, spacing, level})); + + if (gridVerts.Full()) { // Resize HashTable + const glm::vec3 lastVert = vertPos[index[0] - 1]; + const Uint64 lastMorton = + MortonCode(glm::ivec4((lastVert - bounds.min) / spacing, 1)); + const float ratio = static_cast(maxMorton) / lastMorton; + if (ratio > 1000) // do not trust the ratio if it is too large + tableSize *= 2; + else + tableSize *= ratio; + gridVerts = HashTable(tableSize); + vertPos = Vec(gridVerts.Size() * 7); + } else { // Success + vertPos.resize(index[0]); + break; + } + } + + Vec triVerts(gridVerts.Entries() * 12); // worst case + + Vec index(1, 0); + for_each_n(pol, countAt(0), gridVerts.Size(), + BuildTris({triVerts, index, gridVerts.D()})); + triVerts.resize(index[0]); + + out.vertPos.insert(out.vertPos.end(), vertPos.begin(), vertPos.end()); + out.triVerts.insert(out.triVerts.end(), triVerts.begin(), triVerts.end()); + return out; +} +/** @} */ +} // namespace manifold diff --git a/src/utilities/CMakeLists.txt b/src/utilities/CMakeLists.txt index 36c870f82..74807babf 100644 --- a/src/utilities/CMakeLists.txt +++ b/src/utilities/CMakeLists.txt @@ -80,3 +80,10 @@ if(MANIFOLD_DEBUG) endif() target_compile_features(${PROJECT_NAME} INTERFACE cxx_std_17) + +install( + FILES include/public.h + DESTINATION ${CMAKE_INSTALL_FULL_INCLUDEDIR}/manifold + COMPONENT devel +) + diff --git a/src/utilities/include/hashtable.h b/src/utilities/include/hashtable.h index ca80a50ff..7a1f28ccb 100644 --- a/src/utilities/include/hashtable.h +++ b/src/utilities/include/hashtable.h @@ -18,12 +18,12 @@ #include "public.h" #include "utils.h" -#include "vec_dh.h" +#include "vec.h" namespace { typedef unsigned long long int Uint64; typedef Uint64 (*hash_fun_t)(Uint64); -constexpr Uint64 kOpen = std::numeric_limits::max(); +inline constexpr Uint64 kOpen = std::numeric_limits::max(); template T AtomicCAS(T& target, T compare, T val) { @@ -63,7 +63,7 @@ namespace manifold { template class HashTableD { public: - HashTableD(VecDH& keys, VecDH& values, VecDH& used, + HashTableD(Vec& keys, Vec& values, Vec& used, uint32_t step = 1) : step_{step}, keys_{keys}, values_{values}, used_{used} {} @@ -87,7 +87,18 @@ class HashTableD { } } - V& operator[](Uint64 key) const { + V& operator[](Uint64 key) { + uint32_t idx = H(key) & (Size() - 1); + while (1) { + const Uint64 k = AtomicLoad(keys_[idx]); + if (k == key || k == kOpen) { + return values_[idx]; + } + idx = (idx + step_) & (Size() - 1); + } + } + + const V& operator[](Uint64 key) const { uint32_t idx = H(key) & (Size() - 1); while (1) { const Uint64 k = AtomicLoad(keys_[idx]); @@ -99,13 +110,14 @@ class HashTableD { } Uint64 KeyAt(int idx) const { return AtomicLoad(keys_[idx]); } - V& At(int idx) const { return values_[idx]; } + V& At(int idx) { return values_[idx]; } + const V& At(int idx) const { return values_[idx]; } private: uint32_t step_; - VecD keys_; - VecD values_; - VecD used_; + VecView keys_; + VecView values_; + VecView used_; }; template @@ -128,14 +140,14 @@ class HashTable { return static_cast(AtomicLoad(used_[0])) / Size(); } - VecDH& GetValueStore() { return values_; } + Vec& GetValueStore() { return values_; } static Uint64 Open() { return kOpen; } private: - VecDH keys_; - VecDH values_; - VecDH used_ = VecDH(1, 0); + Vec keys_; + Vec values_; + Vec used_ = Vec(1, 0); HashTableD table_; }; diff --git a/src/utilities/include/public.h b/src/utilities/include/public.h index 31a66a985..06a76d071 100644 --- a/src/utilities/include/public.h +++ b/src/utilities/include/public.h @@ -401,7 +401,7 @@ class Quality { int nSegL = 2.0f * radius * glm::pi() / circularEdgeLength_; int nSeg = fmin(nSegA, nSegL) + 3; nSeg -= nSeg % 4; - return nSeg; + return std::max(nSeg, 3); } }; /** @} */ diff --git a/src/utilities/include/sparse.h b/src/utilities/include/sparse.h index 9291f11df..b43fc343e 100644 --- a/src/utilities/include/sparse.h +++ b/src/utilities/include/sparse.h @@ -19,7 +19,7 @@ #include "par.h" #include "public.h" #include "utils.h" -#include "vec_dh.h" +#include "vec.h" namespace manifold { @@ -52,8 +52,8 @@ class SparseIndices { int size() const { return data.size(); } - VecDH Copy(bool use_q) const { - VecDH out(data.size()); + Vec Copy(bool use_q) const { + Vec out(data.size()); int offset = pOffset; if (use_q) offset = 1 - offset; const int* p = ptr(); @@ -86,7 +86,7 @@ class SparseIndices { inline void SetPQ(int i, int64_t pq) { data[i] = pq; } - const VecDH& AsVec64() const { return data; } + const Vec& AsVec64() const { return data; } inline void Add(int p, int q) { data.push_back(EncodePQ(p, q)); } @@ -98,7 +98,7 @@ class SparseIndices { Resize(newSize); } - size_t RemoveZeros(VecDH& S) { + size_t RemoveZeros(Vec& S) { ASSERT(S.size() == data.size(), userErr, "Different number of values than indicies!"); auto zBegin = zip(S.begin(), data.begin()); @@ -128,7 +128,7 @@ class SparseIndices { }; template - size_t KeepFinite(VecDH& v, VecDH& x) { + size_t KeepFinite(Vec& v, Vec& x) { ASSERT(x.size() == data.size(), userErr, "Different number of values than indicies!"); auto zBegin = zip(v.begin(), x.begin(), data.begin()); @@ -155,10 +155,10 @@ class SparseIndices { #endif private: - VecDH data; - inline int* ptr() { return reinterpret_cast(data.ptrD()); } + Vec data; + inline int* ptr() { return reinterpret_cast(data.data()); } inline const int* ptr() const { - return reinterpret_cast(data.ptrD()); + return reinterpret_cast(data.data()); } }; diff --git a/src/utilities/include/vec.h b/src/utilities/include/vec.h new file mode 100644 index 000000000..a48672bce --- /dev/null +++ b/src/utilities/include/vec.h @@ -0,0 +1,320 @@ +// Copyright 2021 The Manifold Authors. +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#pragma once +#include +#if TRACY_ENABLE && TRACY_MEMORY_USAGE +#include "tracy/Tracy.hpp" +#else +#define TracyAllocS(ptr, size, n) (void)0 +#define TracyFreeS(ptr, n) (void)0 +#endif + +// #include "optional_assert.h" +#include "par.h" +#include "public.h" + +namespace manifold { + +/** @addtogroup Private + * @{ + */ +template +class Vec; + +/** + * View for Vec, can perform offset operation. + * This will be invalidated when the original vector is dropped or changes + * length. + */ +template +class VecView { + public: + using Iter = T *; + using IterC = const T *; + + VecView(const VecView &other) { + ptr_ = other.ptr_; + size_ = other.size_; + } + + VecView &operator=(const VecView &other) { + ptr_ = other.ptr_; + size_ = other.size_; + return *this; + } + + // allows conversion to a const VecView + operator VecView() const { return {ptr_, size_}; } + + inline const T &operator[](int i) const { + if (i < 0 || i >= size_) throw std::out_of_range("Vec out of range"); + return ptr_[i]; + } + + inline T &operator[](int i) { + if (i < 0 || i >= size_) throw std::out_of_range("Vec out of range"); + return ptr_[i]; + } + + IterC cbegin() const { return ptr_; } + IterC cend() const { return ptr_ + size_; } + + IterC begin() const { return cbegin(); } + IterC end() const { return cend(); } + + Iter begin() { return ptr_; } + Iter end() { return ptr_ + size_; } + + const T &front() const { + if (size_ == 0) + throw std::out_of_range("attempt to take the front of an empty vector"); + return ptr_[0]; + } + + const T &back() const { + if (size_ == 0) + throw std::out_of_range("attempt to take the back of an empty vector"); + return ptr_[size_ - 1]; + } + + T &front() { + if (size_ == 0) + throw std::out_of_range("attempt to take the front of an empty vector"); + return ptr_[0]; + } + + T &back() { + if (size_ == 0) + throw std::out_of_range("attempt to take the back of an empty vector"); + return ptr_[size_ - 1]; + } + + int size() const { return size_; } + + bool empty() const { return size_ == 0; } + + protected: + T *ptr_ = nullptr; + int size_ = 0; + + VecView() = default; + VecView(T *ptr_, int size_) : ptr_(ptr_), size_(size_) {} + friend class Vec; + friend class Vec::type>; + friend class VecView::type>; +}; + +/* + * Specialized vector implementation with multithreaded fill and uninitialized + * memory optimizations. + * Note that the constructor and resize function will not perform initialization + * if the parameter val is not set. Also, this implementation is a toy + * implementation that did not consider things like non-trivial + * constructor/destructor, please keep T trivial. + */ +template +class Vec : public VecView { + public: + Vec() {} + + // Note that the vector constructed with this constructor will contain + // uninitialized memory. Please specify `val` if you need to make sure that + // the data is initialized. + Vec(int size) { + reserve(size); + this->size_ = size; + } + + Vec(int size, T val) { resize(size, val); } + + Vec(const Vec &vec) { + this->size_ = vec.size(); + this->capacity_ = this->size_; + auto policy = autoPolicy(this->size_); + if (this->size_ != 0) { + this->ptr_ = reinterpret_cast(malloc(this->size_ * sizeof(T))); + if (this->ptr_ == nullptr) throw std::bad_alloc(); + TracyAllocS(ptr_, size_ * sizeof(T), 3); + uninitialized_copy(policy, vec.begin(), vec.end(), this->ptr_); + } + } + + Vec(const std::vector &vec) { + this->size_ = vec.size(); + this->capacity_ = this->size_; + auto policy = autoPolicy(this->size_); + if (this->size_ != 0) { + this->ptr_ = reinterpret_cast(malloc(this->size_ * sizeof(T))); + if (this->ptr_ == nullptr) throw std::bad_alloc(); + TracyAllocS(ptr_, size_ * sizeof(T), 3); + uninitialized_copy(policy, vec.begin(), vec.end(), this->ptr_); + } + } + + Vec(Vec &&vec) { + this->ptr_ = vec.ptr_; + this->size_ = vec.size_; + capacity_ = vec.capacity_; + vec.ptr_ = nullptr; + vec.size_ = 0; + vec.capacity_ = 0; + } + + operator VecView() { return {this->ptr_, this->size_}; } + operator VecView() const { return {this->ptr_, this->size_}; } + + ~Vec() { + if (this->ptr_ != nullptr) { + TracyFreeS(ptr_, 3); + free(this->ptr_); + } + this->ptr_ = nullptr; + this->size_ = 0; + capacity_ = 0; + } + + Vec &operator=(const Vec &other) { + if (&other == this) return *this; + if (this->ptr_ != nullptr) { + TracyFreeS(this->ptr_, 3); + free(this->ptr_); + } + this->size_ = other.size_; + capacity_ = other.size_; + auto policy = autoPolicy(this->size_); + if (this->size_ != 0) { + this->ptr_ = reinterpret_cast(malloc(this->size_ * sizeof(T))); + if (this->ptr_ == nullptr) throw std::bad_alloc(); + TracyAllocS(ptr_, size_ * sizeof(T), 3); + uninitialized_copy(policy, other.begin(), other.end(), this->ptr_); + } + return *this; + } + + Vec &operator=(Vec &&other) { + if (&other == this) return *this; + if (this->ptr_ != nullptr) { + TracyFreeS(ptr_, 3); + free(this->ptr_); + } + this->size_ = other.size_; + capacity_ = other.capacity_; + this->ptr_ = other.ptr_; + other.ptr_ = nullptr; + other.size_ = 0; + other.capacity_ = 0; + return *this; + } + + operator VecView() const { return {this->ptr_, this->size_}; } + + void swap(Vec &other) { + std::swap(this->ptr_, other.ptr_); + std::swap(this->size_, other.size_); + std::swap(capacity_, other.capacity_); + } + + inline void push_back(const T &val) { + if (this->size_ >= capacity_) { + // avoid dangling pointer in case val is a reference of our array + T val_copy = val; + reserve(capacity_ == 0 ? 128 : capacity_ * 2); + this->ptr_[this->size_++] = val_copy; + return; + } + this->ptr_[this->size_++] = val; + } + + void reserve(int n) { + if (n > capacity_) { + T *newBuffer = reinterpret_cast(malloc(n * sizeof(T))); + if (newBuffer == nullptr) throw std::bad_alloc(); + TracyAllocS(newBuffer, n * sizeof(T), 3); + if (this->size_ > 0) + uninitialized_copy(autoPolicy(this->size_), this->ptr_, + this->ptr_ + this->size_, newBuffer); + if (this->ptr_ != nullptr) { + TracyFreeS(ptr_, 3); + free(this->ptr_); + } + this->ptr_ = newBuffer; + capacity_ = n; + } + } + + void resize(int newSize, T val = T()) { + bool shrink = this->size_ > 2 * newSize; + reserve(newSize); + if (this->size_ < newSize) { + uninitialized_fill(autoPolicy(newSize - this->size_), + this->ptr_ + this->size_, this->ptr_ + newSize, val); + } + this->size_ = newSize; + if (shrink) shrink_to_fit(); + } + + void shrink_to_fit() { + T *newBuffer = nullptr; + if (this->size_ > 0) { + newBuffer = reinterpret_cast(malloc(this->size_ * sizeof(T))); + if (newBuffer == nullptr) throw std::bad_alloc(); + TracyAllocS(newBuffer, size_ * sizeof(T), 3); + uninitialized_copy(autoPolicy(this->size_), this->ptr_, + this->ptr_ + this->size_, newBuffer); + if (this->ptr_ != nullptr) { + TracyFreeS(ptr_, 3); + free(this->ptr_); + } + } + this->ptr_ = newBuffer; + capacity_ = this->size_; + } + + VecView view(int offset = 0, int length = -1) { + if (length == -1) { + length = this->size_ - offset; + if (length < 0) throw std::out_of_range("Vec::view out of range"); + } else if (offset + length > this->size_ || offset < 0) { + throw std::out_of_range("Vec::view out of range"); + } else if (length < 0) { + throw std::out_of_range("Vec::view negative length is not allowed"); + } + return VecView(this->ptr_ + offset, length); + } + + VecView cview(int offset = 0, int length = -1) const { + if (length == -1) { + length = this->size_ - offset; + if (length < 0) throw std::out_of_range("Vec::cview out of range"); + } else if (offset + length > this->size_ || offset < 0) { + throw std::out_of_range("Vec::cview out of range"); + } else if (length < 0) { + throw std::out_of_range("Vec::cview negative length is not allowed"); + } + return VecView(this->ptr_ + offset, length); + } + + VecView view(int offset = 0, int length = -1) const { + return cview(offset, length); + } + + T *data() { return this->ptr_; } + const T *data() const { return this->ptr_; } + + private: + int capacity_ = 0; +}; +/** @} */ +} // namespace manifold diff --git a/src/utilities/include/vec_dh.h b/src/utilities/include/vec_dh.h deleted file mode 100644 index 7d86633ba..000000000 --- a/src/utilities/include/vec_dh.h +++ /dev/null @@ -1,268 +0,0 @@ -// Copyright 2021 The Manifold Authors. -// -// Licensed under the Apache License, Version 2.0 (the "License"); -// you may not use this file except in compliance with the License. -// You may obtain a copy of the License at -// -// http://www.apache.org/licenses/LICENSE-2.0 -// -// Unless required by applicable law or agreed to in writing, software -// distributed under the License is distributed on an "AS IS" BASIS, -// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. -// See the License for the specific language governing permissions and -// limitations under the License. - -#pragma once -#include -#if TRACY_ENABLE && TRACY_MEMORY_USAGE -#include -#else -#define TracyAllocS(ptr, size, n) (void)0 -#define TracyFreeS(ptr, n) (void)0 -#endif - -// #include "optional_assert.h" -#include "par.h" -#include "public.h" - -namespace manifold { - -/** @addtogroup Private - * @{ - */ - -/* - * Specialized vector implementation with multithreaded fill and uninitialized - * memory optimizations. - * Note that the constructor and resize function will not perform initialization - * if the parameter val is not set. Also, this implementation is a toy - * implementation that did not consider things like non-trivial - * constructor/destructor, please keep T trivial. - */ -template -class VecDH { - public: - using Iter = T *; - using IterC = const T *; - - VecDH() {} - - // Note that the vector constructed with this constructor will contain - // uninitialized memory. Please specify `val` if you need to make sure that - // the data is initialized. - VecDH(int size) { resize(size); } - - VecDH(int size, T val) { resize(size, val); } - - VecDH(const std::vector &vec) { - size_ = vec.size(); - capacity_ = size_; - auto policy = autoPolicy(size_); - if (size_ != 0) { - ptr_ = reinterpret_cast(malloc(size_ * sizeof(T))); - if (ptr_ == nullptr) throw std::bad_alloc(); - TracyAllocS(ptr_, size_ * sizeof(T), 3); - uninitialized_copy(policy, vec.begin(), vec.end(), ptr_); - } - } - - VecDH(const VecDH &vec) { - size_ = vec.size_; - capacity_ = size_; - auto policy = autoPolicy(size_); - if (size_ != 0) { - ptr_ = reinterpret_cast(malloc(size_ * sizeof(T))); - if (ptr_ == nullptr) throw std::bad_alloc(); - TracyAllocS(ptr_, size_ * sizeof(T), 3); - uninitialized_copy(policy, vec.begin(), vec.end(), ptr_); - } - } - - VecDH(VecDH &&vec) { - ptr_ = vec.ptr_; - size_ = vec.size_; - capacity_ = vec.capacity_; - vec.ptr_ = nullptr; - vec.size_ = 0; - vec.capacity_ = 0; - } - - ~VecDH() { - if (ptr_ != nullptr) { - TracyFreeS(ptr_, 3); - free(ptr_); - } - ptr_ = nullptr; - size_ = 0; - capacity_ = 0; - } - - VecDH &operator=(const VecDH &other) { - if (&other == this) return *this; - if (ptr_ != nullptr) { - TracyFreeS(ptr_, 3); - free(ptr_); - } - size_ = other.size_; - capacity_ = other.size_; - auto policy = autoPolicy(size_); - if (size_ != 0) { - ptr_ = reinterpret_cast(malloc(size_ * sizeof(T))); - if (ptr_ == nullptr) throw std::bad_alloc(); - TracyAllocS(ptr_, size_ * sizeof(T), 3); - uninitialized_copy(policy, other.begin(), other.end(), ptr_); - } - return *this; - } - - VecDH &operator=(VecDH &&other) { - if (&other == this) return *this; - if (ptr_ != nullptr) { - TracyFreeS(ptr_, 3); - free(ptr_); - } - size_ = other.size_; - capacity_ = other.capacity_; - ptr_ = other.ptr_; - other.ptr_ = nullptr; - other.size_ = 0; - other.capacity_ = 0; - return *this; - } - - int size() const { return size_; } - - void swap(VecDH &other) { - std::swap(ptr_, other.ptr_); - std::swap(size_, other.size_); - std::swap(capacity_, other.capacity_); - } - - Iter begin() { return ptr_; } - - Iter end() { return ptr_ + size_; } - - IterC cbegin() const { return ptr_; } - - IterC cend() const { return ptr_ + size_; } - - IterC begin() const { return cbegin(); } - IterC end() const { return cend(); } - - T *ptrD() { return ptr_; } - - const T *cptrD() const { return ptr_; } - - const T *ptrD() const { return cptrD(); } - - T *ptrH() { return ptr_; } - - const T *cptrH() const { return ptr_; } - - const T *ptrH() const { return cptrH(); } - - inline T &operator[](int i) { return ptr_[i]; } - - inline const T &operator[](int i) const { return ptr_[i]; } - - T &back() { return ptr_[size_ - 1]; } - - const T &back() const { return ptr_[size_ - 1]; } - - inline void push_back(const T &val) { - if (size_ >= capacity_) { - // avoid dangling pointer in case val is a reference of our array - T val_copy = val; - reserve(capacity_ == 0 ? 128 : capacity_ * 2); - ptr_[size_++] = val_copy; - return; - } - ptr_[size_++] = val; - } - - void reserve(int n) { - if (n > capacity_) { - T *newBuffer = reinterpret_cast(malloc(n * sizeof(T))); - if (newBuffer == nullptr) throw std::bad_alloc(); - TracyAllocS(newBuffer, n * sizeof(T), 3); - if (size_ > 0) - uninitialized_copy(autoPolicy(size_), ptr_, ptr_ + size_, newBuffer); - if (ptr_ != nullptr) { - TracyFreeS(ptr_, 3); - free(ptr_); - } - ptr_ = newBuffer; - capacity_ = n; - } - } - - void resize(int newSize, T val = T()) { - bool shrink = size_ > 2 * newSize; - reserve(newSize); - if (size_ < newSize) { - uninitialized_fill(autoPolicy(newSize - size_), ptr_ + size_, - ptr_ + newSize, val); - } - size_ = newSize; - if (shrink) shrink_to_fit(); - } - - void shrink_to_fit() { - T *newBuffer = nullptr; - if (size_ > 0) { - newBuffer = reinterpret_cast(malloc(size_ * sizeof(T))); - if (newBuffer == nullptr) throw std::bad_alloc(); - TracyAllocS(newBuffer, size_ * sizeof(T), 3); - uninitialized_copy(autoPolicy(size_), ptr_, ptr_ + size_, newBuffer); - if (ptr_ != nullptr) { - TracyFreeS(ptr_, 3); - free(ptr_); - } - } - ptr_ = newBuffer; - capacity_ = size_; - } - -#ifdef MANIFOLD_DEBUG - void Dump() const { - std::cout << "VecDH = " << std::endl; - for (int i = 0; i < size_; ++i) { - std::cout << i << ", " << ptr_[i] << ", " << std::endl; - } - std::cout << std::endl; - } -#endif - - private: - int size_ = 0; - int capacity_ = 0; - T *ptr_ = nullptr; -}; - -template -class VecDc { - public: - VecDc(const VecDH &vec) : ptr_(vec.ptrD()), size_(vec.size()) {} - - __host__ __device__ const T &operator[](int i) const { return ptr_[i]; } - __host__ __device__ int size() const { return size_; } - - private: - T const *const ptr_; - const int size_; -}; - -template -class VecD { - public: - VecD(VecDH &vec) : ptr_(vec.ptrD()), size_(vec.size()) {} - - __host__ __device__ T &operator[](int i) const { return ptr_[i]; } - __host__ __device__ int size() const { return size_; } - - private: - T *ptr_; - int size_; -}; -/** @} */ -} // namespace manifold diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index d612dc289..8e3d595e9 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -47,7 +47,7 @@ if(NOT TARGET GTest::GTest) endif() add_executable(${PROJECT_NAME} ${SOURCE_FILES}) -target_link_libraries(${PROJECT_NAME} polygon GTest::GTest manifold sdf samples samplesGPU cross_section) +target_link_libraries(${PROJECT_NAME} polygon GTest::GTest manifold sdf samples cross_section) if(MANIFOLD_CBIND) target_link_libraries(${PROJECT_NAME} manifoldc) diff --git a/test/manifold_test.cpp b/test/manifold_test.cpp index f4514be5d..96f6dfde8 100644 --- a/test/manifold_test.cpp +++ b/test/manifold_test.cpp @@ -14,6 +14,8 @@ #include "manifold.h" +#include + #include "cross_section.h" #include "test.h" diff --git a/test/sdf_test.cpp b/test/sdf_test.cpp index 90139e705..00f486250 100644 --- a/test/sdf_test.cpp +++ b/test/sdf_test.cpp @@ -24,7 +24,7 @@ using namespace manifold; struct CubeVoid { - __host__ __device__ float operator()(glm::vec3 p) const { + float operator()(glm::vec3 p) const { const glm::vec3 min = p + glm::vec3(1); const glm::vec3 max = glm::vec3(1) - p; const float min3 = glm::min(min.x, glm::min(min.y, min.z)); @@ -34,7 +34,7 @@ struct CubeVoid { }; struct Layers { - __host__ __device__ float operator()(glm::vec3 p) const { + float operator()(glm::vec3 p) const { int a = glm::mod(glm::round(2 * p.z), 4.0f); return a == 0 ? 1 : (a == 2 ? -1 : 0); } @@ -114,4 +114,4 @@ TEST(SDF, Resize) { EXPECT_EQ(layers.Status(), Manifold::Error::NoError); EXPECT_EQ(layers.Genus(), -8); -} \ No newline at end of file +} diff --git a/test/test_main.cpp b/test/test_main.cpp index 964295416..a9afc6a7f 100644 --- a/test/test_main.cpp +++ b/test/test_main.cpp @@ -11,6 +11,7 @@ // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. // See the License for the specific language governing permissions and // limitations under the License. +#include #include "manifold.h" #include "polygon.h" @@ -126,7 +127,7 @@ Mesh Csaszar() { } struct GyroidSDF { - __host__ __device__ float operator()(glm::vec3 p) const { + float operator()(glm::vec3 p) const { const glm::vec3 min = p; const glm::vec3 max = glm::vec3(glm::two_pi()) - p; const float min3 = glm::min(min.x, glm::min(min.y, min.z));