diff --git a/Simulations/MassSpringSystemSimulator.cpp b/Simulations/MassSpringSystemSimulator.cpp index 0ea12ac..6ac432f 100644 --- a/Simulations/MassSpringSystemSimulator.cpp +++ b/Simulations/MassSpringSystemSimulator.cpp @@ -10,12 +10,12 @@ MassSpringSystemSimulator::MassSpringSystemSimulator() auto first = addMassPoint(Vec3(0, 0, 0), Vec3(-1, 0, 0), true); auto second = addMassPoint(Vec3(0, 2, 0), Vec3(1, 0, 0), true); addSpring(first, second, 1.0); - m_fStiffness = 40; + } const char* MassSpringSystemSimulator::getTestCasesStr() { - return "Demo1,Demo2,Demo3,Demo4"; + return "Euler,LeapFrog,Midpoint"; } void MassSpringSystemSimulator::initUI(DrawingUtilitiesClass* DUC) @@ -71,21 +71,19 @@ void MassSpringSystemSimulator::notifyCaseChanged(int testCase) switch (m_iTestCase) { case 0: - cout << "Demo1 !\n"; + cout << "Euler !\n"; //simulateTimestep(1); break; case 1: - cout << "Demo2 \n"; + break; + case 2: + cout << "Midpoint \n"; //m_iNumSpheres = 100; //m_fSphereSize = 0.05f; break; - case 2: - cout << "Demo3 !\n"; - break; - case 3: - cout << "Demo 4 !\n"; + default: //cout << "Demo4 !\n"; break; @@ -124,19 +122,20 @@ void MassSpringSystemSimulator::simulateTimestep(float timeStep) case 0: //update the masspoint cout << "Euler \n"; - - cout << "Midpoint\n"; + Euler(0, 1, 0, timeStep); + break; - case 1:cout << "demo 2 \n"; + case 1: break; - case 2:cout << "demo 3\n"; - break; - case 3: cout << "demo 4\n"; + + case 2:cout << "midpoint \n"; + Midpoint(0, 1, timeStep); break; + default: break; } - Euler(0, 1, 0, timeStep); + //Euler(0, 1, 0, timeStep); } void MassSpringSystemSimulator::onClick(int x, int y) @@ -216,15 +215,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; } @@ -234,6 +244,51 @@ Vec3 MassSpringSystemSimulator::calculateAcceleration(Vec3 force, float mass) return force / mass; } +void MassSpringSystemSimulator::Midpoint(int index1, int index2, float timestep) { + //here some implementation about Midpoint + auto massPoint1 = masspoints.at(index1); + auto massPoint2 = masspoints.at(index2); + //old position + auto mp = massPoint1->position; + auto mp2 = massPoint2 ->position; + //old Velocity + auto mOld_v = massPoint1->velocity; + auto m2Old_v = massPoint1->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(int index1, int index2, int indexSpring, float timestep) { //take old position and send to calculatePositionTimestepEuler @@ -247,11 +302,11 @@ void MassSpringSystemSimulator::Euler(int index1, int index2, int indexSpring, f // Force of spring is -k * (l - L) * normalizedVector [for P2 we can take -F1) auto force = -m_fStiffness * (lengthVector - springs.at(0).initialLength) * normalized; auto foreP2 = -1 * force; - auto veloc = calcualteVelocityTimestepEuler(mp->velocity, timestep, calculateAcceleration(force, 10.)); - auto pos = calcualtePositionTimestepEuler(mp->position, timestep, veloc); + auto veloc = calculateVelocityTimestepEuler(mp->velocity, timestep, calculateAcceleration(force, 10.)); + auto pos = calculatePositionTimestepEuler(mp->position, timestep, veloc); - auto veloc2 = calcualteVelocityTimestepEuler(mp2->velocity, timestep, calculateAcceleration(foreP2, 10.)); - auto pos2 = calcualtePositionTimestepEuler(mp2->position, timestep, veloc2); + auto veloc2 = calculateVelocityTimestepEuler(mp2->velocity, timestep, calculateAcceleration(foreP2, 10.)); + auto pos2 = calculatePositionTimestepEuler(mp2->position, timestep, veloc2); // Update Positions and Velocity mp->position = pos; diff --git a/Simulations/MassSpringSystemSimulator.h b/Simulations/MassSpringSystemSimulator.h index d40aee9..ecdc8d3 100644 --- a/Simulations/MassSpringSystemSimulator.h +++ b/Simulations/MassSpringSystemSimulator.h @@ -39,8 +39,11 @@ 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(int index1, int index2, 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(int index1, int index2, int indexSpring, float timestep);