Алгоритм нахождения собственных векторов по заданным собственным значениям матрицы 3х3 в C # - PullRequest
0 голосов
/ 25 ноября 2018

Я пытаюсь найти лучшие хитбоксы OOBB для моих мешей, используя PCA.Для того, чтобы сделать это, мне нужны собственные векторы, но я отчасти потерян, как вычислять их, не используя огромную библиотеку.

Я реализовал алгоритм, который вычисляет три собственных значения с учетом матрицы 3x3.Код для этого изначально взят из Википедии :

private Vector3 CalculateEigenvalues(ref Matrix3 A)
{
    Vector3 val = new Vector3(0, 0, 0);

    float p1 = A.M12 * A.M12 + A.M13 * A.M13 + A.M23 * A.M23;
    if (p1 == 0)
    {
        val.X = A.M11;
        val.Y = A.M22;
        val.Z = A.M33;
    }
    else
    {
        float q = A.Trace / 3f;
        float p2 = (float)(Math.Pow(A.M11 - q, 2) + Math.Pow(A.M22 - q, 2) + Math.Pow(A.M33 - q, 2)) + 2 * p1;
        float p = (float)Math.Sqrt(p2 / 6);
        Matrix4 I = Matrix4.Identity;
        Matrix4.Mult(ref I, q, out Matrix4 tmp);
        Matrix4 tmp2 = Matrix4.Subtract(new Matrix4(A), tmp);
        Matrix4 B = Matrix4.Mult(tmp2, 1 / p);
        float r = new Matrix3(B).Determinant / 2;

        float phi = 0;
        if (r <= -1)
            phi = (float)Math.PI / 3;
        else if (r >= 1)
            phi = 0;
        else
            phi = (float)Math.Acos(r) / 3;
        val.X = q + 2 * p * (float)Math.Cos(phi);
        val.Z = q + 2 * p * (float)Math.Cos(phi + (2 * Math.PI / 3));
        val.Y = 3 * q - val.X - val.Z;
    }
    return val;
}

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

Поэтому мой вопрос:

Если у меня есть матрица 3x3 и три собственных значения, есть ли какой-нибудь простой способвычислить соответствующие собственные векторы без использования внешних библиотек?

1 Ответ

0 голосов
/ 25 ноября 2018

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

// Observe that the function doesn't use rZ,
// it is expected that it will become zero vector in triangular form
static Vector3 EigenVector(Vector3 rX, Vector3 rY, Vector3 rZ, float lambda)
{
    // Move RHS to LHS
    rX.X -= lambda;
    rY.Y -= lambda;
    // Transform to upper triangle
    rY -= rX * (rY.X / rX.X);
    // Backsubstitute
    var res = new Vector3(1f);
    res.Y = -rY.Z / rY.Y;
    res.X = -(rX.Y * res.Y + rX.Z * res.Z) / rX.X;
    return res;
}

// Case of eigenvalue with algebraic multiplicity two
static (Vector3, Vector3) EigenVector2(Vector3 rX, Vector3 rY, Vector3 rZ, float lambda)
{
    // Move RHS to LHS
    rX.X -= lambda;
    float x2 = rX.Y / rX.X;
    float x3 = rX.Z / rX.X;
    return (new Vector3(x2, 1, 0), new Vector3(x3, 0, 1));
}

static void Main(string[] args)
{
    var rX = new Vector3(1, -3, 3);
    var rY = new Vector3(3, -5, 3);
    var rZ = new Vector3(6, -6, 4);
    var e = EigenVector(rX, rY, rZ, 4);
    var e2 = EigenVector2(rX, rY, rZ, 2);

    System.Diagnostics.Debug.WriteLine(e.ToString());
    System.Diagnostics.Debug.WriteLine(e2.Item1.ToString());
    System.Diagnostics.Debug.WriteLine(e2.Item2.ToString());
}

<0,5 0,5 1>
<3 1 0>
<-3 0 1>

В реальной жизни нужно много проверять на ошибки.Входные данные взяты из этого документа .

...