Physics: PGS resolution implementation issues

I’m making an attempt to implement a PGS resolution for the collision response, however I’m at present caught.

The way in which wherein the issue is manifested is, for instance, a field that rests on the bottom, if the field has a small rotation, then it may get uncontrolled after some time.

As well as, the containers in a stack will not be steady, however they separate in a short time. (With out pgs, my money batteries are decently steady)

Loads grateful for any recommendation! 🙂

int maxIterations = 10;
for (int i = 0; i  0)
           proceed;
       // Calculate J
       float e = std::min(objA.restitution, objB.restitution);
       float J = -(1.0f + e) * glm::dot(relativeVelocity, contact.regular);
       J /= objA.invMass + objB.invMass +
           glm::dot(glm::cross(rA, contact.regular), objA.inverseInertia * glm::cross(rA, contact.regular)) +
           glm::dot(glm::cross(rB, contact.regular), objB.inverseInertia * glm::cross(rB, contact.regular));
       // Clamp
       float temp = contactPoint.accumulatedImpulse;
       contactPoint.accumulatedImpulse = glm::max(temp + J, 0.0f);
       float deltaImpulse = contactPoint.accumulatedImpulse - temp;
       // Add regular
       glm::vec3 deltaNormalImpulse = deltaImpulse * contact.regular;
       // Apply impulse
       if (glm::size(deltaNormalImpulse) > 1e-6f) {
           objA.linearVelocity -= deltaNormalImpulse * objA.invMass;
           objA.angularVelocity -= objA.inverseInertia * glm::cross(rA, deltaNormalImpulse);
           objB.linearVelocity += deltaNormalImpulse * objB.invMass;
           objB.angularVelocity += objB.inverseInertia * glm::cross(rB, deltaNormalImpulse);
           converged = false;
       }
   }
   if (converged)
       break;
}

// place correction
if (objA.isStatic)
   objB.place += (collisionNormal * depth);
else if (objB.isStatic)
   objA.place += (-collisionNormal * depth);
else
{
   float totalInverseMass = objA.invMass + objB.invMass;
   float displacementA = (objA.invMass / totalInverseMass) * depth * energy;
   float displacementB = (objB.invMass / totalInverseMass) * depth * energy;

   objA.place -= collisionNormal * displacementA;
   objB.place += collisionNormal * displacementB;
}