diff --git a/Simulations/MassSpringSystemSimulator.cpp b/Simulations/MassSpringSystemSimulator.cpp index 4ce7d5f..c9e8026 100644 --- a/Simulations/MassSpringSystemSimulator.cpp +++ b/Simulations/MassSpringSystemSimulator.cpp @@ -101,7 +101,7 @@ void MassSpringSystemSimulator::notifyCaseChanged(int testCase) //calculate Midpoint for one step and print results setIntegrator(MIDPOINT); cout << "\n\n\t -- MIDPOINT RESULT --\n"; - simulateTimestep(0.005); + simulateTimestep(1); printSpring(springs.at(0)); break; } @@ -188,16 +188,13 @@ void MassSpringSystemSimulator::simulateTimestep(float timeStep) springs.erase(springs.begin() + i); continue; } - auto mp1 = sp.mp1.lock(); - auto mp2 = sp.mp2.lock(); - if (m_iIntegrator == EULER) { - Euler(*mp1.get(), *mp2.get(), sp, timeStep); + Euler(sp, timeStep); } else if (m_iIntegrator == MIDPOINT) { - //TODO: Add Midpoint + Midpoint(sp, timeStep); } else if (m_iIntegrator == LEAPFROG) { //TODO: Add Leapfrog @@ -283,15 +280,26 @@ Vec3 MassSpringSystemSimulator::getVelocityOfMassPoint(int index) void MassSpringSystemSimulator::applyExternalForce(Vec3 force) { - + //eine Vorstellung: for all Masspoints, update force } -Vec3 MassSpringSystemSimulator::calcualtePositionTimestepEuler(Vec3 oldPosition, float timestep, Vec3 velocity) +int MassSpringSystemSimulator::LengthCalculator(Vec3 position1, Vec3 position2) { + + Vec3 PosVector = position1 - position2; + //wurzel aus Vektor + int length = sqrt(pow(PosVector.x, 2) + pow(PosVector.y, 2) + pow(PosVector.z, 2)); + + return length; +} + + + +Vec3 MassSpringSystemSimulator::calculatePositionTimestepEuler(Vec3 oldPosition, float timestep, Vec3 velocity) { return oldPosition + timestep * velocity; } -Vec3 MassSpringSystemSimulator::calcualteVelocityTimestepEuler(Vec3 oldVelocity, float timestep, Vec3 acceleration) +Vec3 MassSpringSystemSimulator::calculateVelocityTimestepEuler(Vec3 oldVelocity, float timestep, Vec3 acceleration) { return oldVelocity + acceleration * timestep; } @@ -301,10 +309,59 @@ Vec3 MassSpringSystemSimulator::calculateAcceleration(Vec3 force, float mass) return force / mass; } -void MassSpringSystemSimulator::Euler(MassPoint& masspoint1, MassPoint& masspoint2, Spring& spring, float timestep) +void MassSpringSystemSimulator::Midpoint(Spring& spring, float timestep) { + //here some implementation about Midpoint + auto massPoint1 = spring.mp1.lock(); + auto massPoint2 = spring.mp2.lock(); + + //old position + auto mp = massPoint1->position; + auto mp2 = massPoint2->position; + //old Velocity + auto mOld_v = massPoint1->velocity; + auto m2Old_v = massPoint2->velocity; + + Vec3 PosVector = mp - mp2; + //Abstand ausrechnen + int d = LengthCalculator(mp, mp2); + //normalize + Vec3 PosNorm1 = PosVector / d; + Vec3 PosNorm2 = -1 * PosNorm1; + + Vec3 Force = -m_fStiffness * (d - springs.at(0).initialLength) * PosNorm1; + Vec3 Force2 = -1 * Force; + + Vec3 oldAcc = calculateAcceleration(Force,m_fMass); + Vec3 oldAcc2 = calculateAcceleration(Force2, m_fMass); + + //Midpoint calculator + //Pos of Midstep + Vec3 PosOfMidstep = mp + 0.5 * timestep * mOld_v; + Vec3 PosOfMidstep2 = mp2 + 0.5 * timestep * m2Old_v; + //Vel at Midstep + Vec3 VelAtMidstep = mOld_v + 0.5 * timestep * oldAcc; + Vec3 VelAtMidstep2 = m2Old_v + 0.5 * timestep * oldAcc2; + + Vec3 NewPos = mp + timestep * VelAtMidstep; + Vec3 NewPos2 = mp2 + timestep * VelAtMidstep2; + + Vec3 NewVel = mOld_v + timestep * oldAcc; + Vec3 NewVel2 = m2Old_v + timestep * oldAcc2; + //cout << NewPos; + //cout << NewVel; + massPoint1->position = NewPos; + massPoint1->velocity = NewVel; + massPoint2->position = NewPos2; + massPoint2->velocity = NewVel2; +} + +void MassSpringSystemSimulator::Euler(Spring& spring, float timestep) { + auto massPoint1 = spring.mp1.lock(); + auto massPoint2 = spring.mp2.lock(); + //take old position and send to calculatePositionTimestepEuler - auto PosVector = masspoint1.position - masspoint2.position; + auto PosVector = massPoint1->position - massPoint2->position; auto lengthVector = sqrt(PosVector.x * PosVector.x + PosVector.y * PosVector.y + PosVector.z * PosVector.z); auto normalized = PosVector / lengthVector; @@ -312,18 +369,18 @@ void MassSpringSystemSimulator::Euler(MassPoint& masspoint1, MassPoint& masspoin // Force of spring is -k * (l - L) * normalizedVector [for P2 we can take -F1) auto force = -m_fStiffness * (lengthVector - spring.initialLength) * normalized; auto foreP2 = -1 * force; - auto veloc = calcualteVelocityTimestepEuler(masspoint1.velocity, timestep, calculateAcceleration(force, 10.)); - auto pos = calcualtePositionTimestepEuler(masspoint1.position, timestep, veloc); + auto veloc = calculatePositionTimestepEuler(massPoint1->velocity, timestep, calculateAcceleration(force, 10.)); + auto pos = calculatePositionTimestepEuler(massPoint1->position, timestep, veloc); - auto veloc2 = calcualteVelocityTimestepEuler(masspoint2.velocity, timestep, calculateAcceleration(foreP2, 10.)); - auto pos2 = calcualtePositionTimestepEuler(masspoint2.position, timestep, veloc2); + auto veloc2 = calculateVelocityTimestepEuler(massPoint2->velocity, timestep, calculateAcceleration(foreP2, 10.)); + auto pos2 = calculatePositionTimestepEuler(massPoint2->position, timestep, veloc2); // Update Positions and Velocity - masspoint1.position = pos; - masspoint1.velocity = veloc; + massPoint1->position = pos; + massPoint1->velocity = veloc; - masspoint2.position = pos2; - masspoint2.velocity = veloc2; + massPoint2->position = pos2; + massPoint2->velocity = veloc2; } void MassSpringSystemSimulator::printSpring(const Spring& spring) diff --git a/Simulations/MassSpringSystemSimulator.h b/Simulations/MassSpringSystemSimulator.h index 520ed02..c2f4878 100644 --- a/Simulations/MassSpringSystemSimulator.h +++ b/Simulations/MassSpringSystemSimulator.h @@ -39,10 +39,13 @@ public: Vec3 getVelocityOfMassPoint(int index); void applyExternalForce(Vec3 force); - Vec3 calcualtePositionTimestepEuler(Vec3 oldPosition, float timestep, Vec3 veloctiy); - Vec3 calcualteVelocityTimestepEuler(Vec3 oldVelocity, float timestep, Vec3 acceleration); + void Midpoint(Spring& spring, float timestep); + int LengthCalculator(Vec3 position1, Vec3 position2); + + Vec3 calculatePositionTimestepEuler(Vec3 oldPosition, float timestep, Vec3 veloctiy); + Vec3 calculateVelocityTimestepEuler(Vec3 oldVelocity, float timestep, Vec3 acceleration); Vec3 calculateAcceleration(Vec3 acceleration, float mass); - void Euler(MassPoint& masspoint1, MassPoint& masspoint2, Spring& spring, float timestep); + void Euler(Spring& spring, float timestep); void printSpring(const Spring& spring);