Расчет обратной матрицы с использованием итерационного метода

Мой профессор назначил нам лабораторию, в которой я застрял и не могу найти ничего в Интернете. Это специально для матрицы с собственным значением ‹ 1.

Учитывая матрицу 3x3 A = {(.5, -1, 0), (0, .6666, 0), (.5, -1, .6666)}, вы должны вычислить обратную, используя A^- 1 = I + B + B ^ 2 + B ^ 3 ..... B ^ 20, где I — единичная матрица, а B = I — A. Есть два выхода, один из которых A^-1, а другой будучи A * A^-1, чтобы проверить результат.

До сих пор я пробовал две разные функции для умножения. Извините за мое ужасное имя. i, j и k — счетчик строк, столбцов и циклов соответственно.

void mPower(double x[][3], double y[][3], double z[][3])
{
  int i, j, k;

  for (i = 0; i < 3; i++) 
    for (j = 0; j < 3; j++)   {
      for (k = 0; k < 3; k++)
        z[i][j] = z[i][j] + (x[i][k] * y[k][j]);
}
}

mPower - моя, наверное, неправильная функция для умножения. В качестве альтернативы я пробовал

void multiplyMatrix(double x[][3], double y[][3], double z[][3])
{
  int i, j, k;
   int sum = 0;

for (i = 0; i <= 2; i++) {
  for (j = 0; j <= 2; j++) {
     sum = 0;
     for (k = 0; k <= 2; k++) {
        sum = sum + x[i][k] * y[k][j];
     }
     z[i][j] = sum;
  }
  }

}

Однако я получаю только нули.

void printIt(double a[][3], double b[][3])
{
 int i, j;

 for (i = 0; i < 3; i++)
 {
     printf ("\n\t\t\t");
     for (j = 0; j < 3; j++)
         printf ("%7.4f", a[i][j]);
     printf ("\t");
     for (j = 0; j < 3; j++)
         printf ("%7.4f", b[i][j]);
         printf ("\n\n\t");
 }
}

int main()
{
  int i, j, k;
  double b1[][3] = {.5, 1., 0., 0., .6666, 0., -.5, -1., .6666};
  double b3[3][3] = {0};
  double addOn[3][3];
  double aInverse[3][3];
  double identity[][3] = {1., 0., 0., 0., 1., 0., 0., 0., 1.};
  mPower (b1, b1, b3);
  for (k = 0; k < 20; k++)     {
      addOn[i][j] = addOn[i][j] + b3[i][j];
      mPower (b1, addOn, b3);
}

for (i = 0; i < 3; i++)  {
    for (j = 0; j < 3; j++)
        aInverse[i][j] = identity[i][j] + addOn[i][j];
}

printIt(b3, aInverse);
system("pause");

return 0;
}

PrintIt, по крайней мере, сделал форматирование дисплея приятным, но меня беспокоят циклы for в моей основной проблеме. Я просто не могу уложить в голове, что делать. Я получаю значения повсюду, и мне действительно трудно понять математику, связанную с такой нишевой проблемой. Любая помощь будет оценена по достоинству.


c matrix matrix-inverse
person Dub    schedule 10.12.2014    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Я попытался исправить вашу логику (и получившийся код очень неэффективен и может быть написан намного лучше). Вот некоторые проблемы, которые я нашел:

  1. Переменная addOn не была инициализирована
  2. Ваш b был неправильным (плохое вычитание), учитывая определение B = I-A в вашем предисловии.
  3. Ваше накопление addOn = addOn + b3 не было в двойном цикле с использованием неинициализированных i и j.
  4. Кажется, у вас нет обратной проверки (a * a^-1)=I

(Обратите внимание, я получаю -0,000, потому что мы используем числа двойной точности IEEE, поэтому математика не точна, особенно когда речь идет о третях.)

#include <stdio.h>

void mPower(double x[3][3], double y[3][3], double z[3][3])
{
  int i, j, k;

  for (i = 0; i < 3; i++) 
    for (j = 0; j < 3; j++)   {
      z[i][j] = 0;
      for (k = 0; k < 3; k++)
        z[i][j] += (x[i][k] * y[k][j]);
}
}

void mAdd(double x[3][3], double y[3][3]) 
{
  int i, j;
  for (i = 0; i < 3; i++)  
    for (j = 0; j < 3; j++)
       x[i][j] += y[i][j];
}

void mCopy(double x[3][3], double y[3][3]) 
{
  int i, j;
  for (i = 0; i < 3; i++)  
    for (j = 0; j < 3; j++)
       x[i][j] = y[i][j];
}

void printIt(double a[3][3], double b[3][3])
{
 int i, j;

 for (i = 0; i < 3; i++)
 {
     printf ("\n\t\t\t");
     for (j = 0; j < 3; j++)
         printf ("%7.4f", a[i][j]);
     printf ("\t");
     for (j = 0; j < 3; j++)
         printf ("%7.4f", b[i][j]);
         printf ("\n\n\t");
 }
}

int main()
{
  int i, j, k;
  double a[3][3] = {.5, -1., 0., 0., .6666, 0., .5, -1., .6666};
  double b1[3][3] = {.5, 1., 0., 0., 1-.6666, 0., -.5, 1., 1-.6666};
  double b3[3][3] = {0};
  double addOn[3][3] = {0};
  double aInverse[3][3] = {0};
  double identity[][3] = {1., 0., 0., 0., 1., 0., 0., 0., 1.};
printIt(aInverse, b1);
  mAdd( aInverse, identity );
printIt(aInverse, b1);
  mCopy( b3, b1 );
  mCopy( addOn, b1 );
  for (k = 0; k < 20; k++)     {
    mAdd( aInverse, addOn );
    mPower(b3, b1, addOn);
    mCopy( b3, addOn );
printf( "loop iteration %d\n", k );
printIt(aInverse, addOn );
}
printf( "aInverse\n" );
printIt(b1, aInverse);
printf( "verify\n");
mPower(a,aInverse,b3);
printIt(a,b3);

return 0;
}

Результаты, которые я получил:

aInverse

         0.5000 1.0000 0.0000    2.0000 3.0003 0.0000


         0.0000 0.3334 0.0000    0.0000 1.5002 0.0000


        -0.5000 1.0000 0.3334   -1.5001 0.0000 1.5002

verify

         0.5000-1.0000 0.0000    1.0000-0.0000 0.0000


         0.0000 0.6666 0.0000    0.0000 1.0000 0.0000


         0.5000-1.0000 0.6666    0.0000 0.0000 1.0000
person George Houpis    schedule 10.12.2014
comment
Ух ты, ты выдул математику из воды. Я случайно разместил более старую версию с моей функцией mPower, поэтому addOn не был инициализирован, но помимо этого вы мне многое показали. Спасибо. Но как вы получаете -1,5001 и положительные 3,003? - person Dub; 10.12.2014
comment
Собственно, нашел! Первый 1. в b1 должен быть отрицательным, а -.5 должен быть положительным. - person Dub; 10.12.2014