Вычислительный тензор инерции в 2D

Я исследую, как найти инерцию для 2D фигуры. Контур этой формы зацеплен с несколькими точками, координаты x и y каждой точки уже известны.

Я знаю выражение Ixx, Iyy а также Ixy но тело не имеет массы. Как мне поступить?

2 ответа

Решение

Для любой фигуры вам нужно разбить ее на треугольники и обрабатывать каждый треугольник отдельно. Затем в итоге объединили результаты, используя следующие правила

В общем и целом

% Combined total area of all triangles
total_area = SUM( area(i), i=1:n )
total_mass = SUM( mass(i), i=1:n )
% Combined centroid (center of mass) coordinates
combined_centroid_x = SUM( mass(i)*centroid_x(i), i=1:n)/total_mass
combined_centroid_y = SUM( mass(i)*centroid_y(i), i=1:n)/total_mass
% Each distance to triangle (squared)
centroid_distance_sq(i) = centroid_x(i)*centroid_x(i)+centroid_y(i)*centroid_y(i)
% Combined mass moment of inertia
combined_mmoi = SUM(mmoi(i)+mass(i)*centroid_distance_sq(i), i=1:n)

Теперь для каждого треугольника.

Рассмотрим три угловые вершины с векторными координатами, точки A, B и C

a=[ax,ay]
b=[bx,by]
c=[cx,cy]

и следующие комбинации точек и перекрестных произведений (скалярных)

a·a = ax*ax+ay*ay
b·b = bx*bx+by*by
c·c = cx*cx+cy*cy
a·b = ax*bx+ay*by
b·c = bx*cx+by*cy
c·a = cx*ax+cy*ay
a×b = ax*by-ay*bx
b×c = bx*cy-by*cx
c×a = cx*ay-cy*ax

Свойства треугольника (с t(i) толщина и rho массовая плотность)

area(i) = 1/2*ABS( a×b + b×c + c×a )
mass(i) = rho*t(i)*area(i)
centroid_x(i) = 1/3*(ax + bx + cx)
centroid_y(i) = 1/3*(ay + by + cy)
mmoi(i) = 1/6*mass(i)*( a·a + b·b + c·c + a·b + b·c + c·a )

По компоненту выше

area(i) = 1/2*ABS( ax*(by-cy)+ay*(cx-bx)+bx*cy-by*cx)
mmoi(i) = mass(i)/6*(ax^2+ax*(bx+cx)+bx^2+bx*cx+cx^2+ay^2+ay*(by+cy)+by^2+by*cy+cy^2)

аппендикс

Немного теории здесь. Площадь каждого треугольника определяется с помощью

Area = 1/2 * || (b-a) × (c-b) ||

где × является векторным перекрестным произведением, и || .. || является векторной нормой (функцией длины).

Треугольник параметризован двумя переменными t а также s такой, что двойной интеграл A = INT(INT(1,dx),dy) дает общую площадь

% position r(s,t) = [x,y]
[x,y] = [ax,ay] + t*[bx-ax, by-zy] + t*s*[cx-bx,cy-by]

% gradient directions along s and t
(dr/dt) = [bx-ax,by-ay] + s*[cx-bx,cy-by]
(dr/ds) = t*[cx-bx,cy-by]

% Integration area element
dA = || (dr/ds)×(dr/dt) || = (2*A*t)*ds*dt
%
%   where A = 1/2*||(b-a)×(c-b)||

% Check that the integral returns the area
Area = INT( INT( 2*A*t,s=0..1), t=0..1) = 2*A*(1/2) = A

% Mass moment of inertia components

         /  /  /  | y^2+z^2  -x*y    -x*z   |
I = 2*m*|  |  | t*|  -x*y   x^2+z^2  -y*z   | dz ds dt
        /  /  /   |  -x*z    -y*z   x^2+y^2 |

% where [x,y] are defined from the parametrization

Я хочу немного исправить отличный ответ Джона Алексиу:

  1. Алгоритм треугольника mmoi предполагает, что углы будут определены относительно центра масс треугольника (центроида). Так что вычтите центроид из углов перед вычислением mmoi
  2. Алгоритм Shape mmoi предполагает, что центроиды каждого отдельного треугольника будут определены относительно центра масс формы (комбинированный центроид). Поэтому перед вычислением mmoi вычтите комбинированный центр тяжести из центра тяжести каждого треугольника.

Таким образом, код результата будет выглядеть так:

      public static float CalculateMMOI(Triangle[] triangles, float thickness, float density)
{
    float[] mass = new float[triangles.Length];
    float[] mmoi = new float[triangles.Length];
    Vector2[] centroid = new Vector2[triangles.Length];
    
    float combinedMass = 0;
    float combinedMMOI = 0;
    Vector2 combinedCentroid = new Vector2(0, 0);
    
    for (var i = 0; i < triangles.Length; ++i) 
    {
        mass[i] = triangles[i].CalculateMass(thickness, density);
        mmoi[i] = triangles[i].CalculateMMOI(mass[i]);
        centroid[i] = triangles[i].CalculateCentroid();

        combinedMass += mass[i];
        combinedCentroid += mass[i] * centroid[i];
    }
    combinedCentroid /= combinedMass;


    for (var i = 0; i < triangles.Length; ++i) {
        combinedMMOI += mmoi[i] + mass[i] * (centroid[i] - combinedCentroid).LengthSquared();
    }

    return combinedMMOI;
}


public struct Triangle
{
    public Vector2 A, B, C;
    
    public float CalculateMass(float thickness, float density)
    {
        var area = CalculateArea();
        return area * thickness * density;
    }
    
    public float CalculateArea()
    {
        return 0.5f * Math.Abs(Vector2.Cross(A - B, A - C));
    }

    public float CalculateMMOI(float mass)
    {
        var centroid = CalculateCentroid()
        var a = A - centroid;
        var b = B - centroid;
        var c = C - centroid;
        
        var aa = Vector2.Dot(a, a);
        var bb = Vector2.Dot(b, b);
        var cc = Vector2.Dot(c, c);
        var ab = Vector2.Dot(a, b);
        var bc = Vector2.Dot(b, c);
        var ca = Vector2.Dot(c, a);
        return (aa + bb + cc + ab + bc + ca) * mass / 6f;
    }

    public Vector2 CalculateCentroid()
    {
        return (A + B + C) / 3f;
    }
}


public struct Vector2
{
    public float X, Y;
    
    public float LengthSquared()
    {
        return X * X + Y * Y;
    }
    
    public static float Dot(Vector2 a, Vector2 b)
    {
        return a.X * b.X + a.Y * b.Y;
    }

    public static float Cross(Vector2 a, Vector2 b)
    {
        return a.X * b.Y - a.Y * b.X;
    }
}
Другие вопросы по тегам