2010-10-14 4 views
1

В настоящее время я нахожусь в процессе написания функции, чтобы найти «точную» ограничительную сферу для набора точек в трехмерном пространстве. Я думаю, что у меня есть приличное понимание процесса до сих пор, но я застрял.Поиск собственных векторов ковариационной матрицы для создания трехмерной ограничивающей сферы

Вот что я работаю с: A) точек в 3D-пространстве B) 3x3 ковариационной матрицы хранятся в классе 4x4 матрицы (на который ссылается ячеек m0, m1, m2, m3, m4, ЭСТ, вместо строк и cols)

Я нашел 3 собственных значения для ковариационной матрицы точек, и я создал функцию преобразования матрицы в уменьшенную форму эшелона строки (rref) с помощью исключения Гаусса.

Я тестировал обе эти функции против цифр в примерах, которые я нашел в Интернете, и они, похоже, работают правильно.

Следующим шагом является нахождение собственных векторов, используя уравнение: (M - & лямбда; * I) * V

..., где М представляет собой матрицу ковариации, & Lambda; - одно из собственных значений, I - единичная матрица, а V - собственный вектор.

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

Вот функция до сих пор:

 
Vect eigenVector(const Matrix &M, const float eval) { 
    Matrix A = Matrix(M); 
    A -= Matrix(IDENTITY)*eval; 
    A.rref(); 
    return Vect(A[m3],A[m7],A[m11]); 
} 

3х3 ковариационная матрица передается как М, и собственное значение, как Eval. Матрица (IDENTITY) возвращает идентификационную матрицу. m3, m7 и m11 соответствуют крайнему правому столбцу матрицы 4x3.

Вот матрица пример 3x3 (хранится в классе 4x4 матрицы) Я использую для тестирования функций:

 
Matrix(1.5f, 0.5f, 0.75f, 0, 
     0.5f, 0.5f, 0.25f, 0, 
     0.75f, 0.25f, 0.5f, 0, 
      0,  0, 0, 0); 

Я правильно получать собственные значения 2,097, 0,3055, 0.09756 из (?) моя другая функция.

собственного вектор() выше правильно вычитает прошел собственное значение от диагонали (0,0 1,1 2,2)

матрицы А после того, как RREF():

 
[(1, 0, 0, -0), 
(-0, 1, 0, -0), 
(-0, -0, 1, -0), 
(0, 0, 0, -2.09694)] 

Для RREF() функция, я использую переведенную функцию python, которую можно найти здесь: http://elonen.iki.fi/code/misc-notes/python-gaussj/index.html

Что должна передать матрица, передаваемая rref(), чтобы получить собственный вектор?

Благодаря

ответ

1

(M - & Lambda; I), V не является уравнением, это просто выражение. Однако (M - & lambda; I) V = 0. И это уравнение связывает собственные векторы с собственными значениями.

Так что, предполагая, что функция rref работает, я бы предположил, что вы создаете дополненную матрицу как [(M - λI) | 0], где 0 обозначает нулевой вектор. Это похоже на то, что вы делаете уже, поэтому я должен был предположить, что ваша функция rref сломана. Или, альтернативно, он не знает, как обрабатывать матрицы 4x4 (в отличие от матриц 4x3, что и ожидалось бы для расширенной матрицы).

+0

Так я и думал, и я потратил довольно много времени на различные варианты. На самом деле, я думаю, что <0,0,0> является одним из бесконечных решений, которые решают эту систему линейных уравнений. Сначала мне нужно было выбросить систему произвольным значением, прежде чем я смог получить значимые результаты. См. Мой ответ ниже для того, что я сделал. Спасибо, что нашли время, чтобы протянуть руку, хотя :) –

1

А, еще несколько часов изнурительных исследований, мне удалось решить мою проблему.

Проблема в том, что нет «одного» набора собственных векторов, а есть бесконечное число с различными величинами.

Метод, который я выбрал, состоял в том, чтобы вместо RREF использовать форму REF (строка эшелона), оставляя достаточную информацию в матрице, чтобы я мог заменить произвольное значение для z и работать назад для решения для y и x. Затем я нормализовал вектор, чтобы получить единичный собственный вектор, который должен работать для моих целей.

Мой окончательный код:

 
Vect eigenVector(const Matrix &M, const float eVal) { 
    Matrix A = Matrix(M); 
    A -= Matrix(IDENTITY)*eVal; 
    A.ref(); 
    float K = 16; // Arbitrary value 
    float J = -K*A[m6]; // Substitute in K to find J 
    float I = -K*A[m2]-J*A[m1]; // Substitute in K and J to find I 

    Vect eVec = Vect(I,J,K); 
    eVec.norm(); // Normalize eigenvector 

    return eVec; 
} 

Единственная странность в том, что собственные векторы выходят лицом в противоположном направлении, чем я ожидал (они были сведены на нет!), Но это спорный вопрос.