Движок физики 2D с упругостью полигонов с нарушениями

Я пытался создать простой двухмерный физический движок на C#. Следовал уроку от “Two-Bit Coding” на YouTube, но заметил, что квадраты с двумя контактными точками вообще не отскакивают. Причину этой проблемы я выяснил: импульсы были просто слишком маленькими при разрешении столкновений.

Я добавил немного упругости, умножив импульс на значение. Однако это не решение, так как если столкновение происходит очень близко к центру тела (в середине края квадрата), импульс будет огромным, что приведет к нарушению законов сохранения энергии.

Я попробовал другую формулу для разрешения столкновения, на этот раз из Википедии. Формула предназначена для 3D, но я предположил, что она будет работать и для 2D.
Формула из Википедии

Вот реализация на C#:

        public void Resolve_collision_wiki(in Physics_manifold contact)
        {
            Physics_body body_a = contact.body_a;
            Physics_body body_b = contact.body_b;
            Physics_vector normal = contact.normal;

            int contact_count = contact.contact_count;

            float e = 1;

            contact_list[0] = contact.contact_point_1;
            contact_list[1] = contact.contact_point_2;

            for (int i = 0; i < contact_count; i++)
            {
                Physics_vector ra = contact_list[i] - body_a.Position;
                Physics_vector rb = contact_list[i] - body_b.Position;

                center_a_contact_list[i] = ra;
                center_b_contact_list[i] = rb;

                Physics_vector relative_velocity = body_b.Linear_velocity - body_a.Linear_velocity + Physics_math.Cross_product(body_b.Angular_velocity, rb) - Physics_math.Cross_product(body_a.Angular_velocity, ra);

                //() сначала, затем скалярное произведение?
                float numerator = Physics_math.Dot_product(-(1 + e) * relative_velocity, normal);

                if (numerator < 0)
                    return;

                //https://www.cs.cmu.edu/~baraff/sigcourse/notesd2.pdf#page=16

                numerator /= body_a.inverse_mass + body_b.inverse_mass + Physics_math.Dot_product(Physics_math.Cross_product(body_a.inverse_inertia * Physics_math.Cross_product(ra, normal), ra) + Physics_math.Cross_product(body_b.inverse_inertia * Physics_math.Cross_product(rb, normal), rb), normal);
                numerator = MathF.Max(numerator, 0);

                j_list[i] = numerator;
            }

            //Применение в первом цикле сделает второй импульс немного неправильным
            for (int i = 0; i < contact_count; i++)
            {
                body_a.Linear_velocity -= j_list[i] * body_a.inverse_mass * normal;
                body_a.Angular_velocity -=  Physics_math.Cross_product(center_a_contact_list[i], normal * j_list[i]) * body_a.inverse_inertia;

                body_b.Linear_velocity += j_list[i] * body_b.inverse_mass * normal;
                body_b.Angular_velocity +=  Physics_math.Cross_product(center_b_contact_list[i], normal * j_list[i]) * body_b.inverse_inertia;
            }
        }

Реализация скалярного произведения и двумерного векторного произведения:

        public static float Dot_product(Physics_vector a, Physics_vector b)
        {
            return a.x * b.x + a.y * b.y;
        }

        public static float Cross_product(Physics_vector a, Physics_vector b)
        {
            return a.x * b.y - a.y * b.x;
        }

        public static Physics_vector Cross_product(float a, Physics_vector v)
        {
            return new Physics_vector(-a * v.y, a * v.x);
        }

Но высота отскока все равно драматически изменяется, когда я перемещаю контактные точки ближе друг к другу, и становится ниже, когда контактные точки находятся дальше друг от друга. Более того, значение e равное 1 не обеспечивает идеально упругое столкновение для квадратов, как это происходит для кругов.