diff --git a/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.cpp b/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.cpp index 0450ced2dc..0d831c9dfa 100644 --- a/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.cpp +++ b/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.cpp @@ -374,43 +374,51 @@ namespace SPH pos_n_(particles_->pos_n_), pos_0_(particles_->pos_0_), vel_n_(particles_->vel_n_), - dvel_dt_prior_(particles_->dvel_dt_prior_), - mass_(particles_->mass_) + dvel_dt_prior_(particles_->dvel_dt_prior_) { - stiffness_ = stiffness; - damping_ratio_ = damping_ratio; + // calculate total mass + total_mass_ = 0.0; + for (int i = 0; i < particles_->mass_.size(); i++) + { + total_mass_ += particles_->mass_[i]; + } + // scale stiffness and damping by mass here, so it's not necessary in each iteration + stiffness_ = stiffness / total_mass_; + damping_coeff_ = stiffness * damping_ratio / total_mass_; } //=================================================================================================// + SpringDamperConstraintParticleWise::~SpringDamperConstraintParticleWise() + {} + //=================================================================================================// void SpringDamperConstraintParticleWise::setupDynamics(Real dt) { particles_->total_ghost_particles_ = 0; } //=================================================================================================// - Vecd SpringDamperConstraintParticleWise::getAcceleration(Vecd& disp, Real mass) - { + Vecd SpringDamperConstraintParticleWise::getSpringForce(size_t index_i, Vecd& disp) + { Vecd spring_force(0); for(int i = 0; i < disp.size(); i++) { - spring_force[i] = -stiffness_[i] * disp[i] / mass; - } + spring_force[i] = -stiffness_[i] * disp[i]; + } return spring_force; } //=================================================================================================// - Vecd SpringDamperConstraintParticleWise::getDampingForce(size_t index_i, Real mass) + Vecd SpringDamperConstraintParticleWise::getDampingForce(size_t index_i) { Vecd damping_force(0); for(int i = 0; i < vel_n_[index_i].size(); i++) { - damping_force[i] = -stiffness_[i] * damping_ratio_ * vel_n_[index_i][i] / mass; + damping_force[i] = -damping_coeff_[i] * vel_n_[index_i][i]; } return damping_force; } //=================================================================================================// void SpringDamperConstraintParticleWise::Update(size_t index_i, Real dt) { - Vecd disp_from_0 = pos_n_[index_i] - pos_0_[index_i]; - dvel_dt_prior_[index_i] += getAcceleration(disp_from_0, mass_[index_i]); - dvel_dt_prior_[index_i] += getDampingForce(index_i, mass_[index_i]); + dvel_dt_prior_[index_i] += getSpringForce(index_i, pos_n_[index_i] - pos_0_[index_i]); + dvel_dt_prior_[index_i] += getDampingForce(index_i); } //=================================================================================================// AccelerationForBodyPartInBoundingBox:: diff --git a/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.h b/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.h index 6eb6cbc999..c3fe257cc5 100644 --- a/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.h +++ b/SPHINXsys/src/shared/particle_dynamics/solid_dynamics/solid_dynamics.h @@ -296,23 +296,23 @@ namespace SPH * @brief Exerts spring force and damping force in the form of acceleration to each particle. * The spring force is calculated based on the difference from the particle's initial position. * The damping force is calculated based on the particle's current velocity. + * Only for 3D applications */ class SpringDamperConstraintParticleWise : public ParticleDynamicsSimple, public SolidDataSimple { public: - SpringDamperConstraintParticleWise(SolidBody* body, Vecd stiffness, Real damping_ratio = 0.01); - virtual ~SpringDamperConstraintParticleWise() {}; + SpringDamperConstraintParticleWise(SolidBody* body, Vecd stiffness, Real damping_ratio = 0.05); + ~SpringDamperConstraintParticleWise(); protected: - StdLargeVec& mass_; + Real total_mass_; StdLargeVec& pos_n_,& pos_0_,& vel_n_,& dvel_dt_prior_; Vecd stiffness_; - // damping component parallel to the spring force component - // damping coefficient = stiffness_ * damping_ratio_ - Real damping_ratio_; + Vecd damping_coeff_; // damping component parallel to the spring force component + virtual void setupDynamics(Real dt = 0.0) override; - virtual Vecd getAcceleration(Vecd& disp, Real mass); - virtual Vecd getDampingForce(size_t index_i, Real mass); + virtual Vecd getSpringForce(size_t index_i, Vecd& disp); + virtual Vecd getDampingForce(size_t index_i); virtual void Update(size_t index_i, Real dt = 0.0) override; }; /** @@ -389,12 +389,6 @@ namespace SPH Real ReduceFunction(size_t index_i, Real dt = 0.0) override; }; - /** - * @function getSmallestTimeStepAmongSolidBodies - * @brief computing smallest time step to use in a simulation - */ - Real getSmallestTimeStepAmongSolidBodies(SPHBodyVector solid_bodies); - /** * @class DeformationGradientTensorBySummation * @brief computing deformation gradient tensor by summation