diff --git a/src/shared/geometries/level_set.cpp b/src/shared/geometries/level_set.cpp index 2748d97af5..bc056d9e4b 100644 --- a/src/shared/geometries/level_set.cpp +++ b/src/shared/geometries/level_set.cpp @@ -29,13 +29,13 @@ MultilevelLevelSet::MultilevelLevelSet( finish_data_packages.exec(); probe_signed_distance_set_.push_back( probe_signed_distance_vector_keeper_ - .template createPtr>(*mesh_data_set_[0])); + .template createPtr(*mesh_data_set_[0])); probe_normal_direction_set_.push_back( probe_normal_direction_vector_keeper_ - .template createPtr>(*mesh_data_set_[0])); + .template createPtr(*mesh_data_set_[0])); probe_level_set_gradient_set_.push_back( probe_level_set_gradient_vector_keeper_ - .template createPtr>(*mesh_data_set_[0])); + .template createPtr(*mesh_data_set_[0])); for (size_t level = 1; level != total_levels_; ++level) { @@ -58,13 +58,13 @@ MultilevelLevelSet::MultilevelLevelSet( probe_signed_distance_set_.push_back( probe_signed_distance_vector_keeper_ - .template createPtr>(*mesh_data_set_[level])); + .template createPtr(*mesh_data_set_[level])); probe_normal_direction_set_.push_back( probe_normal_direction_vector_keeper_ - .template createPtr>(*mesh_data_set_[level])); + .template createPtr(*mesh_data_set_[level])); probe_level_set_gradient_set_.push_back( probe_level_set_gradient_vector_keeper_ - .template createPtr>(*mesh_data_set_[level])); + .template createPtr(*mesh_data_set_[level])); } clean_interface = new CleanInterface{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()}; @@ -96,13 +96,13 @@ MultilevelLevelSet::MultilevelLevelSet( probe_signed_distance_set_.push_back( probe_signed_distance_vector_keeper_ - .template createPtr>(*mesh_data_set_[0])); + .template createPtr(*mesh_data_set_[0])); probe_normal_direction_set_.push_back( probe_normal_direction_vector_keeper_ - .template createPtr>(*mesh_data_set_[0])); + .template createPtr(*mesh_data_set_[0])); probe_level_set_gradient_set_.push_back( probe_level_set_gradient_vector_keeper_ - .template createPtr>(*mesh_data_set_[0])); + .template createPtr(*mesh_data_set_[0])); clean_interface = new CleanInterface{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()}; correct_topology = new CorrectTopology{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()}; @@ -132,24 +132,24 @@ void MultilevelLevelSet::correctTopology(Real small_shift_factor) //=============================================================================================// Real MultilevelLevelSet::probeSignedDistance(const Vecd &position) { - return probe_signed_distance_set_[getProbeLevel(position)]->exec(position); + return probe_signed_distance_set_[getProbeLevel(position)]->update(position); } //=============================================================================================// Vecd MultilevelLevelSet::probeNormalDirection(const Vecd &position) { - return probe_normal_direction_set_[getProbeLevel(position)]->exec(position); + return probe_normal_direction_set_[getProbeLevel(position)]->update(position); } //=============================================================================================// Vecd MultilevelLevelSet::probeLevelSetGradient(const Vecd &position) { - return probe_level_set_gradient_set_[getProbeLevel(position)]->exec(position); + return probe_level_set_gradient_set_[getProbeLevel(position)]->update(position); } //=============================================================================================// size_t MultilevelLevelSet::getProbeLevel(const Vecd &position) { for (size_t level = total_levels_; level != 0; --level){ - MeshCalculateDynamics is_within_core_package{*mesh_data_set_[level - 1]}; - if(is_within_core_package.exec(position)) + IsWithinCorePackage is_within_core_package{*mesh_data_set_[level - 1]}; + if(is_within_core_package.update(position)) return level - 1; // jump out of the loop! } return 0; @@ -158,16 +158,16 @@ size_t MultilevelLevelSet::getProbeLevel(const Vecd &position) Real MultilevelLevelSet::probeKernelIntegral(const Vecd &position, Real h_ratio) { if(mesh_data_set_.size() == 1){ - MeshCalculateDynamics refine_probe{*mesh_data_set_[0]}; - return refine_probe.exec(position); + ProbeKernelIntegral refine_probe{*mesh_data_set_[0]}; + return refine_probe.update(position); } size_t coarse_level = getCoarseLevel(h_ratio); Real alpha = (global_h_ratio_vec_[coarse_level + 1] - h_ratio) / (global_h_ratio_vec_[coarse_level + 1] - global_h_ratio_vec_[coarse_level]); - MeshCalculateDynamics coarse_probe{*mesh_data_set_[coarse_level]}; - MeshCalculateDynamics fine_probe{*mesh_data_set_[coarse_level + 1]}; - Real coarse_level_value = coarse_probe.exec(position); - Real fine_level_value = fine_probe.exec(position); + ProbeKernelIntegral coarse_probe{*mesh_data_set_[coarse_level]}; + ProbeKernelIntegral fine_probe{*mesh_data_set_[coarse_level + 1]}; + Real coarse_level_value = coarse_probe.update(position); + Real fine_level_value = fine_probe.update(position); return alpha * coarse_level_value + (1.0 - alpha) * fine_level_value; } @@ -175,16 +175,16 @@ Real MultilevelLevelSet::probeKernelIntegral(const Vecd &position, Real h_ratio) Vecd MultilevelLevelSet::probeKernelGradientIntegral(const Vecd &position, Real h_ratio) { if(mesh_data_set_.size() == 1){ - MeshCalculateDynamics refine_probe{*mesh_data_set_[0]}; - return refine_probe.exec(position); + ProbeKernelGradientIntegral refine_probe{*mesh_data_set_[0]}; + return refine_probe.update(position); } size_t coarse_level = getCoarseLevel(h_ratio); Real alpha = (global_h_ratio_vec_[coarse_level + 1] - h_ratio) / (global_h_ratio_vec_[coarse_level + 1] - global_h_ratio_vec_[coarse_level]); - MeshCalculateDynamics coarse_probe{*mesh_data_set_[coarse_level]}; - MeshCalculateDynamics fine_probe{*mesh_data_set_[coarse_level + 1]}; - Vecd coarse_level_value = coarse_probe.exec(position); - Vecd fine_level_value = fine_probe.exec(position); + ProbeKernelGradientIntegral coarse_probe{*mesh_data_set_[coarse_level]}; + ProbeKernelGradientIntegral fine_probe{*mesh_data_set_[coarse_level + 1]}; + Vecd coarse_level_value = coarse_probe.update(position); + Vecd fine_level_value = fine_probe.update(position); return alpha * coarse_level_value + (1.0 - alpha) * fine_level_value; } @@ -194,8 +194,8 @@ bool MultilevelLevelSet::probeIsWithinMeshBound(const Vecd &position) bool is_bounded = true; for (size_t l = 0; l != total_levels_; ++l) { - MeshCalculateDynamics probe_is_within_mesh_bound{*mesh_data_set_[l]}; - if (!probe_is_within_mesh_bound.exec(position)) + ProbeIsWithinMeshBound probe_is_within_mesh_bound{*mesh_data_set_[l]}; + if (!probe_is_within_mesh_bound.update(position)) { is_bounded = false; break; diff --git a/src/shared/geometries/level_set.h b/src/shared/geometries/level_set.h index eab2c61a1f..d596c295bd 100644 --- a/src/shared/geometries/level_set.h +++ b/src/shared/geometries/level_set.h @@ -65,8 +65,8 @@ class MultilevelLevelSet : public BaseMeshField { for (size_t l = 0; l != total_levels_; ++l) { - MeshCalculateDynamics write_mesh_field_to_plt(*mesh_data_set_[l]); - write_mesh_field_to_plt.exec(output_file); + WriteMeshFieldToPlt write_mesh_field_to_plt(*mesh_data_set_[l]); + write_mesh_field_to_plt.update(output_file); } } @@ -79,13 +79,13 @@ class MultilevelLevelSet : public BaseMeshField size_t total_levels_; /**< level 0 is the coarsest */ StdVec global_h_ratio_vec_; StdVec mesh_data_set_; - StdVec *> probe_signed_distance_set_; - StdVec *> probe_normal_direction_set_; - StdVec *> probe_level_set_gradient_set_; + StdVec probe_signed_distance_set_; + StdVec probe_normal_direction_set_; + StdVec probe_level_set_gradient_set_; UniquePtrsKeeper mesh_data_ptr_vector_keeper_; - UniquePtrsKeeper> probe_signed_distance_vector_keeper_; - UniquePtrsKeeper> probe_normal_direction_vector_keeper_; - UniquePtrsKeeper> probe_level_set_gradient_vector_keeper_; + UniquePtrsKeeper probe_signed_distance_vector_keeper_; + UniquePtrsKeeper probe_normal_direction_vector_keeper_; + UniquePtrsKeeper probe_level_set_gradient_vector_keeper_; CleanInterface *clean_interface; CorrectTopology *correct_topology; diff --git a/src/shared/meshes/all_mesh_dynamics.cpp b/src/shared/meshes/all_mesh_dynamics.cpp deleted file mode 100644 index 36fa08b76a..0000000000 --- a/src/shared/meshes/all_mesh_dynamics.cpp +++ /dev/null @@ -1,69 +0,0 @@ -/* ------------------------------------------------------------------------- * - * SPHinXsys * - * ------------------------------------------------------------------------- * - * SPHinXsys (pronunciation: s'finksis) is an acronym from Smoothed Particle * - * Hydrodynamics for industrial compleX systems. It provides C++ APIs for * - * physical accurate simulation and aims to model coupled industrial dynamic * - * systems including fluid, solid, multi-body dynamics and beyond with SPH * - * (smoothed particle hydrodynamics), a meshless computational method using * - * particle discretization. * - * * - * SPHinXsys is partially funded by German Research Foundation * - * (Deutsche Forschungsgemeinschaft) DFG HU1527/6-1, HU1527/10-1, * - * HU1527/12-1 and HU1527/12-4. * - * * - * Portions copyright (c) 2017-2023 Technical University of Munich and * - * the authors' affiliations. * - * * - * Licensed under the Apache License, Version 2.0 (the "License"); you may * - * not use this file except in compliance with the License. You may obtain a * - * copy of the License at http://www.apache.org/licenses/LICENSE-2.0. * - * * - * ------------------------------------------------------------------------- */ -/** - * @file base_relax_dynamics.h - * @brief This is the classes of particle relaxation in order to produce body fitted - * initial particle distribution. - * @author Chi Zhang and Xiangyu Hu - */ - -#include "all_mesh_dynamics.h" - -#include "all_body_relations.h" -#include "all_particle_dynamics.h" -#include "base_kernel.h" -#include "base_particles.hpp" -#include "cell_linked_list.h" - -namespace SPH -{ -//=================================================================================================// -// void DefineShapeOnMeshWithGridDataPackages::exec() -// { -// initialize_data_in_a_cell.exec(); -// tag_a_cell_is_inner_package.exec(); -// initialize_index_mesh.exec(); -// initialize_cell_neighborhood.exec(); -// mesh_data_.resizeMeshVariableData(); -// } -//=================================================================================================// -// void CleanInterface::exec() -// { -// mark_near_interface.exec(); -// redistance_interface.exec(); -// reinitialize_level_set.exec(); -// update_level_set_gradient.exec(); -// update_kernel_integrals.exec(); -// } -//=================================================================================================// -// void CorrectTopology::exec() -// { -// mark_near_interface.exec(); -// for (size_t i = 0; i != 10; ++i) -// diffuse_level_set_sign.exec(); -// update_level_set_gradient.exec(); -// update_kernel_integrals.exec(); -// } -//=================================================================================================// -} // namespace SPH -//=================================================================================================// diff --git a/src/shared/meshes/all_mesh_dynamics.h b/src/shared/meshes/all_mesh_dynamics.h index 304abbe933..b94092c545 100644 --- a/src/shared/meshes/all_mesh_dynamics.h +++ b/src/shared/meshes/all_mesh_dynamics.h @@ -57,8 +57,8 @@ class FinishDataPackages : public BaseMeshDynamics mesh_data_.resizeMeshVariableData(); Real far_field_distance = grid_spacing_ * (Real)buffer_width_; - initialize_data_for_singular_package.exec(0, -far_field_distance); - initialize_data_for_singular_package.exec(1, far_field_distance); + initialize_data_for_singular_package.update(0, -far_field_distance); + initialize_data_for_singular_package.update(1, far_field_distance); initialize_basic_data_for_a_package.exec(); update_level_set_gradient.exec(); @@ -72,7 +72,7 @@ class FinishDataPackages : public BaseMeshDynamics Real grid_spacing_; size_t buffer_width_; - MeshSingleDynamics initialize_data_for_singular_package{mesh_data_}; + InitializeDataForSingularPackage initialize_data_for_singular_package{mesh_data_}; MeshAllDynamics tag_a_cell_is_inner_package{mesh_data_}; MeshInnerDynamics initialize_index_mesh{mesh_data_}; MeshInnerDynamics initialize_cell_neighborhood{mesh_data_}; @@ -90,7 +90,7 @@ class ProbeNormalDirection : public BaseMeshLocalDynamics Vecd update(const Vecd &position) { - Vecd probed_value = probe_level_set_gradient.exec(position); + Vecd probed_value = probe_level_set_gradient.update(position); Real threshold = 1.0e-2 * data_spacing_; while (probed_value.norm() < threshold) @@ -98,13 +98,13 @@ class ProbeNormalDirection : public BaseMeshLocalDynamics Vecd jittered = position; // jittering for (int l = 0; l != position.size(); ++l) jittered[l] += rand_uniform(-0.5, 0.5) * 0.5 * data_spacing_; - probed_value = probe_level_set_gradient.exec(jittered); + probed_value = probe_level_set_gradient.update(jittered); } return probed_value.normalized(); } private: - MeshCalculateDynamics probe_level_set_gradient{mesh_data_}; + ProbeLevelSetGradient probe_level_set_gradient{mesh_data_}; }; class CleanInterface : public BaseMeshDynamics diff --git a/src/shared/meshes/mesh_dynamics.h b/src/shared/meshes/mesh_dynamics.h index f32f54261c..596facd55b 100644 --- a/src/shared/meshes/mesh_dynamics.h +++ b/src/shared/meshes/mesh_dynamics.h @@ -52,12 +52,16 @@ class BaseMeshDynamics public: BaseMeshDynamics(MeshWithGridDataPackages<4> &mesh_data) : mesh_data_(mesh_data), - all_cells_(mesh_data.AllCells()){}; + all_cells_(mesh_data.AllCells()), + num_grid_pkgs_(mesh_data.num_grid_pkgs_), + meta_data_cell_(mesh_data.meta_data_cell_){}; virtual ~BaseMeshDynamics(){}; protected: MeshWithGridDataPackages<4> &mesh_data_; Arrayi all_cells_; + size_t &num_grid_pkgs_; + std::pair* &meta_data_cell_; template void grid_parallel_for(const FunctionOnData &function) @@ -68,6 +72,21 @@ class BaseMeshDynamics function(cell_index); }); } + + /** Iterator on a collection of mesh data packages. parallel computing. */ + template + void package_parallel_for(const FunctionOnData &function) + { + parallel_for(IndexRange(2, num_grid_pkgs_), + [&](const IndexRange &r) + { + for (size_t i = r.begin(); i != r.end(); ++i) + { + function(i); + } + }, + ap); + } }; /** @@ -111,7 +130,7 @@ class MeshInnerDynamics : public LocalDynamicsType, public BaseMeshDynamics void exec() { - mesh_data_.package_parallel_for( + package_parallel_for( [&](size_t package_index) { this->update(package_index); @@ -136,10 +155,10 @@ class MeshCoreDynamics : public LocalDynamicsType, public BaseMeshDynamics void exec() { - mesh_data_.package_parallel_for( + package_parallel_for( [&](size_t package_index) { - std::pair &metadata = mesh_data_.meta_data_cell_[package_index]; + std::pair &metadata = meta_data_cell_[package_index]; if (metadata.second == 1) { this->update(package_index); @@ -148,43 +167,5 @@ class MeshCoreDynamics : public LocalDynamicsType, public BaseMeshDynamics ); }; }; - -/** - * @class MeshSingleDynamics - * @brief Mesh dynamics for only core cells on the mesh - */ -template -class MeshSingleDynamics : public LocalDynamicsType, public BaseMeshDynamics -{ - public: - template - MeshSingleDynamics(DynamicsIdentifier &identifier, Args &&...args) - : LocalDynamicsType(identifier, std::forward(args)...), - BaseMeshDynamics(identifier){}; - virtual ~MeshSingleDynamics(){}; - - template - void exec(Args &&...args) - { - this->update(std::forward(args)...); - }; -}; - -template -class MeshCalculateDynamics : public LocalDynamicsType, public BaseMeshDynamics -{ - public: - template - MeshCalculateDynamics(DynamicsIdentifier &identifier, Args &&...args) - : LocalDynamicsType(identifier, std::forward(args)...), - BaseMeshDynamics(identifier){}; - virtual ~MeshCalculateDynamics(){}; - - template - ReturnType exec(Args &&...args) - { - return this->update(std::forward(args)...); - }; -}; } // namespace SPH #endif // MESH_DYNAMICS_H \ No newline at end of file diff --git a/src/shared/meshes/mesh_with_data_packages.h b/src/shared/meshes/mesh_with_data_packages.h index 04ddf25684..00e5ce98cc 100644 --- a/src/shared/meshes/mesh_with_data_packages.h +++ b/src/shared/meshes/mesh_with_data_packages.h @@ -96,8 +96,8 @@ class MeshWithGridDataPackages : public Mesh template explicit MeshWithGridDataPackages(BoundingBox tentative_bounds, Real data_spacing, size_t buffer_size) : Mesh(tentative_bounds, pkg_size * data_spacing, buffer_size), - data_spacing_(data_spacing), - global_mesh_(mesh_lower_bound_ + 0.5 * data_spacing * Vecd::Ones(), data_spacing, all_cells_ * pkg_size) + global_mesh_(mesh_lower_bound_ + 0.5 * data_spacing * Vecd::Ones(), data_spacing, all_cells_ * pkg_size), + data_spacing_(data_spacing) { allocateIndexDataMatrix(); }; @@ -117,12 +117,12 @@ class MeshWithGridDataPackages : public Mesh CellNeighborhood *cell_neighborhood_; /**< 3*3(*3) array to store indicies of neighborhood cells. */ std::pair *meta_data_cell_; /**< metadata for each occupied cell: (arrayi)cell index, (int)core1/inner0. */ BaseMesh global_mesh_; /**< the mesh for the locations of all possible data points. */ + size_t num_grid_pkgs_ = 2; /**< the number of all distinct packages, initially only 2 singular packages. */ protected: MeshVariableAssemble all_mesh_variables_; /**< all mesh variables on this mesh. */ static constexpr int pkg_size = PKG_SIZE; /**< the size of the data package matrix*/ const Real data_spacing_; /**< spacing of data in the data packages*/ - size_t num_grid_pkgs_ = 2; /**< the number of all distinct packages, initially only 2 singular packages. */ using MetaData = std::pair; /**< stores the metadata for each cell: (int)singular0/inner1/core2, (size_t)package data index*/ MeshDataMatrix index_data_mesh_; /**< metadata for all cells. */ using NeighbourIndex = std::pair; /**< stores shifted neighbour info: (size_t)package index, (arrayi)local grid index. */ @@ -218,21 +218,6 @@ class MeshWithGridDataPackages : public Mesh void computeGradient(MeshVariable &in_variable, MeshVariable &out_variable, const size_t package_index); - /** Iterator on a collection of mesh data packages. parallel computing. */ - template - void package_parallel_for(const FunctionOnData &function) - { - parallel_for( - IndexRange(2, num_grid_pkgs_), - [&](const IndexRange &r) - { - for (size_t i = r.begin(); i != r.end(); ++i) - { - function(i); - } - }, - ap); - } void registerOccupied(std::pair &occupied) {