From 41774c6a9a2749cae49117459bab0ae3f463f34d Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 23 Aug 2024 16:32:39 +0200 Subject: [PATCH 01/18] delete SplitCellLinkedList --- src/shared/common/sph_data_containers.h | 2 - src/shared/meshes/cell_linked_list.cpp | 47 +---------- src/shared/meshes/cell_linked_list.h | 20 ----- .../particle_dynamics/particle_iterators.h | 80 ------------------- 4 files changed, 1 insertion(+), 148 deletions(-) diff --git a/src/shared/common/sph_data_containers.h b/src/shared/common/sph_data_containers.h index d597b3a324..a51e60dcd4 100644 --- a/src/shared/common/sph_data_containers.h +++ b/src/shared/common/sph_data_containers.h @@ -62,8 +62,6 @@ using ListData = std::pair; using ListDataVector = StdLargeVec; using DataListsInCells = StdLargeVec; using ConcurrentCellLists = ConcurrentVec; -/** Cell list for splitting algorithms. */ -using SplitCellLists = StdVec; /** Cell list for periodic boundary condition algorithms. */ using CellLists = std::pair; diff --git a/src/shared/meshes/cell_linked_list.cpp b/src/shared/meshes/cell_linked_list.cpp index ec5d2b2091..837ef78401 100644 --- a/src/shared/meshes/cell_linked_list.cpp +++ b/src/shared/meshes/cell_linked_list.cpp @@ -14,36 +14,13 @@ BaseCellLinkedList:: : BaseMeshField("CellLinkedList"), kernel_(*sph_adaptation.getKernel()) {} //=================================================================================================// -SplitCellLists *BaseCellLinkedList::getSplitCellLists() -{ - std::cout << "\n Error: SplitCellList not defined!" << std::endl; - std::cout << __FILE__ << ':' << __LINE__ << std::endl; - exit(1); - return nullptr; -} -//=================================================================================================// -void BaseCellLinkedList::setUseSplitCellLists() -{ - std::cout << "\n Error: SplitCellList not defined!" << std::endl; - std::cout << __FILE__ << ':' << __LINE__ << std::endl; - exit(1); -}; -//=================================================================================================// -void BaseCellLinkedList::clearSplitCellLists(SplitCellLists &split_cell_lists) -{ - for (size_t i = 0; i < split_cell_lists.size(); i++) - split_cell_lists[i].clear(); -} -//=================================================================================================// CellLinkedList::CellLinkedList(BoundingBox tentative_bounds, Real grid_spacing, SPHAdaptation &sph_adaptation) : BaseCellLinkedList(sph_adaptation), Mesh(tentative_bounds, grid_spacing, 2), - use_split_cell_lists_(false), cell_index_lists_(nullptr), cell_data_lists_(nullptr) + cell_index_lists_(nullptr), cell_data_lists_(nullptr) { allocateMeshDataMatrix(); single_cell_linked_list_level_.push_back(this); - size_t number_of_split_cell_lists = pow(3, Dimensions); - split_cell_lists_.resize(number_of_split_cell_lists); } //=================================================================================================// void CellLinkedList ::allocateMeshDataMatrix() @@ -86,23 +63,6 @@ void CellLinkedList::UpdateCellListData(BaseParticles &base_particles) }); } //=================================================================================================// -void CellLinkedList::updateSplitCellLists(SplitCellLists &split_cell_lists) -{ - clearSplitCellLists(split_cell_lists); - mesh_parallel_for( - MeshRange(Arrayi::Zero(), all_cells_), - [&](const Arrayi &cell_index) - { - ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); - size_t real_particles_in_cell = cell_list.size(); - if (real_particles_in_cell != 0) - { - split_cell_lists[transferMeshIndexTo1D(3 * Arrayi::Ones(), mod(cell_index, 3))] - .push_back(&cell_list); - } - }); -} -//=================================================================================================// void CellLinkedList::UpdateCellLists(BaseParticles &base_particles) { clearCellLists(); @@ -120,11 +80,6 @@ void CellLinkedList::UpdateCellLists(BaseParticles &base_particles) ap); UpdateCellListData(base_particles); - - if (use_split_cell_lists_) - { - updateSplitCellLists(split_cell_lists_); - } } //=================================================================================================// void CellLinkedList ::insertParticleIndex(size_t particle_index, const Vecd &particle_position) diff --git a/src/shared/meshes/cell_linked_list.h b/src/shared/meshes/cell_linked_list.h index f1f243cfac..a6fe1bb5f6 100644 --- a/src/shared/meshes/cell_linked_list.h +++ b/src/shared/meshes/cell_linked_list.h @@ -52,11 +52,6 @@ class BaseCellLinkedList : public BaseMeshField protected: Kernel &kernel_; - /** clear split cell lists in this mesh*/ - virtual void clearSplitCellLists(SplitCellLists &split_cell_lists); - /** update split particle list in this mesh */ - virtual void updateSplitCellLists(SplitCellLists &split_cell_lists) = 0; - public: BaseCellLinkedList(SPHAdaptation &sph_adaptation); virtual ~BaseCellLinkedList(){}; @@ -64,8 +59,6 @@ class BaseCellLinkedList : public BaseMeshField /** access concrete cell linked list levels*/ virtual StdVec CellLinkedListLevels() = 0; virtual void UpdateCellLists(BaseParticles &base_particles) = 0; - virtual SplitCellLists *getSplitCellLists(); - virtual void setUseSplitCellLists(); /** Insert a cell-linked_list entry to the concurrent index list. */ virtual void insertParticleIndex(size_t particle_index, const Vecd &particle_position) = 0; /** Insert a cell-linked_list entry of the index and particle position pair. */ @@ -88,14 +81,6 @@ class BaseCellLinkedList : public BaseMeshField class CellLinkedList : public BaseCellLinkedList, public Mesh { StdVec single_cell_linked_list_level_; - /** - * @brief particle by cells lists is for parallel splitting algorithm. - * All particles in each cell are collected together. - * If two particles each belongs two different cell entries, - * they have no interaction because they are too far. - */ - SplitCellLists split_cell_lists_; - bool use_split_cell_lists_; protected: /** using concurrent vectors due to writing conflicts when building the list */ @@ -105,7 +90,6 @@ class CellLinkedList : public BaseCellLinkedList, public Mesh void allocateMeshDataMatrix(); /**< allocate memories for addresses of data packages. */ void deleteMeshDataMatrix(); /**< delete memories for addresses of data packages. */ - virtual void updateSplitCellLists(SplitCellLists &split_cell_lists) override; template DataListsType &getCellDataList(DataListsType *data_lists, const Arrayi &cell_index) { @@ -117,8 +101,6 @@ class CellLinkedList : public BaseCellLinkedList, public Mesh ~CellLinkedList() { deleteMeshDataMatrix(); }; void clearCellLists(); - virtual SplitCellLists *getSplitCellLists() override { return &split_cell_lists_; }; - virtual void setUseSplitCellLists() override { use_split_cell_lists_ = true; }; void UpdateCellListData(BaseParticles &base_particles); virtual void UpdateCellLists(BaseParticles &base_particles) override; void insertParticleIndex(size_t particle_index, const Vecd &particle_position) override; @@ -153,8 +135,6 @@ class MultilevelCellLinkedList : public MultilevelMesh &h_ratio_; /**< Smoothing length for each level. */ - /** Update split cell list. */ - virtual void updateSplitCellLists(SplitCellLists &split_cell_lists) override{}; /** determine mesh level from particle cutoff radius */ inline size_t getMeshLevel(Real particle_cutoff_radius); diff --git a/src/shared/particle_dynamics/particle_iterators.h b/src/shared/particle_dynamics/particle_iterators.h index f101d56a22..f3c8ab08ce 100644 --- a/src/shared/particle_dynamics/particle_iterators.h +++ b/src/shared/particle_dynamics/particle_iterators.h @@ -161,86 +161,6 @@ inline void particle_for(const ParallelPolicy &par, const DataListsInCells &body }, ap); }; -/** - * Splitting algorithm (for sequential and parallel computing). - */ -template -inline void particle_for(const SequencedPolicy &seq, const SplitCellLists &split_cell_lists, - const LocalDynamicsFunction &local_dynamics_function) -{ - // forward sweeping - for (size_t k = 0; k != split_cell_lists.size(); ++k) - { - const ConcurrentCellLists &cell_lists = split_cell_lists[k]; - for (size_t l = 0; l != cell_lists.size(); ++l) - { - const ConcurrentIndexVector &particle_indexes = *cell_lists[l]; - for (size_t i = 0; i != particle_indexes.size(); ++i) - { - local_dynamics_function(particle_indexes[i]); - } - } - } - - // backward sweeping - for (size_t k = split_cell_lists.size(); k != 0; --k) - { - const ConcurrentCellLists &cell_lists = split_cell_lists[k - 1]; - for (size_t l = 0; l != cell_lists.size(); ++l) - { - const ConcurrentIndexVector &particle_indexes = *cell_lists[l]; - for (size_t i = particle_indexes.size(); i != 0; --i) - { - local_dynamics_function(particle_indexes[i - 1]); - } - } - } -} - -template -inline void particle_for(const ParallelPolicy &par, const SplitCellLists &split_cell_lists, - const LocalDynamicsFunction &local_dynamics_function) -{ - // forward sweeping - for (size_t k = 0; k != split_cell_lists.size(); ++k) - { - const ConcurrentCellLists &cell_lists = split_cell_lists[k]; - parallel_for( - IndexRange(0, cell_lists.size()), - [&](const IndexRange &r) - { - for (size_t l = r.begin(); l < r.end(); ++l) - { - const ConcurrentIndexVector &particle_indexes = *cell_lists[l]; - for (size_t i = 0; i < particle_indexes.size(); ++i) - { - local_dynamics_function(particle_indexes[i]); - } - } - }, - ap); - } - - // backward sweeping - for (size_t k = split_cell_lists.size(); k != 0; --k) - { - const ConcurrentCellLists &cell_lists = split_cell_lists[k - 1]; - parallel_for( - IndexRange(0, cell_lists.size()), - [&](const IndexRange &r) - { - for (size_t l = r.begin(); l < r.end(); ++l) - { - const ConcurrentIndexVector &particle_indexes = *cell_lists[l]; - for (size_t i = particle_indexes.size(); i != 0; --i) - { - local_dynamics_function(particle_indexes[i - 1]); - } - } - }, - ap); - } -} template From 4a7723ec529e16b6944e6205a8c3c2f6f5749b0a Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 23 Aug 2024 16:44:36 +0200 Subject: [PATCH 02/18] add split method in cell linked list --- src/for_2D_build/meshes/mesh_iterators.hpp | 42 ++++++++++++++++ src/for_3D_build/meshes/mesh_iterators.hpp | 48 +++++++++++++++++++ src/shared/meshes/cell_linked_list.h | 4 ++ src/shared/meshes/cell_linked_list.hpp | 16 +++++++ src/shared/meshes/mesh_iterators.h | 3 ++ .../particle_dynamics_algorithms.h | 11 ++--- 6 files changed, 118 insertions(+), 6 deletions(-) diff --git a/src/for_2D_build/meshes/mesh_iterators.hpp b/src/for_2D_build/meshes/mesh_iterators.hpp index 57bd4d993c..6a08e8dd7f 100644 --- a/src/for_2D_build/meshes/mesh_iterators.hpp +++ b/src/for_2D_build/meshes/mesh_iterators.hpp @@ -84,4 +84,46 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f ap); } //=================================================================================================// +template +void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ + // forward sweeping + for (int m = 0; m < stride[0]; m++) + for (int n = 0; n < stride[1]; n++) + { + parallel_for( + IndexRange2d((mesh_range.first)[0] + m, (mesh_range.second)[0], + (mesh_range.first)[1] + n, (mesh_range.second)[1]), + [&](const IndexRange2d &r) + { + for (size_t i = r.rows().begin(); i != r.rows().end(); i++) + for (size_t j = r.cols().begin(); j != r.cols().end(); j++) + { + if ((i - m) % stride[0] == 0 && (j - n) % stride[1] == 0) + local_function(Array2i(i, j)); + } + }, + ap); + } + + // backward sweeping + for (int m = stride[0] - 1; m >= 0; m--) + for (int n = stride[1] - 1; n >= 0; n--) + { + parallel_for( + IndexRange2d((mesh_range.first)[0] + m, (mesh_range.second)[0], + (mesh_range.first)[1] + n, (mesh_range.second)[1]), + [&](const IndexRange2d &r) + { + for (size_t i = r.rows().end(); i != r.rows().begin(); i--) + for (size_t j = r.cols().end(); j != r.cols().begin(); j--) + { + if ((i - 1 - m) % stride[0] == 0 && (j - 1 - n) % stride[1] == 0) + local_function(Array2i(i - 1, j - 1)); + } + }, + ap); + } +} +//=================================================================================================// } // namespace SPH diff --git a/src/for_3D_build/meshes/mesh_iterators.hpp b/src/for_3D_build/meshes/mesh_iterators.hpp index 4221a3728c..c7f9f8e8a4 100644 --- a/src/for_3D_build/meshes/mesh_iterators.hpp +++ b/src/for_3D_build/meshes/mesh_iterators.hpp @@ -94,4 +94,52 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f ap); } //=================================================================================================// +template +void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ + // forward sweeping + for (int m = 0; m < stride[0]; m++) + for (int n = 0; n < stride[1]; n++) + for (int p = 0; p < stride[2]; p++) + { + parallel_for( + IndexRange3d((mesh_range.first)[0] + m, (mesh_range.second)[0], + (mesh_range.first)[1] + n, (mesh_range.second)[1], + (mesh_range.first)[2] + p, (mesh_range.second)[2]), + [&](const IndexRange3d &r) + { + for (size_t i = r.pages().begin(); i != r.pages().end(); ++i) + for (size_t j = r.rows().begin(); j != r.rows().end(); ++j) + for (size_t k = r.cols().begin(); k != r.cols().end(); ++k) + { + if ((i - m) % stride[0] == 0 && (j - n) % stride[1] == 0 && (k - p) % stride[2] == 0) + local_function(Array3i(i, j, k)); + } + }, + ap); + } + + // backward sweeping + for (int m = stride[0] - 1; m >= 0; m--) + for (int n = stride[1] - 1; n >= 0; n--) + for (int p = stride[2] - 1; p >= 0; p--) + { + parallel_for( + IndexRange3d((mesh_range.first)[0] + m, (mesh_range.second)[0], + (mesh_range.first)[1] + n, (mesh_range.second)[1], + (mesh_range.first)[2] + p, (mesh_range.second)[2]), + [&](const IndexRange3d &r) + { + for (size_t i = r.pages().end(); i != r.pages().begin(); i--) + for (size_t j = r.rows().end(); j != r.rows().begin(); j--) + for (size_t k = r.cols().end(); k != r.cols().begin(); k--) + { + if ((i - 1 - m) % stride[0] == 0 && (j - 1 - n) % stride[1] == 0 && (k - 1 - p) % stride[2] == 0) + local_function(Array3i(i - 1, j - 1, k - 1)); + } + }, + ap); + } +} +//=================================================================================================// } // namespace SPH diff --git a/src/shared/meshes/cell_linked_list.h b/src/shared/meshes/cell_linked_list.h index a6fe1bb5f6..838d2cfe85 100644 --- a/src/shared/meshes/cell_linked_list.h +++ b/src/shared/meshes/cell_linked_list.h @@ -116,6 +116,10 @@ class CellLinkedList : public BaseCellLinkedList, public Mesh template void searchNeighborsByParticles(DynamicsRange &dynamics_range, ParticleConfiguration &particle_configuration, GetSearchDepth &get_search_depth, GetNeighborRelation &get_neighbor_relation); + + /** split algorithm */; + template + void particle_for_split(const LocalDynamicsFunction &local_dynamics_function); }; template <> diff --git a/src/shared/meshes/cell_linked_list.hpp b/src/shared/meshes/cell_linked_list.hpp index 5bf984f0d4..5f3eaf195e 100644 --- a/src/shared/meshes/cell_linked_list.hpp +++ b/src/shared/meshes/cell_linked_list.hpp @@ -44,4 +44,20 @@ void CellLinkedList::searchNeighborsByParticles( }); } //=================================================================================================// +template +void CellLinkedList::particle_for_split(const LocalDynamicsFunction &local_dynamics_function) +{ + mesh_split_parallel_for( + MeshRange(Arrayi::Zero(), all_cells_), + 3 * Arrayi::Ones(), + [&](const Arrayi &cell_index) + { + const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + for (size_t index_i : cell_list) + { + local_dynamics_function(index_i); + } + }); +} +//=================================================================================================// } // namespace SPH diff --git a/src/shared/meshes/mesh_iterators.h b/src/shared/meshes/mesh_iterators.h index d9adb1a3bd..53de29613b 100644 --- a/src/shared/meshes/mesh_iterators.h +++ b/src/shared/meshes/mesh_iterators.h @@ -107,5 +107,8 @@ void mesh_for(const MeshRange &mesh_range, const LocalFunction &local_function, /** Iterator on the mesh by looping index. parallel computing. */ template void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_function, Args &&...args); +/** Iterator on the mesh by looping index with a stride. parallel computing. */ +template +void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); } // namespace SPH #endif // MESH_ITERATORS_H \ No newline at end of file diff --git a/src/shared/particle_dynamics/particle_dynamics_algorithms.h b/src/shared/particle_dynamics/particle_dynamics_algorithms.h index ddcaddfc54..8964f5ff50 100644 --- a/src/shared/particle_dynamics/particle_dynamics_algorithms.h +++ b/src/shared/particle_dynamics/particle_dynamics_algorithms.h @@ -52,6 +52,7 @@ #include "base_local_dynamics.h" #include "base_particle_dynamics.h" +#include "cell_linked_list.hpp" #include "particle_iterators.h" #include @@ -198,16 +199,15 @@ class InteractionSplit : public BaseInteractionDynamics InteractionSplit(Args &&... args) : BaseInteractionDynamics(std::forward(args)...), real_body_(DynamicCast(this, this->getSPHBody())), - split_cell_lists_(*real_body_.getCellLinkedList().getSplitCellLists()) + cell_linked_list_(DynamicCast(this, real_body_.getCellLinkedList())) { - real_body_.getCellLinkedList().setUseSplitCellLists(); static_assert(!has_initialize::value && !has_update::value, "LocalDynamicsType does not fulfill InteractionSplit requirements"); @@ -217,9 +217,8 @@ class InteractionSplit : public BaseInteractionDynamicsinteraction(i, dt * 0.5); }); + cell_linked_list_.particle_for_split([&](size_t i) + { this->interaction(i, dt * 0.5); }); } }; From 91a4092197d7ba92a2e189a6144e42af12ec4c9d Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 23 Aug 2024 16:50:43 +0200 Subject: [PATCH 03/18] fix bug in old test cases --- .../test_2d_three_ring_impact/test_2d_three_ring_impact.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/2d_examples/test_2d_three_ring_impact/test_2d_three_ring_impact.cpp b/tests/2d_examples/test_2d_three_ring_impact/test_2d_three_ring_impact.cpp index 9ade8b2e0a..0f747b5bc0 100644 --- a/tests/2d_examples/test_2d_three_ring_impact/test_2d_three_ring_impact.cpp +++ b/tests/2d_examples/test_2d_three_ring_impact/test_2d_three_ring_impact.cpp @@ -299,7 +299,7 @@ void three_ring_impact(int resolution_factor_l, int resolution_factor_m, int res { const auto &pos_ = particles->ParticlePositions(); for (const auto &pos : pos_) - if (std::isnan(pos[0]) || std::isnan(pos[1]) || std::isnan(pos[2])) + if (std::isnan(pos[0]) || std::isnan(pos[1])) throw std::runtime_error("position has become nan"); }; From 0caef4c953ad8f0fba81008384b1ccee628ed624 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Tue, 10 Sep 2024 13:18:51 +0200 Subject: [PATCH 04/18] add sequential loop --- src/for_2D_build/meshes/mesh_iterators.hpp | 36 +++++++++++++++-- src/for_3D_build/meshes/mesh_iterators.hpp | 46 +++++++++++++++++++--- src/shared/meshes/mesh_iterators.h | 3 ++ 3 files changed, 75 insertions(+), 10 deletions(-) diff --git a/src/for_2D_build/meshes/mesh_iterators.hpp b/src/for_2D_build/meshes/mesh_iterators.hpp index 6a08e8dd7f..a972db5293 100644 --- a/src/for_2D_build/meshes/mesh_iterators.hpp +++ b/src/for_2D_build/meshes/mesh_iterators.hpp @@ -85,6 +85,34 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f } //=================================================================================================// template +void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ + // forward sweeping + for (int m = 0; m < stride[0]; m++) + for (int n = 0; n < stride[1]; n++) + for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) + for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) + { + local_function(Array2i(i, j)); + } + + // backward sweeping + for (int m = stride[0] - 1; m >= 0; m--) + for (int n = stride[1] - 1; n >= 0; n--) + { + const auto index_total_i = (mesh_range.second)[0] - 1 - (mesh_range.first)[0] - m; + const auto index_total_j = (mesh_range.second)[1] - 1 - (mesh_range.first)[1] - n; + const auto i_max = (mesh_range.second)[0] - 1 - index_total_i % stride[0]; + const auto j_max = (mesh_range.second)[1] - 1 - index_total_j % stride[1]; + for (auto i = i_max; i >= (mesh_range.first)[0] + m; i -= stride[0]) + for (auto j = j_max; j >= (mesh_range.first)[1] + n; j -= stride[1]) + { + local_function(Array2i(i, j)); + } + } +} +//=================================================================================================// +template void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { // forward sweeping @@ -96,8 +124,8 @@ void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, (mesh_range.first)[1] + n, (mesh_range.second)[1]), [&](const IndexRange2d &r) { - for (size_t i = r.rows().begin(); i != r.rows().end(); i++) - for (size_t j = r.cols().begin(); j != r.cols().end(); j++) + for (auto i = r.rows().begin(); i != r.rows().end(); i++) + for (auto j = r.cols().begin(); j != r.cols().end(); j++) { if ((i - m) % stride[0] == 0 && (j - n) % stride[1] == 0) local_function(Array2i(i, j)); @@ -115,8 +143,8 @@ void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, (mesh_range.first)[1] + n, (mesh_range.second)[1]), [&](const IndexRange2d &r) { - for (size_t i = r.rows().end(); i != r.rows().begin(); i--) - for (size_t j = r.cols().end(); j != r.cols().begin(); j--) + for (auto i = r.rows().end(); i != r.rows().begin(); i--) + for (auto j = r.cols().end(); j != r.cols().begin(); j--) { if ((i - 1 - m) % stride[0] == 0 && (j - 1 - n) % stride[1] == 0) local_function(Array2i(i - 1, j - 1)); diff --git a/src/for_3D_build/meshes/mesh_iterators.hpp b/src/for_3D_build/meshes/mesh_iterators.hpp index c7f9f8e8a4..af82229449 100644 --- a/src/for_3D_build/meshes/mesh_iterators.hpp +++ b/src/for_3D_build/meshes/mesh_iterators.hpp @@ -95,6 +95,40 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f } //=================================================================================================// template +void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ + // forward sweeping + for (int m = 0; m < stride[0]; m++) + for (int n = 0; n < stride[1]; n++) + for (int p = 0; p < stride[2]; p++) + for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) + for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) + for (auto k = (mesh_range.first)[2] + p; k < (mesh_range.second)[2]; k += stride[2]) + { + local_function(Array3i(i, j, k)); + } + + // backward sweeping + for (int m = stride[0] - 1; m >= 0; m--) + for (int n = stride[1] - 1; n >= 0; n--) + for (int p = stride[2] - 1; p >= 0; p--) + { + const auto index_total_i = (mesh_range.second)[0] - 1 - (mesh_range.first)[0] - m; + const auto index_total_j = (mesh_range.second)[1] - 1 - (mesh_range.first)[1] - n; + const auto index_total_k = (mesh_range.second)[2] - 1 - (mesh_range.first)[2] - p; + const auto i_max = (mesh_range.second)[0] - 1 - index_total_i % stride[0]; + const auto j_max = (mesh_range.second)[1] - 1 - index_total_j % stride[1]; + const auto k_max = (mesh_range.second)[2] - 1 - index_total_k % stride[2]; + for (auto i = i_max; i >= (mesh_range.first)[0] + m; i -= stride[0]) + for (auto j = j_max; j >= (mesh_range.first)[1] + n; j -= stride[1]) + for (auto k = k_max; k >= (mesh_range.first)[2] + p; k -= stride[2]) + { + local_function(Array3i(i, j, k)); + } + } +} +//=================================================================================================// +template void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { // forward sweeping @@ -108,9 +142,9 @@ void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, (mesh_range.first)[2] + p, (mesh_range.second)[2]), [&](const IndexRange3d &r) { - for (size_t i = r.pages().begin(); i != r.pages().end(); ++i) - for (size_t j = r.rows().begin(); j != r.rows().end(); ++j) - for (size_t k = r.cols().begin(); k != r.cols().end(); ++k) + for (auto i = r.pages().begin(); i != r.pages().end(); ++i) + for (auto j = r.rows().begin(); j != r.rows().end(); ++j) + for (auto k = r.cols().begin(); k != r.cols().end(); ++k) { if ((i - m) % stride[0] == 0 && (j - n) % stride[1] == 0 && (k - p) % stride[2] == 0) local_function(Array3i(i, j, k)); @@ -130,9 +164,9 @@ void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, (mesh_range.first)[2] + p, (mesh_range.second)[2]), [&](const IndexRange3d &r) { - for (size_t i = r.pages().end(); i != r.pages().begin(); i--) - for (size_t j = r.rows().end(); j != r.rows().begin(); j--) - for (size_t k = r.cols().end(); k != r.cols().begin(); k--) + for (auto i = r.pages().end(); i != r.pages().begin(); i--) + for (auto j = r.rows().end(); j != r.rows().begin(); j--) + for (auto k = r.cols().end(); k != r.cols().begin(); k--) { if ((i - 1 - m) % stride[0] == 0 && (j - 1 - n) % stride[1] == 0 && (k - 1 - p) % stride[2] == 0) local_function(Array3i(i - 1, j - 1, k - 1)); diff --git a/src/shared/meshes/mesh_iterators.h b/src/shared/meshes/mesh_iterators.h index 53de29613b..1da19ee956 100644 --- a/src/shared/meshes/mesh_iterators.h +++ b/src/shared/meshes/mesh_iterators.h @@ -107,6 +107,9 @@ void mesh_for(const MeshRange &mesh_range, const LocalFunction &local_function, /** Iterator on the mesh by looping index. parallel computing. */ template void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_function, Args &&...args); +/** Iterator on the mesh by looping index with a stride. sequential computing. */ +template +void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); /** Iterator on the mesh by looping index with a stride. parallel computing. */ template void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); From 050ab2a4d1c2dcfec7553ea683d91f2f6fa6b478 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Tue, 10 Sep 2024 13:19:03 +0200 Subject: [PATCH 05/18] add unit tests --- .../for_2D_build/iterators/CMakeLists.txt | 7 +++ .../test_2d_mesh_iterator/CMakeLists.txt | 15 ++++++ .../test_mesh_iterator.cpp | 53 +++++++++++++++++++ .../for_3D_build/iterators/CMakeLists.txt | 7 +++ .../test_3d_mesh_iterator/CMakeLists.txt | 16 ++++++ .../test_mesh_iterator.cpp | 53 +++++++++++++++++++ 6 files changed, 151 insertions(+) create mode 100644 tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp create mode 100644 tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp diff --git a/tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt b/tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt new file mode 100644 index 0000000000..8eb174fa72 --- /dev/null +++ b/tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt @@ -0,0 +1,7 @@ +SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR}) + +foreach(subdir ${SUBDIRS}) + if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt) + add_subdirectory(${subdir}) + endif() +endforeach() \ No newline at end of file diff --git a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt new file mode 100644 index 0000000000..c5b091dac7 --- /dev/null +++ b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt @@ -0,0 +1,15 @@ +STRING( REGEX REPLACE ".*/(.*)" "\\1" CURRENT_FOLDER ${CMAKE_CURRENT_SOURCE_DIR} ) +PROJECT("${CURRENT_FOLDER}") + +SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib) +SET(EXECUTABLE_OUTPUT_PATH "${PROJECT_BINARY_DIR}/bin/") +SET(BUILD_INPUT_PATH "${EXECUTABLE_OUTPUT_PATH}/input") +SET(BUILD_RELOAD_PATH "${EXECUTABLE_OUTPUT_PATH}/reload") + +aux_source_directory(. DIR_SRCS) +ADD_EXECUTABLE(${PROJECT_NAME} ${EXECUTABLE_OUTPUT_PATH} ${DIR_SRCS}) +target_link_libraries(${PROJECT_NAME} sphinxsys_2d GTest::gtest GTest::gtest_main) +set_target_properties(${PROJECT_NAME} PROPERTIES VS_DEBUGGER_WORKING_DIRECTORY "${EXECUTABLE_OUTPUT_PATH}") + +add_test(NAME ${PROJECT_NAME} COMMAND ${PROJECT_NAME} + WORKING_DIRECTORY ${EXECUTABLE_OUTPUT_PATH}) diff --git a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp new file mode 100644 index 0000000000..152e23f371 --- /dev/null +++ b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp @@ -0,0 +1,53 @@ +#include "base_data_type.h" +#include "mesh_iterators.hpp" +#include + +using namespace SPH; + +size_t get_cell_index(const Array2i &index) +{ + return 3 * (index[0] % 3) + (index[1] % 3); +} + +size_t get_1d_index(const Array2i &all_cells, const Array2i &index) +{ + return all_cells[1] * index[0] + index[1]; +} + +void function(const Array2i &index, std::vector &vec, const Array2i &all_cells) +{ + size_t index_cell = get_cell_index(index); + size_t index_1d = get_1d_index(all_cells, index); + vec[index_1d] += index_cell; +} + +TEST(test_mesh_iterator, split_for) +{ + Array2i all_cells(9, 9); + MeshRange mesh_range(Array2i::Zero(), all_cells); + Array2i stride = 3 * Array2i::Ones(); + + const size_t total_cells = all_cells[0] * all_cells[1]; + + std::vector vec_1(total_cells); + std::vector vec_2(total_cells); + + mesh_split_for( + mesh_range, stride, + [&](const Array2i &index) + { + function(index, vec_1, all_cells); + }); + + mesh_split_parallel_for( + mesh_range, stride, + [&](const Array2i &index) + { + function(index, vec_2, all_cells); + }); + + for (size_t i = 0; i < total_cells; ++i) + { + ASSERT_EQ(vec_1[i], vec_2[i]); + } +} \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt b/tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt new file mode 100644 index 0000000000..8eb174fa72 --- /dev/null +++ b/tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt @@ -0,0 +1,7 @@ +SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR}) + +foreach(subdir ${SUBDIRS}) + if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt) + add_subdirectory(${subdir}) + endif() +endforeach() \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt new file mode 100644 index 0000000000..d9a298a2e7 --- /dev/null +++ b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt @@ -0,0 +1,16 @@ +STRING( REGEX REPLACE ".*/(.*)" "\\1" CURRENT_FOLDER ${CMAKE_CURRENT_SOURCE_DIR} ) +PROJECT("${CURRENT_FOLDER}") + +SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib) +SET(EXECUTABLE_OUTPUT_PATH "${PROJECT_BINARY_DIR}/bin/") +SET(BUILD_INPUT_PATH "${EXECUTABLE_OUTPUT_PATH}/input") +SET(BUILD_RELOAD_PATH "${EXECUTABLE_OUTPUT_PATH}/reload") + +aux_source_directory(. DIR_SRCS) +ADD_EXECUTABLE(${PROJECT_NAME} ${EXECUTABLE_OUTPUT_PATH} ${DIR_SRCS}) +target_link_libraries(${PROJECT_NAME} sphinxsys_3d GTest::gtest GTest::gtest_main) +set_target_properties(${PROJECT_NAME} PROPERTIES VS_DEBUGGER_WORKING_DIRECTORY "${EXECUTABLE_OUTPUT_PATH}") + +add_test(NAME ${PROJECT_NAME}_particle_relaxation + COMMAND ${PROJECT_NAME} --relax=true + WORKING_DIRECTORY ${EXECUTABLE_OUTPUT_PATH}) diff --git a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp new file mode 100644 index 0000000000..bd5ef27b2b --- /dev/null +++ b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp @@ -0,0 +1,53 @@ +#include "base_data_type.h" +#include "mesh_iterators.hpp" +#include + +using namespace SPH; + +size_t get_cell_index(const Array3i &index) +{ + return 9 * (index[0] % 3) + 3 * (index[1] % 3) + (index[2] % 3); +} + +size_t get_1d_index(const Array3i &all_cells, const Array3i &index) +{ + return all_cells[1] * all_cells[2] * index[0] + all_cells[2] * index[1] + index[2]; +} + +void function(const Array3i &index, std::vector &vec, const Array3i &all_cells) +{ + size_t index_cell = get_cell_index(index); + size_t index_1d = get_1d_index(all_cells, index); + vec[index_1d] += index_cell; +} + +TEST(test_mesh_iterator, split_for) +{ + Array3i all_cells(9, 9, 9); + MeshRange mesh_range(Array3i::Zero(), all_cells); + Array3i stride = 3 * Array3i::Ones(); + + const size_t total_cells = all_cells[0] * all_cells[1] * all_cells[2]; + + std::vector vec_1(total_cells); + std::vector vec_2(total_cells); + + mesh_split_for( + mesh_range, stride, + [&](const Array3i &index) + { + function(index, vec_1, all_cells); + }); + + mesh_split_parallel_for( + mesh_range, stride, + [&](const Array3i &index) + { + function(index, vec_2, all_cells); + }); + + for (size_t i = 0; i < total_cells; ++i) + { + ASSERT_EQ(vec_1[i], vec_2[i]); + } +} \ No newline at end of file From 46baf491fded1579b120dbb6b48323dbf36cff44 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Thu, 12 Sep 2024 11:53:07 +0200 Subject: [PATCH 06/18] separate forward and backward sweeping --- src/for_2D_build/meshes/mesh_iterators.hpp | 72 ++++++---------- src/for_3D_build/meshes/mesh_iterators.hpp | 84 +++++++------------ src/shared/meshes/mesh_iterators.h | 14 +++- .../test_mesh_iterator.cpp | 27 ++++-- .../test_mesh_iterator.cpp | 27 ++++-- 5 files changed, 110 insertions(+), 114 deletions(-) diff --git a/src/for_2D_build/meshes/mesh_iterators.hpp b/src/for_2D_build/meshes/mesh_iterators.hpp index a972db5293..df4dec36cb 100644 --- a/src/for_2D_build/meshes/mesh_iterators.hpp +++ b/src/for_2D_build/meshes/mesh_iterators.hpp @@ -85,9 +85,8 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f } //=================================================================================================// template -void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +void mesh_stride_forward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - // forward sweeping for (int m = 0; m < stride[0]; m++) for (int n = 0; n < stride[1]; n++) for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) @@ -95,62 +94,43 @@ void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const Loc { local_function(Array2i(i, j)); } - - // backward sweeping +} +//=================================================================================================// +template +void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ + for (int m = 0; m < stride[0]; m++) + for (int n = 0; n < stride[1]; n++) + { + parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) + { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) + { local_function(Array2i(i, j)); }, ap); }, + ap); + } +} +//=================================================================================================// +template +void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ for (int m = stride[0] - 1; m >= 0; m--) for (int n = stride[1] - 1; n >= 0; n--) - { - const auto index_total_i = (mesh_range.second)[0] - 1 - (mesh_range.first)[0] - m; - const auto index_total_j = (mesh_range.second)[1] - 1 - (mesh_range.first)[1] - n; - const auto i_max = (mesh_range.second)[0] - 1 - index_total_i % stride[0]; - const auto j_max = (mesh_range.second)[1] - 1 - index_total_j % stride[1]; - for (auto i = i_max; i >= (mesh_range.first)[0] + m; i -= stride[0]) - for (auto j = j_max; j >= (mesh_range.first)[1] + n; j -= stride[1]) + for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) + for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) { local_function(Array2i(i, j)); } - } } //=================================================================================================// template -void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - // forward sweeping - for (int m = 0; m < stride[0]; m++) - for (int n = 0; n < stride[1]; n++) - { - parallel_for( - IndexRange2d((mesh_range.first)[0] + m, (mesh_range.second)[0], - (mesh_range.first)[1] + n, (mesh_range.second)[1]), - [&](const IndexRange2d &r) - { - for (auto i = r.rows().begin(); i != r.rows().end(); i++) - for (auto j = r.cols().begin(); j != r.cols().end(); j++) - { - if ((i - m) % stride[0] == 0 && (j - n) % stride[1] == 0) - local_function(Array2i(i, j)); - } - }, - ap); - } - - // backward sweeping for (int m = stride[0] - 1; m >= 0; m--) for (int n = stride[1] - 1; n >= 0; n--) { - parallel_for( - IndexRange2d((mesh_range.first)[0] + m, (mesh_range.second)[0], - (mesh_range.first)[1] + n, (mesh_range.second)[1]), - [&](const IndexRange2d &r) - { - for (auto i = r.rows().end(); i != r.rows().begin(); i--) - for (auto j = r.cols().end(); j != r.cols().begin(); j--) - { - if ((i - 1 - m) % stride[0] == 0 && (j - 1 - n) % stride[1] == 0) - local_function(Array2i(i - 1, j - 1)); - } - }, - ap); + parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) + { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) + { local_function(Array2i(i, j)); }, ap); }, + ap); } } //=================================================================================================// diff --git a/src/for_3D_build/meshes/mesh_iterators.hpp b/src/for_3D_build/meshes/mesh_iterators.hpp index af82229449..b9843c6bf8 100644 --- a/src/for_3D_build/meshes/mesh_iterators.hpp +++ b/src/for_3D_build/meshes/mesh_iterators.hpp @@ -95,9 +95,8 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f } //=================================================================================================// template -void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +void mesh_stride_forward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - // forward sweeping for (int m = 0; m < stride[0]; m++) for (int n = 0; n < stride[1]; n++) for (int p = 0; p < stride[2]; p++) @@ -107,72 +106,49 @@ void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const Loc { local_function(Array3i(i, j, k)); } - - // backward sweeping +} +//=================================================================================================// +template +void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ + for (int m = 0; m < stride[0]; m++) + for (int n = 0; n < stride[1]; n++) + for (int p = 0; p < stride[2]; p++) + { + parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) + { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) + { parallel_for((mesh_range.first)[2] + p, (mesh_range.second)[2], stride[2], [&](auto k) + { local_function(Array3i(i, j, k)); }, + ap); }, ap); }, ap); + } +} +//=================================================================================================// +template +void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +{ for (int m = stride[0] - 1; m >= 0; m--) for (int n = stride[1] - 1; n >= 0; n--) for (int p = stride[2] - 1; p >= 0; p--) - { - const auto index_total_i = (mesh_range.second)[0] - 1 - (mesh_range.first)[0] - m; - const auto index_total_j = (mesh_range.second)[1] - 1 - (mesh_range.first)[1] - n; - const auto index_total_k = (mesh_range.second)[2] - 1 - (mesh_range.first)[2] - p; - const auto i_max = (mesh_range.second)[0] - 1 - index_total_i % stride[0]; - const auto j_max = (mesh_range.second)[1] - 1 - index_total_j % stride[1]; - const auto k_max = (mesh_range.second)[2] - 1 - index_total_k % stride[2]; - for (auto i = i_max; i >= (mesh_range.first)[0] + m; i -= stride[0]) - for (auto j = j_max; j >= (mesh_range.first)[1] + n; j -= stride[1]) - for (auto k = k_max; k >= (mesh_range.first)[2] + p; k -= stride[2]) + for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) + for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) + for (auto k = (mesh_range.first)[2] + p; k < (mesh_range.second)[2]; k += stride[2]) { local_function(Array3i(i, j, k)); } - } } //=================================================================================================// template -void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) +void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - // forward sweeping - for (int m = 0; m < stride[0]; m++) - for (int n = 0; n < stride[1]; n++) - for (int p = 0; p < stride[2]; p++) - { - parallel_for( - IndexRange3d((mesh_range.first)[0] + m, (mesh_range.second)[0], - (mesh_range.first)[1] + n, (mesh_range.second)[1], - (mesh_range.first)[2] + p, (mesh_range.second)[2]), - [&](const IndexRange3d &r) - { - for (auto i = r.pages().begin(); i != r.pages().end(); ++i) - for (auto j = r.rows().begin(); j != r.rows().end(); ++j) - for (auto k = r.cols().begin(); k != r.cols().end(); ++k) - { - if ((i - m) % stride[0] == 0 && (j - n) % stride[1] == 0 && (k - p) % stride[2] == 0) - local_function(Array3i(i, j, k)); - } - }, - ap); - } - - // backward sweeping for (int m = stride[0] - 1; m >= 0; m--) for (int n = stride[1] - 1; n >= 0; n--) for (int p = stride[2] - 1; p >= 0; p--) { - parallel_for( - IndexRange3d((mesh_range.first)[0] + m, (mesh_range.second)[0], - (mesh_range.first)[1] + n, (mesh_range.second)[1], - (mesh_range.first)[2] + p, (mesh_range.second)[2]), - [&](const IndexRange3d &r) - { - for (auto i = r.pages().end(); i != r.pages().begin(); i--) - for (auto j = r.rows().end(); j != r.rows().begin(); j--) - for (auto k = r.cols().end(); k != r.cols().begin(); k--) - { - if ((i - 1 - m) % stride[0] == 0 && (j - 1 - n) % stride[1] == 0 && (k - 1 - p) % stride[2] == 0) - local_function(Array3i(i - 1, j - 1, k - 1)); - } - }, - ap); + parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) + { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) + { parallel_for((mesh_range.first)[2] + p, (mesh_range.second)[2], stride[2], [&](auto k) + { local_function(Array3i(i, j, k)); }, + ap); }, ap); }, ap); } } //=================================================================================================// diff --git a/src/shared/meshes/mesh_iterators.h b/src/shared/meshes/mesh_iterators.h index 1da19ee956..254d8a93ac 100644 --- a/src/shared/meshes/mesh_iterators.h +++ b/src/shared/meshes/mesh_iterators.h @@ -107,11 +107,17 @@ void mesh_for(const MeshRange &mesh_range, const LocalFunction &local_function, /** Iterator on the mesh by looping index. parallel computing. */ template void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_function, Args &&...args); -/** Iterator on the mesh by looping index with a stride. sequential computing. */ +/** Iterator on the mesh by looping index with a stride in forward direction. sequential computing. */ template -void mesh_split_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); -/** Iterator on the mesh by looping index with a stride. parallel computing. */ +void mesh_stride_forward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); +/** Iterator on the mesh by looping index with a stride in forward direction. parallel computing. */ template -void mesh_split_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); +void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); +/** Iterator on the mesh by looping index with a stride in backward direction. sequential computing. */ +template +void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); +/** Iterator on the mesh by looping index with a stride in backward direction. parallel computing. */ +template +void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); } // namespace SPH #endif // MESH_ITERATORS_H \ No newline at end of file diff --git a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp index 152e23f371..8ed0fed6a6 100644 --- a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp +++ b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp @@ -23,7 +23,7 @@ void function(const Array2i &index, std::vector &vec, const Array2i &all TEST(test_mesh_iterator, split_for) { - Array2i all_cells(9, 9); + Array2i all_cells(10, 10); MeshRange mesh_range(Array2i::Zero(), all_cells); Array2i stride = 3 * Array2i::Ones(); @@ -32,22 +32,39 @@ TEST(test_mesh_iterator, split_for) std::vector vec_1(total_cells); std::vector vec_2(total_cells); - mesh_split_for( + // forward test + mesh_stride_forward_for( + mesh_range, stride, + [&](const Array2i &index) + { + function(index, vec_2, all_cells); + }); + + mesh_stride_forward_parallel_for( mesh_range, stride, [&](const Array2i &index) { function(index, vec_1, all_cells); }); - mesh_split_parallel_for( + for (size_t i = 0; i < total_cells; ++i) + ASSERT_EQ(vec_1[i], vec_2[i]); + + // backward test + mesh_stride_backward_for( mesh_range, stride, [&](const Array2i &index) { function(index, vec_2, all_cells); }); + mesh_stride_backward_parallel_for( + mesh_range, stride, + [&](const Array2i &index) + { + function(index, vec_1, all_cells); + }); + for (size_t i = 0; i < total_cells; ++i) - { ASSERT_EQ(vec_1[i], vec_2[i]); - } } \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp index bd5ef27b2b..d590b7f60b 100644 --- a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp +++ b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp @@ -23,7 +23,7 @@ void function(const Array3i &index, std::vector &vec, const Array3i &all TEST(test_mesh_iterator, split_for) { - Array3i all_cells(9, 9, 9); + Array3i all_cells(10, 10, 10); MeshRange mesh_range(Array3i::Zero(), all_cells); Array3i stride = 3 * Array3i::Ones(); @@ -32,22 +32,39 @@ TEST(test_mesh_iterator, split_for) std::vector vec_1(total_cells); std::vector vec_2(total_cells); - mesh_split_for( + // forward test + mesh_stride_forward_for( + mesh_range, stride, + [&](const Array3i &index) + { + function(index, vec_2, all_cells); + }); + + mesh_stride_forward_parallel_for( mesh_range, stride, [&](const Array3i &index) { function(index, vec_1, all_cells); }); - mesh_split_parallel_for( + for (size_t i = 0; i < total_cells; ++i) + ASSERT_EQ(vec_1[i], vec_2[i]); + + // backward test + mesh_stride_backward_for( mesh_range, stride, [&](const Array3i &index) { function(index, vec_2, all_cells); }); + mesh_stride_backward_parallel_for( + mesh_range, stride, + [&](const Array3i &index) + { + function(index, vec_1, all_cells); + }); + for (size_t i = 0; i < total_cells; ++i) - { ASSERT_EQ(vec_1[i], vec_2[i]); - } } \ No newline at end of file From 226112d5266562ab7b947f372b760ab98aa763d5 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Thu, 12 Sep 2024 13:07:51 +0200 Subject: [PATCH 07/18] backward starts from end --- src/for_2D_build/meshes/mesh_iterators.hpp | 16 +++++++-------- src/for_3D_build/meshes/mesh_iterators.hpp | 24 +++++++++++----------- 2 files changed, 20 insertions(+), 20 deletions(-) diff --git a/src/for_2D_build/meshes/mesh_iterators.hpp b/src/for_2D_build/meshes/mesh_iterators.hpp index df4dec36cb..61a8424a57 100644 --- a/src/for_2D_build/meshes/mesh_iterators.hpp +++ b/src/for_2D_build/meshes/mesh_iterators.hpp @@ -112,10 +112,10 @@ void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi template void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - for (int m = stride[0] - 1; m >= 0; m--) - for (int n = stride[1] - 1; n >= 0; n--) - for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) - for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) + for (int m = stride[0]; m > 0; m--) + for (int n = stride[1]; n > 0; n--) + for (auto i = (mesh_range.first)[0] + m - 1; i < (mesh_range.second)[0]; i += stride[0]) + for (auto j = (mesh_range.first)[1] + n - 1; j < (mesh_range.second)[1]; j += stride[1]) { local_function(Array2i(i, j)); } @@ -124,11 +124,11 @@ void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, template void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - for (int m = stride[0] - 1; m >= 0; m--) - for (int n = stride[1] - 1; n >= 0; n--) + for (int m = stride[0]; m > 0; m--) + for (int n = stride[1]; n > 0; n--) { - parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) - { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) + parallel_for((mesh_range.first)[0] + m - 1, (mesh_range.second)[0], stride[0], [&](auto i) + { parallel_for((mesh_range.first)[1] + n - 1, (mesh_range.second)[1], stride[1], [&](auto j) { local_function(Array2i(i, j)); }, ap); }, ap); } diff --git a/src/for_3D_build/meshes/mesh_iterators.hpp b/src/for_3D_build/meshes/mesh_iterators.hpp index b9843c6bf8..391372470b 100644 --- a/src/for_3D_build/meshes/mesh_iterators.hpp +++ b/src/for_3D_build/meshes/mesh_iterators.hpp @@ -126,12 +126,12 @@ void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi template void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - for (int m = stride[0] - 1; m >= 0; m--) - for (int n = stride[1] - 1; n >= 0; n--) - for (int p = stride[2] - 1; p >= 0; p--) - for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) - for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) - for (auto k = (mesh_range.first)[2] + p; k < (mesh_range.second)[2]; k += stride[2]) + for (int m = stride[0]; m > 0; m--) + for (int n = stride[1]; n > 0; n--) + for (int p = stride[2]; p > 0; p--) + for (auto i = (mesh_range.first)[0] + m - 1; i < (mesh_range.second)[0]; i += stride[0]) + for (auto j = (mesh_range.first)[1] + n - 1; j < (mesh_range.second)[1]; j += stride[1]) + for (auto k = (mesh_range.first)[2] + p - 1; k < (mesh_range.second)[2]; k += stride[2]) { local_function(Array3i(i, j, k)); } @@ -140,13 +140,13 @@ void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, template void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) { - for (int m = stride[0] - 1; m >= 0; m--) - for (int n = stride[1] - 1; n >= 0; n--) - for (int p = stride[2] - 1; p >= 0; p--) + for (int m = stride[0]; m > 0; m--) + for (int n = stride[1]; n > 0; n--) + for (int p = stride[2]; p > 0; p--) { - parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) - { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) - { parallel_for((mesh_range.first)[2] + p, (mesh_range.second)[2], stride[2], [&](auto k) + parallel_for((mesh_range.first)[0] + m - 1, (mesh_range.second)[0], stride[0], [&](auto i) + { parallel_for((mesh_range.first)[1] + n - 1, (mesh_range.second)[1], stride[1], [&](auto j) + { parallel_for((mesh_range.first)[2] + p - 1, (mesh_range.second)[2], stride[2], [&](auto k) { local_function(Array3i(i, j, k)); }, ap); }, ap); }, ap); } From dd42e69b406ef0d82ed46da759b6c4d321d61b81 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Thu, 12 Sep 2024 13:08:11 +0200 Subject: [PATCH 08/18] modify particle_for_split --- src/shared/meshes/cell_linked_list.h | 7 ++-- src/shared/meshes/cell_linked_list.hpp | 50 ++++++++++++++++++++++++-- 2 files changed, 52 insertions(+), 5 deletions(-) diff --git a/src/shared/meshes/cell_linked_list.h b/src/shared/meshes/cell_linked_list.h index 838d2cfe85..e73c2760d2 100644 --- a/src/shared/meshes/cell_linked_list.h +++ b/src/shared/meshes/cell_linked_list.h @@ -33,6 +33,7 @@ #define MESH_CELL_LINKED_LIST_H #include "base_mesh.h" +#include "execution_policy.h" #include "neighborhood.h" namespace SPH @@ -119,7 +120,9 @@ class CellLinkedList : public BaseCellLinkedList, public Mesh /** split algorithm */; template - void particle_for_split(const LocalDynamicsFunction &local_dynamics_function); + void particle_for_split(const execution::SequencedPolicy &, const LocalDynamicsFunction &local_dynamics_function); + template + void particle_for_split(const execution::ParallelPolicy &, const LocalDynamicsFunction &local_dynamics_function); }; template <> @@ -153,7 +156,7 @@ class MultilevelCellLinkedList : public MultilevelMesh &computingSequence(BaseParticles &base_particles) override; virtual void tagBodyPartByCell(ConcurrentCellLists &cell_lists, std::function &check_included) override; - virtual void tagBoundingCells(StdVec &cell_data_lists, const BoundingBox &bounding_bounds, int axis) override{}; + virtual void tagBoundingCells(StdVec &cell_data_lists, const BoundingBox &bounding_bounds, int axis) override {}; virtual StdVec CellLinkedListLevels() override { return getMeshLevels(); }; }; } // namespace SPH diff --git a/src/shared/meshes/cell_linked_list.hpp b/src/shared/meshes/cell_linked_list.hpp index 5f3eaf195e..7a369e99f0 100644 --- a/src/shared/meshes/cell_linked_list.hpp +++ b/src/shared/meshes/cell_linked_list.hpp @@ -45,19 +45,63 @@ void CellLinkedList::searchNeighborsByParticles( } //=================================================================================================// template -void CellLinkedList::particle_for_split(const LocalDynamicsFunction &local_dynamics_function) +void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, const LocalDynamicsFunction &local_dynamics_function) { - mesh_split_parallel_for( + // fowward sweeping + mesh_stride_forward_for( MeshRange(Arrayi::Zero(), all_cells_), 3 * Arrayi::Ones(), [&](const Arrayi &cell_index) { const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); - for (size_t index_i : cell_list) + for (const size_t index_i : cell_list) { local_dynamics_function(index_i); } }); + + // backward sweeping + mesh_stride_backward_for( + MeshRange(Arrayi::Zero(), all_cells_), + 3 * Arrayi::Ones(), + [&](const Arrayi &cell_index) + { + const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + for (size_t i = cell_list.size(); i != 0; --i) + { + local_dynamics_function(cell_list[i - 1]); + } + }); +} +//=================================================================================================// +template +void CellLinkedList::particle_for_split(const execution::ParallelPolicy &, const LocalDynamicsFunction &local_dynamics_function) +{ + // foward sweeping + mesh_stride_forward_parallel_for( + MeshRange(Arrayi::Zero(), all_cells_), + 3 * Arrayi::Ones(), + [&](const Arrayi &cell_index) + { + const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + for (const size_t index_i : cell_list) + { + local_dynamics_function(index_i); + } + }); + + // backward sweeping + mesh_stride_backward_parallel_for( + MeshRange(Arrayi::Zero(), all_cells_), + 3 * Arrayi::Ones(), + [&](const Arrayi &cell_index) + { + const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + for (size_t i = cell_list.size(); i != 0; --i) + { + local_dynamics_function(cell_list[i - 1]); + } + }); } //=================================================================================================// } // namespace SPH From 26311ae6f0fc5f89ac02b0b48c850c6781626f07 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Thu, 12 Sep 2024 13:08:51 +0200 Subject: [PATCH 09/18] modify interactionsplit --- .../particle_dynamics/particle_dynamics_algorithms.h | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/src/shared/particle_dynamics/particle_dynamics_algorithms.h b/src/shared/particle_dynamics/particle_dynamics_algorithms.h index 8964f5ff50..2856fb5d05 100644 --- a/src/shared/particle_dynamics/particle_dynamics_algorithms.h +++ b/src/shared/particle_dynamics/particle_dynamics_algorithms.h @@ -203,7 +203,7 @@ class InteractionSplit : public BaseInteractionDynamics - InteractionSplit(Args &&... args) + explicit InteractionSplit(Args &&...args) : BaseInteractionDynamics(std::forward(args)...), real_body_(DynamicCast(this, this->getSPHBody())), cell_linked_list_(DynamicCast(this, real_body_.getCellLinkedList())) @@ -212,12 +212,11 @@ class InteractionSplit : public BaseInteractionDynamics::value, "LocalDynamicsType does not fulfill InteractionSplit requirements"); }; - virtual ~InteractionSplit(){}; /** run the main interaction step between particles. */ - virtual void runMainStep(Real dt) override + void runMainStep(Real dt) override { - cell_linked_list_.particle_for_split([&](size_t i) + cell_linked_list_.particle_for_split(ExecutionPolicy(), [&](size_t i) { this->interaction(i, dt * 0.5); }); } }; From ec60fc06d4167f67be9db79006abdab8ecf34cdc Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 14:31:55 +0200 Subject: [PATCH 10/18] delete tests --- .../for_2D_build/iterators/CMakeLists.txt | 7 -- .../test_2d_mesh_iterator/CMakeLists.txt | 15 ---- .../test_mesh_iterator.cpp | 70 ------------------- .../for_3D_build/iterators/CMakeLists.txt | 7 -- .../test_3d_mesh_iterator/CMakeLists.txt | 16 ----- .../test_mesh_iterator.cpp | 70 ------------------- 6 files changed, 185 deletions(-) delete mode 100644 tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt delete mode 100644 tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt delete mode 100644 tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp delete mode 100644 tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt delete mode 100644 tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt delete mode 100644 tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp diff --git a/tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt b/tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt deleted file mode 100644 index 8eb174fa72..0000000000 --- a/tests/unit_tests_src/for_2D_build/iterators/CMakeLists.txt +++ /dev/null @@ -1,7 +0,0 @@ -SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR}) - -foreach(subdir ${SUBDIRS}) - if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt) - add_subdirectory(${subdir}) - endif() -endforeach() \ No newline at end of file diff --git a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt deleted file mode 100644 index c5b091dac7..0000000000 --- a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/CMakeLists.txt +++ /dev/null @@ -1,15 +0,0 @@ -STRING( REGEX REPLACE ".*/(.*)" "\\1" CURRENT_FOLDER ${CMAKE_CURRENT_SOURCE_DIR} ) -PROJECT("${CURRENT_FOLDER}") - -SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib) -SET(EXECUTABLE_OUTPUT_PATH "${PROJECT_BINARY_DIR}/bin/") -SET(BUILD_INPUT_PATH "${EXECUTABLE_OUTPUT_PATH}/input") -SET(BUILD_RELOAD_PATH "${EXECUTABLE_OUTPUT_PATH}/reload") - -aux_source_directory(. DIR_SRCS) -ADD_EXECUTABLE(${PROJECT_NAME} ${EXECUTABLE_OUTPUT_PATH} ${DIR_SRCS}) -target_link_libraries(${PROJECT_NAME} sphinxsys_2d GTest::gtest GTest::gtest_main) -set_target_properties(${PROJECT_NAME} PROPERTIES VS_DEBUGGER_WORKING_DIRECTORY "${EXECUTABLE_OUTPUT_PATH}") - -add_test(NAME ${PROJECT_NAME} COMMAND ${PROJECT_NAME} - WORKING_DIRECTORY ${EXECUTABLE_OUTPUT_PATH}) diff --git a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp b/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp deleted file mode 100644 index 8ed0fed6a6..0000000000 --- a/tests/unit_tests_src/for_2D_build/iterators/test_2d_mesh_iterator/test_mesh_iterator.cpp +++ /dev/null @@ -1,70 +0,0 @@ -#include "base_data_type.h" -#include "mesh_iterators.hpp" -#include - -using namespace SPH; - -size_t get_cell_index(const Array2i &index) -{ - return 3 * (index[0] % 3) + (index[1] % 3); -} - -size_t get_1d_index(const Array2i &all_cells, const Array2i &index) -{ - return all_cells[1] * index[0] + index[1]; -} - -void function(const Array2i &index, std::vector &vec, const Array2i &all_cells) -{ - size_t index_cell = get_cell_index(index); - size_t index_1d = get_1d_index(all_cells, index); - vec[index_1d] += index_cell; -} - -TEST(test_mesh_iterator, split_for) -{ - Array2i all_cells(10, 10); - MeshRange mesh_range(Array2i::Zero(), all_cells); - Array2i stride = 3 * Array2i::Ones(); - - const size_t total_cells = all_cells[0] * all_cells[1]; - - std::vector vec_1(total_cells); - std::vector vec_2(total_cells); - - // forward test - mesh_stride_forward_for( - mesh_range, stride, - [&](const Array2i &index) - { - function(index, vec_2, all_cells); - }); - - mesh_stride_forward_parallel_for( - mesh_range, stride, - [&](const Array2i &index) - { - function(index, vec_1, all_cells); - }); - - for (size_t i = 0; i < total_cells; ++i) - ASSERT_EQ(vec_1[i], vec_2[i]); - - // backward test - mesh_stride_backward_for( - mesh_range, stride, - [&](const Array2i &index) - { - function(index, vec_2, all_cells); - }); - - mesh_stride_backward_parallel_for( - mesh_range, stride, - [&](const Array2i &index) - { - function(index, vec_1, all_cells); - }); - - for (size_t i = 0; i < total_cells; ++i) - ASSERT_EQ(vec_1[i], vec_2[i]); -} \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt b/tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt deleted file mode 100644 index 8eb174fa72..0000000000 --- a/tests/unit_tests_src/for_3D_build/iterators/CMakeLists.txt +++ /dev/null @@ -1,7 +0,0 @@ -SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR}) - -foreach(subdir ${SUBDIRS}) - if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt) - add_subdirectory(${subdir}) - endif() -endforeach() \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt deleted file mode 100644 index d9a298a2e7..0000000000 --- a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/CMakeLists.txt +++ /dev/null @@ -1,16 +0,0 @@ -STRING( REGEX REPLACE ".*/(.*)" "\\1" CURRENT_FOLDER ${CMAKE_CURRENT_SOURCE_DIR} ) -PROJECT("${CURRENT_FOLDER}") - -SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib) -SET(EXECUTABLE_OUTPUT_PATH "${PROJECT_BINARY_DIR}/bin/") -SET(BUILD_INPUT_PATH "${EXECUTABLE_OUTPUT_PATH}/input") -SET(BUILD_RELOAD_PATH "${EXECUTABLE_OUTPUT_PATH}/reload") - -aux_source_directory(. DIR_SRCS) -ADD_EXECUTABLE(${PROJECT_NAME} ${EXECUTABLE_OUTPUT_PATH} ${DIR_SRCS}) -target_link_libraries(${PROJECT_NAME} sphinxsys_3d GTest::gtest GTest::gtest_main) -set_target_properties(${PROJECT_NAME} PROPERTIES VS_DEBUGGER_WORKING_DIRECTORY "${EXECUTABLE_OUTPUT_PATH}") - -add_test(NAME ${PROJECT_NAME}_particle_relaxation - COMMAND ${PROJECT_NAME} --relax=true - WORKING_DIRECTORY ${EXECUTABLE_OUTPUT_PATH}) diff --git a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp b/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp deleted file mode 100644 index d590b7f60b..0000000000 --- a/tests/unit_tests_src/for_3D_build/iterators/test_3d_mesh_iterator/test_mesh_iterator.cpp +++ /dev/null @@ -1,70 +0,0 @@ -#include "base_data_type.h" -#include "mesh_iterators.hpp" -#include - -using namespace SPH; - -size_t get_cell_index(const Array3i &index) -{ - return 9 * (index[0] % 3) + 3 * (index[1] % 3) + (index[2] % 3); -} - -size_t get_1d_index(const Array3i &all_cells, const Array3i &index) -{ - return all_cells[1] * all_cells[2] * index[0] + all_cells[2] * index[1] + index[2]; -} - -void function(const Array3i &index, std::vector &vec, const Array3i &all_cells) -{ - size_t index_cell = get_cell_index(index); - size_t index_1d = get_1d_index(all_cells, index); - vec[index_1d] += index_cell; -} - -TEST(test_mesh_iterator, split_for) -{ - Array3i all_cells(10, 10, 10); - MeshRange mesh_range(Array3i::Zero(), all_cells); - Array3i stride = 3 * Array3i::Ones(); - - const size_t total_cells = all_cells[0] * all_cells[1] * all_cells[2]; - - std::vector vec_1(total_cells); - std::vector vec_2(total_cells); - - // forward test - mesh_stride_forward_for( - mesh_range, stride, - [&](const Array3i &index) - { - function(index, vec_2, all_cells); - }); - - mesh_stride_forward_parallel_for( - mesh_range, stride, - [&](const Array3i &index) - { - function(index, vec_1, all_cells); - }); - - for (size_t i = 0; i < total_cells; ++i) - ASSERT_EQ(vec_1[i], vec_2[i]); - - // backward test - mesh_stride_backward_for( - mesh_range, stride, - [&](const Array3i &index) - { - function(index, vec_2, all_cells); - }); - - mesh_stride_backward_parallel_for( - mesh_range, stride, - [&](const Array3i &index) - { - function(index, vec_1, all_cells); - }); - - for (size_t i = 0; i < total_cells; ++i) - ASSERT_EQ(vec_1[i], vec_2[i]); -} \ No newline at end of file From b398a94c2446b2e97e39c246d05346791778e91a Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 14:32:03 +0200 Subject: [PATCH 11/18] delete cmake preset --- CMakePresets.json | 27 --------------------------- 1 file changed, 27 deletions(-) delete mode 100644 CMakePresets.json diff --git a/CMakePresets.json b/CMakePresets.json deleted file mode 100644 index 97aa089298..0000000000 --- a/CMakePresets.json +++ /dev/null @@ -1,27 +0,0 @@ -{ - "version": 8, - "configurePresets": [ - { - "name": "default", - "displayName": "SPHinXsys Debug", - "description": "Sets Ninja generator, build and install directory", - "generator": "Ninja", - "binaryDir": "${sourceDir}/build", - "cacheVariables": { - "CMAKE_BUILD_TYPE": "Debug", - "CMAKE_TOOLCHAIN_FILE": "'$/home/xyhu/vcpkg/scripts/buildsystems/vcpkg.cmake'" - } - }, - { - "name": "release", - "displayName": "SPHinXsys Release", - "description": "Sets Ninja generator, build and install directory", - "generator": "Ninja", - "binaryDir": "${sourceDir}/build", - "cacheVariables": { - "CMAKE_BUILD_TYPE": "Release", - "CMAKE_TOOLCHAIN_FILE": "'$/home/xyhu/vcpkg/scripts/buildsystems/vcpkg.cmake'" - } - } - ] -} \ No newline at end of file From 74f2fe78f8f2285193faa83ffc059ad358eae095 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 14:32:20 +0200 Subject: [PATCH 12/18] delete functions in mesh iterator --- src/for_2D_build/meshes/mesh_iterators.hpp | 50 ------------------- src/for_3D_build/meshes/mesh_iterators.hpp | 58 ---------------------- src/shared/meshes/mesh_iterators.h | 13 +---- 3 files changed, 1 insertion(+), 120 deletions(-) diff --git a/src/for_2D_build/meshes/mesh_iterators.hpp b/src/for_2D_build/meshes/mesh_iterators.hpp index 61a8424a57..57bd4d993c 100644 --- a/src/for_2D_build/meshes/mesh_iterators.hpp +++ b/src/for_2D_build/meshes/mesh_iterators.hpp @@ -84,54 +84,4 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f ap); } //=================================================================================================// -template -void mesh_stride_forward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = 0; m < stride[0]; m++) - for (int n = 0; n < stride[1]; n++) - for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) - for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) - { - local_function(Array2i(i, j)); - } -} -//=================================================================================================// -template -void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = 0; m < stride[0]; m++) - for (int n = 0; n < stride[1]; n++) - { - parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) - { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) - { local_function(Array2i(i, j)); }, ap); }, - ap); - } -} -//=================================================================================================// -template -void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = stride[0]; m > 0; m--) - for (int n = stride[1]; n > 0; n--) - for (auto i = (mesh_range.first)[0] + m - 1; i < (mesh_range.second)[0]; i += stride[0]) - for (auto j = (mesh_range.first)[1] + n - 1; j < (mesh_range.second)[1]; j += stride[1]) - { - local_function(Array2i(i, j)); - } -} -//=================================================================================================// -template -void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = stride[0]; m > 0; m--) - for (int n = stride[1]; n > 0; n--) - { - parallel_for((mesh_range.first)[0] + m - 1, (mesh_range.second)[0], stride[0], [&](auto i) - { parallel_for((mesh_range.first)[1] + n - 1, (mesh_range.second)[1], stride[1], [&](auto j) - { local_function(Array2i(i, j)); }, ap); }, - ap); - } -} -//=================================================================================================// } // namespace SPH diff --git a/src/for_3D_build/meshes/mesh_iterators.hpp b/src/for_3D_build/meshes/mesh_iterators.hpp index 391372470b..4221a3728c 100644 --- a/src/for_3D_build/meshes/mesh_iterators.hpp +++ b/src/for_3D_build/meshes/mesh_iterators.hpp @@ -94,62 +94,4 @@ void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_f ap); } //=================================================================================================// -template -void mesh_stride_forward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = 0; m < stride[0]; m++) - for (int n = 0; n < stride[1]; n++) - for (int p = 0; p < stride[2]; p++) - for (auto i = (mesh_range.first)[0] + m; i < (mesh_range.second)[0]; i += stride[0]) - for (auto j = (mesh_range.first)[1] + n; j < (mesh_range.second)[1]; j += stride[1]) - for (auto k = (mesh_range.first)[2] + p; k < (mesh_range.second)[2]; k += stride[2]) - { - local_function(Array3i(i, j, k)); - } -} -//=================================================================================================// -template -void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = 0; m < stride[0]; m++) - for (int n = 0; n < stride[1]; n++) - for (int p = 0; p < stride[2]; p++) - { - parallel_for((mesh_range.first)[0] + m, (mesh_range.second)[0], stride[0], [&](auto i) - { parallel_for((mesh_range.first)[1] + n, (mesh_range.second)[1], stride[1], [&](auto j) - { parallel_for((mesh_range.first)[2] + p, (mesh_range.second)[2], stride[2], [&](auto k) - { local_function(Array3i(i, j, k)); }, - ap); }, ap); }, ap); - } -} -//=================================================================================================// -template -void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = stride[0]; m > 0; m--) - for (int n = stride[1]; n > 0; n--) - for (int p = stride[2]; p > 0; p--) - for (auto i = (mesh_range.first)[0] + m - 1; i < (mesh_range.second)[0]; i += stride[0]) - for (auto j = (mesh_range.first)[1] + n - 1; j < (mesh_range.second)[1]; j += stride[1]) - for (auto k = (mesh_range.first)[2] + p - 1; k < (mesh_range.second)[2]; k += stride[2]) - { - local_function(Array3i(i, j, k)); - } -} -//=================================================================================================// -template -void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args) -{ - for (int m = stride[0]; m > 0; m--) - for (int n = stride[1]; n > 0; n--) - for (int p = stride[2]; p > 0; p--) - { - parallel_for((mesh_range.first)[0] + m - 1, (mesh_range.second)[0], stride[0], [&](auto i) - { parallel_for((mesh_range.first)[1] + n - 1, (mesh_range.second)[1], stride[1], [&](auto j) - { parallel_for((mesh_range.first)[2] + p - 1, (mesh_range.second)[2], stride[2], [&](auto k) - { local_function(Array3i(i, j, k)); }, - ap); }, ap); }, ap); - } -} -//=================================================================================================// } // namespace SPH diff --git a/src/shared/meshes/mesh_iterators.h b/src/shared/meshes/mesh_iterators.h index 254d8a93ac..5f3c29e869 100644 --- a/src/shared/meshes/mesh_iterators.h +++ b/src/shared/meshes/mesh_iterators.h @@ -33,6 +33,7 @@ #define MESH_ITERATORS_H #include "base_data_package.h" +#include "base_mesh.h" namespace SPH { @@ -107,17 +108,5 @@ void mesh_for(const MeshRange &mesh_range, const LocalFunction &local_function, /** Iterator on the mesh by looping index. parallel computing. */ template void mesh_parallel_for(const MeshRange &mesh_range, const LocalFunction &local_function, Args &&...args); -/** Iterator on the mesh by looping index with a stride in forward direction. sequential computing. */ -template -void mesh_stride_forward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); -/** Iterator on the mesh by looping index with a stride in forward direction. parallel computing. */ -template -void mesh_stride_forward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); -/** Iterator on the mesh by looping index with a stride in backward direction. sequential computing. */ -template -void mesh_stride_backward_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); -/** Iterator on the mesh by looping index with a stride in backward direction. parallel computing. */ -template -void mesh_stride_backward_parallel_for(const MeshRange &mesh_range, const Arrayi &stride, const LocalFunction &local_function, Args &&...args); } // namespace SPH #endif // MESH_ITERATORS_H \ No newline at end of file From 6ca6e57bac2385c16602566e2fb0c7210569858f Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 14:32:48 +0200 Subject: [PATCH 13/18] add function to get M*N --- src/for_2D_build/meshes/base_mesh_2d.cpp | 5 +++++ src/for_3D_build/meshes/base_mesh_3d.cpp | 5 +++++ src/shared/meshes/base_mesh.h | 2 ++ 3 files changed, 12 insertions(+) diff --git a/src/for_2D_build/meshes/base_mesh_2d.cpp b/src/for_2D_build/meshes/base_mesh_2d.cpp index 3985460ab7..78571bc16b 100644 --- a/src/for_2D_build/meshes/base_mesh_2d.cpp +++ b/src/for_2D_build/meshes/base_mesh_2d.cpp @@ -15,6 +15,11 @@ size_t BaseMesh::transferMeshIndexTo1D(const Array2i &mesh_size, const Array2i & return mesh_index[0] * mesh_size[1] + mesh_index[1]; } //=============================================================================================// +size_t BaseMesh::get1DMeshSize(const Array2i &mesh_size) +{ + return mesh_size[0] * mesh_size[1]; +} +//=============================================================================================// size_t BaseMesh::transferMeshIndexToMortonOrder(const Array2i &mesh_index) { return MortonCode(mesh_index[0]) | (MortonCode(mesh_index[1]) << 1); diff --git a/src/for_3D_build/meshes/base_mesh_3d.cpp b/src/for_3D_build/meshes/base_mesh_3d.cpp index d9902fcfe8..6bf2008f93 100644 --- a/src/for_3D_build/meshes/base_mesh_3d.cpp +++ b/src/for_3D_build/meshes/base_mesh_3d.cpp @@ -20,6 +20,11 @@ size_t BaseMesh::transferMeshIndexTo1D(const Array3i &mesh_size, const Array3i & mesh_index[2]; } //=================================================================================================// +size_t BaseMesh::get1DMeshSize(const Array3i &mesh_size) +{ + return mesh_size[0] * mesh_size[1] * mesh_size[2]; +} +//=================================================================================================// size_t BaseMesh::transferMeshIndexToMortonOrder(const Array3i &mesh_index) { return MortonCode(mesh_index[0]) | (MortonCode(mesh_index[1]) << 1) | (MortonCode(mesh_index[2]) << 2); diff --git a/src/shared/meshes/base_mesh.h b/src/shared/meshes/base_mesh.h index a5b3a04681..97a0b54b5e 100644 --- a/src/shared/meshes/base_mesh.h +++ b/src/shared/meshes/base_mesh.h @@ -90,6 +90,8 @@ class BaseMesh Arrayi transfer1DtoMeshIndex(const Arrayi &mesh_size, size_t i); /** Transfer mesh index to 1D int. */ size_t transferMeshIndexTo1D(const Arrayi &mesh_size, const Arrayi &mesh_index); + /** Get 1D mesh size from 2D/3D mesh size */ + size_t get1DMeshSize(const Arrayi &mesh_size); /** converts mesh index into a Morton order. * Interleave a 10 bit number in 32 bits, fill one bit and leave the other 2 as zeros * https://stackoverflow.com/questions/18529057/ From c34a9294ce1384dcdc139884994c50b4243601ea Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 14:33:05 +0200 Subject: [PATCH 14/18] modify split algorithm --- src/shared/meshes/cell_linked_list.cpp | 2 +- src/shared/meshes/cell_linked_list.h | 2 + src/shared/meshes/cell_linked_list.hpp | 94 +++++++++++++++++--------- 3 files changed, 66 insertions(+), 32 deletions(-) diff --git a/src/shared/meshes/cell_linked_list.cpp b/src/shared/meshes/cell_linked_list.cpp index 837ef78401..6099f86bb2 100644 --- a/src/shared/meshes/cell_linked_list.cpp +++ b/src/shared/meshes/cell_linked_list.cpp @@ -17,7 +17,7 @@ BaseCellLinkedList:: CellLinkedList::CellLinkedList(BoundingBox tentative_bounds, Real grid_spacing, SPHAdaptation &sph_adaptation) : BaseCellLinkedList(sph_adaptation), Mesh(tentative_bounds, grid_spacing, 2), - cell_index_lists_(nullptr), cell_data_lists_(nullptr) + cell_index_lists_(nullptr), cell_data_lists_(nullptr), number_of_split_cell_lists_(pow(3, Dimensions)) { allocateMeshDataMatrix(); single_cell_linked_list_level_.push_back(this); diff --git a/src/shared/meshes/cell_linked_list.h b/src/shared/meshes/cell_linked_list.h index e73c2760d2..8349315d07 100644 --- a/src/shared/meshes/cell_linked_list.h +++ b/src/shared/meshes/cell_linked_list.h @@ -88,6 +88,8 @@ class CellLinkedList : public BaseCellLinkedList, public Mesh ConcurrentIndexVector *cell_index_lists_; /** non-concurrent list data rewritten for building neighbor list */ ListDataVector *cell_data_lists_; + /**< number of split cell lists */ + size_t number_of_split_cell_lists_; void allocateMeshDataMatrix(); /**< allocate memories for addresses of data packages. */ void deleteMeshDataMatrix(); /**< delete memories for addresses of data packages. */ diff --git a/src/shared/meshes/cell_linked_list.hpp b/src/shared/meshes/cell_linked_list.hpp index 7a369e99f0..fb832b3001 100644 --- a/src/shared/meshes/cell_linked_list.hpp +++ b/src/shared/meshes/cell_linked_list.hpp @@ -47,61 +47,93 @@ void CellLinkedList::searchNeighborsByParticles( template void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, const LocalDynamicsFunction &local_dynamics_function) { - // fowward sweeping - mesh_stride_forward_for( - MeshRange(Arrayi::Zero(), all_cells_), - 3 * Arrayi::Ones(), - [&](const Arrayi &cell_index) + // foward sweeping + for (size_t k = 0; k < number_of_split_cell_lists_; k++) + { + const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k); + const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); + const size_t number_of_cells = get1DMeshSize(all_cells_k); + + for (size_t l = 0; l < number_of_cells; l++) { + const Arrayi cell_index = split_cell_index + 3 * transfer1DtoMeshIndex(all_cells_k, l); const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); for (const size_t index_i : cell_list) { local_dynamics_function(index_i); } - }); + } + } // backward sweeping - mesh_stride_backward_for( - MeshRange(Arrayi::Zero(), all_cells_), - 3 * Arrayi::Ones(), - [&](const Arrayi &cell_index) + for (size_t k = number_of_split_cell_lists_; k != 0; --k) + { + const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k - 1); + const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); + const size_t number_of_cells = get1DMeshSize(all_cells_k); + + for (size_t l = 0; l < number_of_cells; l++) { + const Arrayi cell_index = split_cell_index + 3 * transfer1DtoMeshIndex(all_cells_k, l); const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); for (size_t i = cell_list.size(); i != 0; --i) { local_dynamics_function(cell_list[i - 1]); } - }); + } + } } //=================================================================================================// template void CellLinkedList::particle_for_split(const execution::ParallelPolicy &, const LocalDynamicsFunction &local_dynamics_function) { // foward sweeping - mesh_stride_forward_parallel_for( - MeshRange(Arrayi::Zero(), all_cells_), - 3 * Arrayi::Ones(), - [&](const Arrayi &cell_index) - { - const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); - for (const size_t index_i : cell_list) + for (size_t k = 0; k < number_of_split_cell_lists_; k++) + { + const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k); + const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); + const size_t number_of_cells = get1DMeshSize(all_cells_k); + + parallel_for( + IndexRange(0, number_of_cells), + [&](const IndexRange &r) { - local_dynamics_function(index_i); - } - }); + for (size_t l = r.begin(); l < r.end(); ++l) + { + const Arrayi cell_index = split_cell_index + 3 * transfer1DtoMeshIndex(all_cells_k, l); + const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + for (const size_t index_i : cell_list) + { + local_dynamics_function(index_i); + } + } + }, + ap); + } // backward sweeping - mesh_stride_backward_parallel_for( - MeshRange(Arrayi::Zero(), all_cells_), - 3 * Arrayi::Ones(), - [&](const Arrayi &cell_index) - { - const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); - for (size_t i = cell_list.size(); i != 0; --i) + for (size_t k = number_of_split_cell_lists_; k != 0; --k) + { + const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k - 1); + const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); + const size_t number_of_cells = get1DMeshSize(all_cells_k); + + parallel_for( + IndexRange(0, number_of_cells), + [&](const IndexRange &r) { - local_dynamics_function(cell_list[i - 1]); - } - }); + for (size_t l = r.begin(); l < r.end(); ++l) + { + const Arrayi cell_index = split_cell_index + 3 * transfer1DtoMeshIndex(all_cells_k, l); + const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + for (size_t i = cell_list.size(); i != 0; --i) + { + local_dynamics_function(cell_list[i - 1]); + } + } + }, + ap); + } } //=================================================================================================// } // namespace SPH From 6c487f20bfaadc8fa50ef556214748022c019bbb Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 16:23:05 +0200 Subject: [PATCH 15/18] add unit test --- .../for_2D_build/meshes/CMakeLists.txt | 7 ++ .../test_2d_particle_for_split/CMakeLists.txt | 15 +++++ .../test_particle_for_split.cpp | 64 +++++++++++++++++++ .../for_3D_build/meshes/CMakeLists.txt | 7 ++ .../test_3d_particle_for_split/CMakeLists.txt | 16 +++++ .../test_particle_for_split.cpp | 62 ++++++++++++++++++ 6 files changed, 171 insertions(+) create mode 100644 tests/unit_tests_src/for_2D_build/meshes/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp create mode 100644 tests/unit_tests_src/for_3D_build/meshes/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/CMakeLists.txt create mode 100644 tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp diff --git a/tests/unit_tests_src/for_2D_build/meshes/CMakeLists.txt b/tests/unit_tests_src/for_2D_build/meshes/CMakeLists.txt new file mode 100644 index 0000000000..8eb174fa72 --- /dev/null +++ b/tests/unit_tests_src/for_2D_build/meshes/CMakeLists.txt @@ -0,0 +1,7 @@ +SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR}) + +foreach(subdir ${SUBDIRS}) + if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt) + add_subdirectory(${subdir}) + endif() +endforeach() \ No newline at end of file diff --git a/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/CMakeLists.txt b/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/CMakeLists.txt new file mode 100644 index 0000000000..c5b091dac7 --- /dev/null +++ b/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/CMakeLists.txt @@ -0,0 +1,15 @@ +STRING( REGEX REPLACE ".*/(.*)" "\\1" CURRENT_FOLDER ${CMAKE_CURRENT_SOURCE_DIR} ) +PROJECT("${CURRENT_FOLDER}") + +SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib) +SET(EXECUTABLE_OUTPUT_PATH "${PROJECT_BINARY_DIR}/bin/") +SET(BUILD_INPUT_PATH "${EXECUTABLE_OUTPUT_PATH}/input") +SET(BUILD_RELOAD_PATH "${EXECUTABLE_OUTPUT_PATH}/reload") + +aux_source_directory(. DIR_SRCS) +ADD_EXECUTABLE(${PROJECT_NAME} ${EXECUTABLE_OUTPUT_PATH} ${DIR_SRCS}) +target_link_libraries(${PROJECT_NAME} sphinxsys_2d GTest::gtest GTest::gtest_main) +set_target_properties(${PROJECT_NAME} PROPERTIES VS_DEBUGGER_WORKING_DIRECTORY "${EXECUTABLE_OUTPUT_PATH}") + +add_test(NAME ${PROJECT_NAME} COMMAND ${PROJECT_NAME} + WORKING_DIRECTORY ${EXECUTABLE_OUTPUT_PATH}) diff --git a/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp b/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp new file mode 100644 index 0000000000..d1fe326df9 --- /dev/null +++ b/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp @@ -0,0 +1,64 @@ +#include "sphinxsys.h" +#include + +using namespace SPH; + +TEST(test_meshes, split_for) +{ + Real length = 10; + Real dp = 1; + + MultiPolygon shape; + shape.addABox(Transform(0.5 * length * Vec2d::Ones()), 0.5 * length * Vec2d::Ones(), ShapeBooleanOps::add); + auto polygon_shape = makeShared(shape, "PolygonShape"); + + BoundingBox bb_system = polygon_shape->getBounds(); + + SPHSystem system(bb_system, dp); + + SolidBody body(system, polygon_shape); + body.defineMaterial(); + body.generateParticles(); + auto &particles = body.getBaseParticles(); + const auto &pos = *particles.registerSharedVariable("Position"); + auto &quantity = *particles.registerSharedVariable("Quantity", [&](size_t i) + { return pos[i].norm(); }); + + InnerRelation inner(body); + + body.updateCellLinkedList(); + inner.updateConfiguration(); + + auto interaction = [&](size_t index_i) + { + const auto &configuration = inner.inner_configuration_; + const auto &neighborhood = configuration[index_i]; + for (size_t n = 0; n < neighborhood.current_size_; n++) + { + size_t q_ij = quantity[index_i] - quantity[neighborhood.j_[n]]; + quantity[index_i] += 0.5 * q_ij; + quantity[neighborhood.j_[n]] -= 0.5 * q_ij; + } + }; + + auto &cell_linked_list = *dynamic_cast(&body.getCellLinkedList()); + + // run the interaction in sequenced policy + cell_linked_list.particle_for_split(execution::SequencedPolicy(), interaction); + // record the result + auto q_seq = quantity; + + // reset the value + for (size_t i = 0; i < particles.TotalRealParticles(); i++) + { + quantity[i] = pos[i].norm(); + } + // run the interaction in parallel policy + cell_linked_list.particle_for_split(execution::ParallelPolicy(), interaction); + + // check the result + for (size_t i = 0; i < particles.TotalRealParticles(); i++) + { + ASSERT_EQ(q_seq[i], quantity[i]); + } +} \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/meshes/CMakeLists.txt b/tests/unit_tests_src/for_3D_build/meshes/CMakeLists.txt new file mode 100644 index 0000000000..8eb174fa72 --- /dev/null +++ b/tests/unit_tests_src/for_3D_build/meshes/CMakeLists.txt @@ -0,0 +1,7 @@ +SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR}) + +foreach(subdir ${SUBDIRS}) + if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt) + add_subdirectory(${subdir}) + endif() +endforeach() \ No newline at end of file diff --git a/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/CMakeLists.txt b/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/CMakeLists.txt new file mode 100644 index 0000000000..d9a298a2e7 --- /dev/null +++ b/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/CMakeLists.txt @@ -0,0 +1,16 @@ +STRING( REGEX REPLACE ".*/(.*)" "\\1" CURRENT_FOLDER ${CMAKE_CURRENT_SOURCE_DIR} ) +PROJECT("${CURRENT_FOLDER}") + +SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib) +SET(EXECUTABLE_OUTPUT_PATH "${PROJECT_BINARY_DIR}/bin/") +SET(BUILD_INPUT_PATH "${EXECUTABLE_OUTPUT_PATH}/input") +SET(BUILD_RELOAD_PATH "${EXECUTABLE_OUTPUT_PATH}/reload") + +aux_source_directory(. DIR_SRCS) +ADD_EXECUTABLE(${PROJECT_NAME} ${EXECUTABLE_OUTPUT_PATH} ${DIR_SRCS}) +target_link_libraries(${PROJECT_NAME} sphinxsys_3d GTest::gtest GTest::gtest_main) +set_target_properties(${PROJECT_NAME} PROPERTIES VS_DEBUGGER_WORKING_DIRECTORY "${EXECUTABLE_OUTPUT_PATH}") + +add_test(NAME ${PROJECT_NAME}_particle_relaxation + COMMAND ${PROJECT_NAME} --relax=true + WORKING_DIRECTORY ${EXECUTABLE_OUTPUT_PATH}) diff --git a/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp b/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp new file mode 100644 index 0000000000..fe608b7699 --- /dev/null +++ b/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp @@ -0,0 +1,62 @@ +#include "sphinxsys.h" +#include + +using namespace SPH; + +TEST(test_meshes, split_for) +{ + Real length = 10; + Real dp = 1; + + auto shape = makeShared(0.5 * length * Vec3d::Ones(), "Shape"); + + BoundingBox bb_system = shape->getBounds(); + + SPHSystem system(bb_system, dp); + + SolidBody body(system, shape); + body.defineMaterial(); + body.generateParticles(); + auto &particles = body.getBaseParticles(); + const auto &pos = *particles.registerSharedVariable("Position"); + auto &quantity = *particles.registerSharedVariable("Quantity", [&](size_t i) + { return pos[i].norm(); }); + + InnerRelation inner(body); + + body.updateCellLinkedList(); + inner.updateConfiguration(); + + auto interaction = [&](size_t index_i) + { + const auto &configuration = inner.inner_configuration_; + const auto &neighborhood = configuration[index_i]; + for (size_t n = 0; n < neighborhood.current_size_; n++) + { + size_t q_ij = quantity[index_i] - quantity[neighborhood.j_[n]]; + quantity[index_i] += 0.5 * q_ij; + quantity[neighborhood.j_[n]] -= 0.5 * q_ij; + } + }; + + auto &cell_linked_list = *dynamic_cast(&body.getCellLinkedList()); + + // run the interaction in sequenced policy + cell_linked_list.particle_for_split(execution::SequencedPolicy(), interaction); + // record the result + auto q_seq = quantity; + + // reset the value + for (size_t i = 0; i < particles.TotalRealParticles(); i++) + { + quantity[i] = pos[i].norm(); + } + // run the interaction in parallel policy + cell_linked_list.particle_for_split(execution::ParallelPolicy(), interaction); + + // check the result + for (size_t i = 0; i < particles.TotalRealParticles(); i++) + { + ASSERT_EQ(q_seq[i], quantity[i]); + } +} \ No newline at end of file From e122e6304667bb0361808867b2714530080e3211 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Fri, 13 Sep 2024 17:47:09 +0200 Subject: [PATCH 16/18] add comments --- src/shared/meshes/cell_linked_list.hpp | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/src/shared/meshes/cell_linked_list.hpp b/src/shared/meshes/cell_linked_list.hpp index fb832b3001..b15148bcd1 100644 --- a/src/shared/meshes/cell_linked_list.hpp +++ b/src/shared/meshes/cell_linked_list.hpp @@ -50,14 +50,26 @@ void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, cons // foward sweeping for (size_t k = 0; k < number_of_split_cell_lists_; k++) { + // get the corresponding 2D/3D split cell index (m, n) + // e.g., for k = 0, split_cell_index = (0,0), for k = 3, split_cell_index = (1,0), etc. const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k); + // get the number of cells belonging to the split cell k + // i_max = (M - m - 1) / 3 + 1, j_max = (N - n - 1) / 3 + 1 + // e.g. all_cells = (M,N) = (6, 9), (m, n) = (1, 1), then i_max = 2, j_max = 3 const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get1DMeshSize(all_cells_k); + const size_t number_of_cells = get1DMeshSize(all_cells_k); // i_max * j_max + // looping over all cells in the split cell k for (size_t l = 0; l < number_of_cells; l++) { + // get the 2D/3D cell index of the l-th cell in the split cell k + // (i , j) = (m + 3 * (l / j_max), n + 3 * l % i_max) + // e.g. all_cells = (M,N) = (6, 9), (m, n) = (1, 1), l = 0, then (i, j) = (1, 1) + // l = 1, then (i, j) = (1, 4), l = 3, then (i, j) = (4, 1), etc. const Arrayi cell_index = split_cell_index + 3 * transfer1DtoMeshIndex(all_cells_k, l); + // get the list of particles in the cell (i, j) const ConcurrentIndexVector &cell_list = getCellDataList(cell_index_lists_, cell_index); + // looping over all particles in the cell (i, j) for (const size_t index_i : cell_list) { local_dynamics_function(index_i); From 28b9e9aed058bb5a044db996c0b94c8c9d63ada5 Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Wed, 18 Sep 2024 10:44:09 +0200 Subject: [PATCH 17/18] fix after merge --- src/shared/meshes/cell_linked_list.cpp | 8 ++++---- src/shared/meshes/cell_linked_list.hpp | 17 +++++++++++++---- src/shared/meshes/mesh_iterators.h | 1 - .../test_particle_for_split.cpp | 6 +++--- .../test_particle_for_split.cpp | 6 +++--- 5 files changed, 23 insertions(+), 15 deletions(-) diff --git a/src/shared/meshes/cell_linked_list.cpp b/src/shared/meshes/cell_linked_list.cpp index 34214dd954..6d10e43950 100644 --- a/src/shared/meshes/cell_linked_list.cpp +++ b/src/shared/meshes/cell_linked_list.cpp @@ -10,9 +10,8 @@ namespace SPH { //=================================================================================================// BaseCellLinkedList:: - BaseCellLinkedList(SPHAdaptation &sph_adaptation) - : BaseMeshField("CellLinkedList"), - kernel_(*sph_adaptation.getKernel()) {} + BaseCellLinkedList(BaseParticles &base_particles, SPHAdaptation &sph_adaptation) + : BaseMeshField("CellLinkedList"), kernel_(*sph_adaptation.getKernel()) {} //=================================================================================================// CellLinkedList::CellLinkedList(BoundingBox tentative_bounds, Real grid_spacing, BaseParticles &base_particles, SPHAdaptation &sph_adaptation) @@ -21,7 +20,8 @@ CellLinkedList::CellLinkedList(BoundingBox tentative_bounds, Real grid_spacing, index_list_size_(SMAX(base_particles.ParticlesBound(), cell_offset_list_size_)), dv_particle_index_(base_particles.registerDiscreteVariableOnly("ParticleIndex", index_list_size_)), dv_cell_offset_(base_particles.registerDiscreteVariableOnly("CellOffset", cell_offset_list_size_)), - cell_index_lists_(nullptr), cell_data_lists_(nullptr), number_of_split_cell_lists_(pow(3, Dimensions)) + cell_index_lists_(nullptr), cell_data_lists_(nullptr), + number_of_split_cell_lists_(static_cast(pow(3, Dimensions))) { allocateMeshDataMatrix(); single_cell_linked_list_level_.push_back(this); diff --git a/src/shared/meshes/cell_linked_list.hpp b/src/shared/meshes/cell_linked_list.hpp index dc89dd19b9..65ad93682c 100644 --- a/src/shared/meshes/cell_linked_list.hpp +++ b/src/shared/meshes/cell_linked_list.hpp @@ -83,6 +83,15 @@ void CellLinkedList::searchNeighborsByParticles( }); } //=================================================================================================// +inline size_t get_1d_number_of_cells(const Array2i &all_cells) +{ + return all_cells[0] * all_cells[1]; +} +inline size_t get_1d_number_of_cells(const Array3i &all_cells) +{ + return all_cells[0] * all_cells[1] * all_cells[2]; +} +//=================================================================================================// template void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, const LocalDynamicsFunction &local_dynamics_function) { @@ -96,7 +105,7 @@ void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, cons // i_max = (M - m - 1) / 3 + 1, j_max = (N - n - 1) / 3 + 1 // e.g. all_cells = (M,N) = (6, 9), (m, n) = (1, 1), then i_max = 2, j_max = 3 const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get1DMeshSize(all_cells_k); // i_max * j_max + const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); // i_max * j_max // looping over all cells in the split cell k for (size_t l = 0; l < number_of_cells; l++) @@ -121,7 +130,7 @@ void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, cons { const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k - 1); const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get1DMeshSize(all_cells_k); + const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); for (size_t l = 0; l < number_of_cells; l++) { @@ -143,7 +152,7 @@ void CellLinkedList::particle_for_split(const execution::ParallelPolicy &, const { const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k); const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get1DMeshSize(all_cells_k); + const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); parallel_for( IndexRange(0, number_of_cells), @@ -167,7 +176,7 @@ void CellLinkedList::particle_for_split(const execution::ParallelPolicy &, const { const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k - 1); const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get1DMeshSize(all_cells_k); + const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); parallel_for( IndexRange(0, number_of_cells), diff --git a/src/shared/meshes/mesh_iterators.h b/src/shared/meshes/mesh_iterators.h index 5f3c29e869..d9adb1a3bd 100644 --- a/src/shared/meshes/mesh_iterators.h +++ b/src/shared/meshes/mesh_iterators.h @@ -33,7 +33,6 @@ #define MESH_ITERATORS_H #include "base_data_package.h" -#include "base_mesh.h" namespace SPH { diff --git a/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp b/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp index d1fe326df9..e0eb37daa0 100644 --- a/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp +++ b/tests/unit_tests_src/for_2D_build/meshes/test_2d_particle_for_split/test_particle_for_split.cpp @@ -20,9 +20,9 @@ TEST(test_meshes, split_for) body.defineMaterial(); body.generateParticles(); auto &particles = body.getBaseParticles(); - const auto &pos = *particles.registerSharedVariable("Position"); - auto &quantity = *particles.registerSharedVariable("Quantity", [&](size_t i) - { return pos[i].norm(); }); + const auto pos = particles.registerStateVariable("Position"); + auto quantity = particles.registerStateVariable("Quantity", [&](size_t i) -> Real + { return pos[i].norm(); }); InnerRelation inner(body); diff --git a/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp b/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp index fe608b7699..ef4c9264e4 100644 --- a/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp +++ b/tests/unit_tests_src/for_3D_build/meshes/test_3d_particle_for_split/test_particle_for_split.cpp @@ -18,9 +18,9 @@ TEST(test_meshes, split_for) body.defineMaterial(); body.generateParticles(); auto &particles = body.getBaseParticles(); - const auto &pos = *particles.registerSharedVariable("Position"); - auto &quantity = *particles.registerSharedVariable("Quantity", [&](size_t i) - { return pos[i].norm(); }); + const auto pos = particles.registerStateVariable("Position"); + auto quantity = particles.registerStateVariable("Quantity", [&](size_t i) -> Real + { return pos[i].norm(); }); InnerRelation inner(body); From 2615d5378ee7a78ecdefb0d643f1ebbb2a284ecf Mon Sep 17 00:00:00 2001 From: WeiyiVirtonomy Date: Wed, 18 Sep 2024 11:23:38 +0200 Subject: [PATCH 18/18] use prod --- src/shared/meshes/cell_linked_list.hpp | 17 ++++------------- 1 file changed, 4 insertions(+), 13 deletions(-) diff --git a/src/shared/meshes/cell_linked_list.hpp b/src/shared/meshes/cell_linked_list.hpp index 65ad93682c..527695667e 100644 --- a/src/shared/meshes/cell_linked_list.hpp +++ b/src/shared/meshes/cell_linked_list.hpp @@ -83,15 +83,6 @@ void CellLinkedList::searchNeighborsByParticles( }); } //=================================================================================================// -inline size_t get_1d_number_of_cells(const Array2i &all_cells) -{ - return all_cells[0] * all_cells[1]; -} -inline size_t get_1d_number_of_cells(const Array3i &all_cells) -{ - return all_cells[0] * all_cells[1] * all_cells[2]; -} -//=================================================================================================// template void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, const LocalDynamicsFunction &local_dynamics_function) { @@ -105,7 +96,7 @@ void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, cons // i_max = (M - m - 1) / 3 + 1, j_max = (N - n - 1) / 3 + 1 // e.g. all_cells = (M,N) = (6, 9), (m, n) = (1, 1), then i_max = 2, j_max = 3 const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); // i_max * j_max + const size_t number_of_cells = all_cells_k.prod(); // i_max * j_max // looping over all cells in the split cell k for (size_t l = 0; l < number_of_cells; l++) @@ -130,7 +121,7 @@ void CellLinkedList::particle_for_split(const execution::SequencedPolicy &, cons { const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k - 1); const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); + const size_t number_of_cells = all_cells_k.prod(); for (size_t l = 0; l < number_of_cells; l++) { @@ -152,7 +143,7 @@ void CellLinkedList::particle_for_split(const execution::ParallelPolicy &, const { const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k); const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); + const size_t number_of_cells = all_cells_k.prod(); parallel_for( IndexRange(0, number_of_cells), @@ -176,7 +167,7 @@ void CellLinkedList::particle_for_split(const execution::ParallelPolicy &, const { const Arrayi split_cell_index = transfer1DtoMeshIndex(3 * Arrayi::Ones(), k - 1); const Arrayi all_cells_k = (all_cells_ - split_cell_index - Arrayi::Ones()) / 3 + Arrayi::Ones(); - const size_t number_of_cells = get_1d_number_of_cells(all_cells_k); + const size_t number_of_cells = all_cells_k.prod(); parallel_for( IndexRange(0, number_of_cells),