From d6630c97134c3d8b5fd471041d9924f774c0994d Mon Sep 17 00:00:00 2001 From: Winterless Date: Tue, 29 Aug 2023 21:36:32 +0200 Subject: [PATCH 1/2] Scale the factor of ShellMidSurfaceBounding to particle spacing; Change to a better function of random number. --- .../particle_generator_lattice_supplementary.cpp | 7 ++++++- .../particle_generator_lattice_supplementary.cpp | 7 ++++++- .../relax_dynamics/relax_dynamics.cpp | 12 ++++-------- .../relax_dynamics/relax_dynamics.h | 8 +++----- .../test_2d_shell_particle_relaxation.cpp | 3 +-- .../test_3d_shell_particle_relaxation.cpp | 2 +- 6 files changed, 21 insertions(+), 18 deletions(-) diff --git a/src/for_2D_build/particle_generator/particle_generator_lattice_supplementary.cpp b/src/for_2D_build/particle_generator/particle_generator_lattice_supplementary.cpp index 86ac057f28..3d6f07af26 100644 --- a/src/for_2D_build/particle_generator/particle_generator_lattice_supplementary.cpp +++ b/src/for_2D_build/particle_generator/particle_generator_lattice_supplementary.cpp @@ -54,6 +54,11 @@ void ThickSurfaceParticleGeneratorLattice::initializeGeometricVariables() Real interval = (Real)planned_number_of_particles_ / (Real)all_cells_; if (interval <= 0) interval = 1; // It has to be lager than 0. + + // initialize a uniform distribution between 0 (inclusive) and 1 (exclusive) + std::mt19937_64 rng; + std::uniform_real_distribution unif(0, 1); + // Add a particle in each interval, randomly. We will skip the last intervals if we already reach the number of particles for (int i = 0; i < number_of_lattices[0]; ++i) for (int j = 0; j < number_of_lattices[1]; ++j) @@ -63,7 +68,7 @@ void ThickSurfaceParticleGeneratorLattice::initializeGeometricVariables() { if (body_shape_.checkContain(particle_position)) { - Real random_real = (Real)rand() / (RAND_MAX); + Real random_real = unif(rng); // If the random_real is smaller than the interval, add a particle, only if we haven't reached the max. number of particles if (random_real <= interval && base_particles_.total_real_particles_ < planned_number_of_particles_) { diff --git a/src/for_3D_build/particle_generator/particle_generator_lattice_supplementary.cpp b/src/for_3D_build/particle_generator/particle_generator_lattice_supplementary.cpp index 58b2013190..269a8f5fa4 100644 --- a/src/for_3D_build/particle_generator/particle_generator_lattice_supplementary.cpp +++ b/src/for_3D_build/particle_generator/particle_generator_lattice_supplementary.cpp @@ -52,10 +52,15 @@ void ThickSurfaceParticleGeneratorLattice::initializeGeometricVariables() Real number_of_particles = total_volume_ / avg_particle_volume_ + 0.5; planned_number_of_particles_ = int(number_of_particles); + // initialize a uniform distribution between 0 (inclusive) and 1 (exclusive) + std::mt19937_64 rng; + std::uniform_real_distribution unif(0, 1); + // Calculate the interval based on the number of particles. Real interval = planned_number_of_particles_ / (all_cells_ + TinyReal); if (interval <= 0) interval = 1; // It has to be lager than 0. + // Add a particle in each interval, randomly. We will skip the last intervals if we already reach the number of particles. for (int i = 0; i < number_of_lattices[0]; ++i) for (int j = 0; j < number_of_lattices[1]; ++j) @@ -66,7 +71,7 @@ void ThickSurfaceParticleGeneratorLattice::initializeGeometricVariables() { if (body_shape_.checkContain(particle_position)) { - Real random_real = (Real)rand() / (RAND_MAX); + Real random_real = unif(rng); // If the random_real is smaller than the interval, add a particle, only if we haven't reached the max. number of particles. if (random_real <= interval && base_particles_.total_real_particles_ < planned_number_of_particles_) { diff --git a/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.cpp b/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.cpp index a67f9b1568..bae559dedf 100644 --- a/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.cpp +++ b/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.cpp @@ -159,19 +159,17 @@ void RelaxationStepComplex::exec(Real dt) } //=================================================================================================// ShellMidSurfaceBounding:: - ShellMidSurfaceBounding(NearShapeSurface &body_part, BaseInnerRelation &inner_relation, - Real thickness, Real level_set_refinement_ratio) + ShellMidSurfaceBounding(NearShapeSurface &body_part, BaseInnerRelation &inner_relation) : BaseLocalDynamics(body_part), RelaxDataDelegateInner(inner_relation), pos_(particles_->pos_), constrained_distance_(0.5 * sph_body_.sph_adaptation_->MinimumSpacing()), particle_spacing_ref_(sph_body_.sph_adaptation_->MinimumSpacing()), - thickness_(thickness), level_set_refinement_ratio_(level_set_refinement_ratio), level_set_shape_(DynamicCast(this, sph_body_.body_shape_)) {} //=================================================================================================// void ShellMidSurfaceBounding::update(size_t index_i, Real dt) { Vecd none_normalized_normal = level_set_shape_->findLevelSetGradient(pos_[index_i]); Vecd normal = level_set_shape_->findNormalDirection(pos_[index_i]); - Real factor = none_normalized_normal.squaredNorm() / level_set_refinement_ratio_; + Real factor = 0.2 * none_normalized_normal.norm(); pos_[index_i] -= factor * constrained_distance_ * normal; } //=================================================================================================// @@ -292,12 +290,10 @@ void ShellNormalDirectionPrediction::SmoothingNormal::update(size_t index_i, Rea } //=================================================================================================// ShellRelaxationStepInner:: - ShellRelaxationStepInner(BaseInnerRelation &inner_relation, Real thickness, - Real level_set_refinement_ratio, bool level_set_correction) + ShellRelaxationStepInner(BaseInnerRelation &inner_relation, bool level_set_correction) : RelaxationStepInner(inner_relation, level_set_correction), update_shell_particle_position_(*real_body_), - mid_surface_bounding_(near_shape_surface_, inner_relation, - thickness, level_set_refinement_ratio) {} + mid_surface_bounding_(near_shape_surface_, inner_relation) {} //=================================================================================================// void ShellRelaxationStepInner::exec(Real ite_p) { diff --git a/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.h b/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.h index ecbac45399..48be2ec033 100644 --- a/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.h +++ b/src/shared/particle_dynamics/relax_dynamics/relax_dynamics.h @@ -302,15 +302,14 @@ class ShellMidSurfaceBounding : public BaseLocalDynamics, public RelaxDataDelegateInner { public: - ShellMidSurfaceBounding(NearShapeSurface &body_part, BaseInnerRelation &inner_relation, - Real thickness, Real level_set_refinement_ratio); + ShellMidSurfaceBounding(NearShapeSurface &body_part, BaseInnerRelation &inner_relation); virtual ~ShellMidSurfaceBounding(){}; void update(size_t index_i, Real dt = 0.0); protected: StdLargeVec &pos_; Real constrained_distance_; - Real particle_spacing_ref_, thickness_, level_set_refinement_ratio_; + Real particle_spacing_ref_; LevelSetShape *level_set_shape_; }; @@ -440,8 +439,7 @@ class ShellNormalDirectionPrediction : public BaseDynamics class ShellRelaxationStepInner : public RelaxationStepInner { public: - explicit ShellRelaxationStepInner(BaseInnerRelation &inner_relation, Real thickness, - Real level_set_refinement_ratio, bool level_set_correction = false); + explicit ShellRelaxationStepInner(BaseInnerRelation &inner_relation, bool level_set_correction = false); virtual ~ShellRelaxationStepInner(){}; SimpleDynamics update_shell_particle_position_; diff --git a/tests/2d_examples/test_2d_shell_particle_relaxation/test_2d_shell_particle_relaxation.cpp b/tests/2d_examples/test_2d_shell_particle_relaxation/test_2d_shell_particle_relaxation.cpp index 1343658e6b..0a9690f2a5 100644 --- a/tests/2d_examples/test_2d_shell_particle_relaxation/test_2d_shell_particle_relaxation.cpp +++ b/tests/2d_examples/test_2d_shell_particle_relaxation/test_2d_shell_particle_relaxation.cpp @@ -62,8 +62,7 @@ int main() /** Random reset the particle position. */ SimpleDynamics random_pipe_body_particles(pipe_body); /** A Physics relaxation step. */ - relax_dynamics::ShellRelaxationStepInner - relaxation_step_pipe_body_inner(pipe_body_inner, thickness, level_set_refinement_ratio); + relax_dynamics::ShellRelaxationStepInner relaxation_step_pipe_body_inner(pipe_body_inner); relax_dynamics::ShellNormalDirectionPrediction shell_normal_prediction(pipe_body_inner, thickness); pipe_body.addBodyStateForRecording("UpdatedIndicator"); /** diff --git a/tests/3d_examples/test_3d_shell_particle_relaxation/test_3d_shell_particle_relaxation.cpp b/tests/3d_examples/test_3d_shell_particle_relaxation/test_3d_shell_particle_relaxation.cpp index cfb172e625..6764c55986 100644 --- a/tests/3d_examples/test_3d_shell_particle_relaxation/test_3d_shell_particle_relaxation.cpp +++ b/tests/3d_examples/test_3d_shell_particle_relaxation/test_3d_shell_particle_relaxation.cpp @@ -70,7 +70,7 @@ int main(int ac, char *av[]) //---------------------------------------------------------------------- SimpleDynamics random_imported_model_particles(imported_model); /** A Physics relaxation step. */ - relax_dynamics::ShellRelaxationStepInner relaxation_step_inner(imported_model_inner, thickness, level_set_refinement_ratio); + relax_dynamics::ShellRelaxationStepInner relaxation_step_inner(imported_model_inner); relax_dynamics::ShellNormalDirectionPrediction shell_normal_prediction(imported_model_inner, thickness); //---------------------------------------------------------------------- // Particle relaxation starts here. From d495f8ee8b2651646ed9a935eda070d637b46e8a Mon Sep 17 00:00:00 2001 From: Winterless Date: Sun, 3 Sep 2023 23:15:30 +0200 Subject: [PATCH 2/2] check cases --- .../test_2d_ball_shell_collision/ball_shell_collision.cpp | 3 +-- .../src/shell_beam_collision.cpp | 7 +++---- .../test_3d_roof_parametric_cvt.cpp | 4 ++-- 3 files changed, 6 insertions(+), 8 deletions(-) diff --git a/tests/2d_examples/test_2d_ball_shell_collision/ball_shell_collision.cpp b/tests/2d_examples/test_2d_ball_shell_collision/ball_shell_collision.cpp index 4399691b14..765623903e 100644 --- a/tests/2d_examples/test_2d_ball_shell_collision/ball_shell_collision.cpp +++ b/tests/2d_examples/test_2d_ball_shell_collision/ball_shell_collision.cpp @@ -112,8 +112,7 @@ int main(int ac, char *av[]) // Define the methods for particle relaxation for wall boundary. //---------------------------------------------------------------------- SimpleDynamics wall_boundary_random_particles(wall_boundary); - relax_dynamics::ShellRelaxationStepInner - relaxation_step_wall_boundary_inner(wall_boundary_inner, thickness, level_set_refinement_ratio); + relax_dynamics::ShellRelaxationStepInner relaxation_step_wall_boundary_inner(wall_boundary_inner); relax_dynamics::ShellNormalDirectionPrediction shell_normal_prediction(wall_boundary_inner, thickness, cos(Pi / 3.75)); wall_boundary.addBodyStateForRecording("UpdatedIndicator"); //---------------------------------------------------------------------- diff --git a/tests/2d_examples/test_2d_shell_beam_collision/src/shell_beam_collision.cpp b/tests/2d_examples/test_2d_shell_beam_collision/src/shell_beam_collision.cpp index d54bdfb510..bb6653f574 100644 --- a/tests/2d_examples/test_2d_shell_beam_collision/src/shell_beam_collision.cpp +++ b/tests/2d_examples/test_2d_shell_beam_collision/src/shell_beam_collision.cpp @@ -83,9 +83,9 @@ int main(int ac, char *av[]) //---------------------------------------------------------------------- SPHSystem sph_system(system_domain_bounds, resolution_ref); /** Tag for running particle relaxation for the initially body-fitted distribution */ - sph_system.setRunParticleRelaxation(true); + sph_system.setRunParticleRelaxation(false); /** Tag for starting with relaxed body-fitted particles distribution */ - sph_system.setReloadParticles(false); + sph_system.setReloadParticles(true); sph_system.handleCommandlineOptions(ac, av); IOEnvironment io_environment(sph_system); //---------------------------------------------------------------------- @@ -135,8 +135,7 @@ int main(int ac, char *av[]) // Define the methods for particle relaxation for wall boundary. //---------------------------------------------------------------------- SimpleDynamics shell_random_particles(shell); - relax_dynamics::ShellRelaxationStepInner - relaxation_step_shell_inner(shell_inner, thickness, level_set_refinement_ratio); + relax_dynamics::ShellRelaxationStepInner relaxation_step_shell_inner(shell_inner); relax_dynamics::ShellNormalDirectionPrediction shell_normal_prediction(shell_inner, thickness, cos(Pi / 3.75)); shell.addBodyStateForRecording("UpdatedIndicator"); //---------------------------------------------------------------------- diff --git a/tests/3d_examples/test_3d_roof_parametric_cvt/test_3d_roof_parametric_cvt.cpp b/tests/3d_examples/test_3d_roof_parametric_cvt/test_3d_roof_parametric_cvt.cpp index 9f4d502d18..a9ea64c9b4 100644 --- a/tests/3d_examples/test_3d_roof_parametric_cvt/test_3d_roof_parametric_cvt.cpp +++ b/tests/3d_examples/test_3d_roof_parametric_cvt/test_3d_roof_parametric_cvt.cpp @@ -14,13 +14,13 @@ using namespace SPH; Real to_rad(Real angle) { return angle * Pi / 180; } -void relax_shell(RealBody &plate_body, Real thickness, Real level_set_refinement_ratio) +void relax_shell(RealBody &plate_body, Real thickness) { // BUG: apparently only works if dp > thickness, otherwise ShellNormalDirectionPrediction::correctNormalDirection() throws error InnerRelation imported_model_inner(plate_body); SimpleDynamics random_imported_model_particles(plate_body); - relax_dynamics::ShellRelaxationStepInner relaxation_step_inner(imported_model_inner, thickness, level_set_refinement_ratio); + relax_dynamics::ShellRelaxationStepInner relaxation_step_inner(imported_model_inner); relax_dynamics::ShellNormalDirectionPrediction shell_normal_prediction(imported_model_inner, thickness); //---------------------------------------------------------------------- // Particle relaxation starts here.