*(F + (j * n + i)) = f (xmp[j], (i + 1));
}
}
// Матрица F готова. Надо вычислить по формуле:
// am = Inverse[Transpose[F].F].Transpose[F].ym значение
// коэффициентов a1, a2, a3, ...
// Транспонируем F
for (j = 0; j < n; j++) // Цикл по строкам TF
{
for (i = 0; i < yr; i++) // И по её столбцам
{
*(TF + (j * yr + i)) = *(F + (i * n + j));
}
}
// Считаем TMP = TF * F
MatrixMultiply (TF, n, yr, F, n, TMP);
// Далее считаем оперделитель от TMP
if ((dt = Determinant (TMP, n)) == 0)
{
printf ("\nТак, как определитель матрицы TF*F равен 0,\n"
"невозможно посчитать обратную к ним матрицу\n");
free (F); abort();
}
// Составляем обратную матрицу.
for (j = 0; j < n; j++)
{
for (i = 0; i < n; i++)
{
// Берём Минор элемента ij
MMinor (TMP, AC2, n, i, j);
// Знак элемента
flag = ((i + j) % 2 == 0) ? 1. : -1.;
// Сразу транспонирование
*(REV + (i * n) + j) = flag * Determinant (AC2, (n - 1)) / dt;
}
}
// Умножаем обратную матрицу на транспонированную к F
// т.е. Inverse (TF*F) * TF
// Такая матрица будет размера yr*n, поэтому вполне хватит памяти для F
MatrixMultiply (REV, n, n, TF, yr, F);
// И, наконец, всё это умножаем на матрицу Y и получаем искомые
// коэффициенты a1, a2, ... aN
// Для такой матрицы (размером 1*n) вполне хватит памяти под REV
MatrixMultiply (F, n, yr, ym, 1, REV);
// Всё, печатаем ответ
printf ("\nВычисления успешны, получен следующие коэффициенты:");
for (i = 0; i < n; i++)
printf ("\na%d = %lg", i, *(REV + i));
// Освободить память
free (F);
printf ("\nНажмите any key");
getch ();
printf ("\nDone.\n");
return 0;
}
После запуска программа сразу же начинает расчёт коэффициентов. На экран выводится следующее:
Программа аппроксимации функции методом наименьших квадратов для модели
y = a1x1 + a2x2^2 + a3x1x2
по заданной таблице эксперимента.
Разработчик: студент группы ПС-146
Нечаев Леонид Владимирович
25.02.2004
Известны результаты наблюдений:
x1 x2 y
1 1 0
-1 -1 -2
2 2 -2
3 -2 29
-2 4 54
Начинаем аппроксимацию...
Вычисления успешны, получены следующие коэффициенты:
a0 = 1
a1 = 2
a2 = -3
Нажмите any key
Done.
Результат верен, так как при подстановке переменных в модель получается верное равенство:
0 = 1 * 1 + 2 * 1 – 3 * 1 * 1
-2 = 1 * (-1) + 2 * (-1) – 3 * (-1) * (-1)
-2 = 1 * 2 + 2 * 2 – 3 * 2 * 2
29 = 1 * 3 + 2 * (-2) – 3 * 3 * (-2)
54 = 1 * (-2) + 2 * 4 – 3 * (-2) * 4
Задача решена верно.