|
Текст подпрограммы и версий afg4r_c.zip , afg4d_c.zip , afg4c_c.zip |
Тексты тестовых примеров tafg4r_c.zip , tafg4d_c.zip , tafg4c_c.zip |
Треугольное разложение вещественной матрицы общего вида методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности матрицы.
Для заданной квадратной вещественной матрицы А порядка N выполняется треугольная факторизация А = LU, где U - верхняя треугольная матрица, а L - нижняя треугольная матрица, и вычисляется величина RСОND, обратная числу обусловленности матрицы А:
rcond = 1 / || A ||1 * || A-1 ||1
Для вычисления || A- 1 ||1 решается система АZ = Y и полагается
|| A-1 ||1 ≈ || Z ||1 / || Y ||1 ,
где Y - решение системы АTY = Е; компоненты вектора Е выбираются таким образом, чтобы максимизировать норму || W ||1, где W - решение системы UTW = Е.
Вычисления проводятся в следующей последовательности:
1) вычисляется || A ||1 = max || aj ||1 , для j = 1, ..., N ,
где аj - j-й вектор-столбец матрицы A, а
|| аj ||1 = | a1 j | + | a2 j | +...+ | aN j | ;
2) выполняется треугольная факторизация А = LU ;
3) решается система UTW = E
одновременно с выбором компонент вектора E;
4) решается система LTY = W;
5) решается система LV = Y;
6) решается система UZ = V и полагается
rcond = || Y ||1 / ( || A ||1 * || Z ||1 )
Дж. Форсайт, М. Малькольм, К. Моулер. Машинные методы математических вычислений. М.: Мир, 1980.
int afg4r_c (real *a, integer *m, integer *n, integer *nlead,
real *rcond, real *z, integer *ierr)
Параметры
| a - | вещественный двумерный массив размера m*n, в котором задается исходная квадратная матрица порядка n, на выходе на соответствующих местах находятся элементы матрицы U и поддиагональные элементы элементарных матриц исключения Li, i = 1, ..., n - 1; |
| m - | первая размерность массива a в вызывающей программе (тип: целый); |
| n - | порядок матрицы A (тип: целый); |
| nlead - | целый вектор длины n, содержащий на выходе информацию о выполненых в процессе факторизации перестановках (см. замечания по использованию); |
| rcond - | вещественная переменная, содержащая на выходе вычисленное значение величины, обратной числу обусловленности матрицы A (см. замечания по использованию); |
| z - | вещественный рабочий вектор длины n (см. замечания по использованию); |
| ierr - | целая переменная, содержащая на выходе информацию о прохождении счета, при этом: |
| ierr=65 - | если значение m ≤ 0 или n ≤ 0; |
| ierr=66 - | если в процессе работы подпрограммы произошло переполнение (это говорит о том, что норма || A ||1 либо некоторые элементы матрицы U превосходят по абсолютной величине максимальное представимое на данной машине число); |
| ierr=-k - | если в результате факторизации диагональный элемент в K - й строке матрицы U равен нулю (это свидетельствует о вырожденности матрицы A). Если таких строк у матрицы несколько, то значение k полагается равным номеру последней из них. |
Версии
| afg4d_c - | треугольное разложение вещественной матрицы, заданной с удвоенной точностью, методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности матрицы. |
| afg4c_c - | треугольное разложение комплексной матрицы методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности матрицы. |
Вызываемые подпрограммы
| utafsi_c - | подпрограмма выдачи диагностических сообщений. |
Замечания по использованию
| 1. |
В подпрограмме afg4c_c массивы a и z имеют тип complex. | |
| 2. |
В подпрограмме afg4d_c массивы a, z и переменная RCOND имеют тип double. | |
| 3. |
На выходе K - й элемент вектора nlead равен номеру строки, переставленной на K - м шаге факторизации с K - й строкой матрицы A. Поскольку факторизация Гаусса требует n - 1 шагов, то nlead (n) = n. | |
| 4. |
Если переменной ierr присвоено значение, отличное от нуля, то величина RCOND полагается равной нулю, выдается соответствующее диагностическое сообщение и происходит выход из подпрограммы. | |
| 5. | Если ierr = 0, то на выходе вектор Z удовлетворяет условию || AZ ||1 = RCOND * || A ||1 * || Z ||1 . |
int main(void)
{
/* Initialized data */
static float a[16] /* was [4][4] */ = { 1.f,.42f,.54f,.66f,.42f,1.f,.32f,
.44f,.54f,.32f,1.f,.22f,.66f,.44f,.22f,1.f };
/* Local variables */
static int ierr;
extern int afg4r_c(float *, int *, int *, int *, float *, float *, int *);
static int i__, m, n, nlead[4];
static float z__[4], rcond;
#define a_ref(a_1,a_2) a[(a_2)*4 + a_1 - 5]
m = 4;
n = m;
for (i__ = 1; i__ <= 4; ++i__) {
printf("\n %16.7e %16.7e \n %16.7e %16.7e \n",
a_ref(i__, 1), a_ref(i__, 2), a_ref(i__, 3), a_ref(i__, 4));
}
afg4r_c(a, &m, &n, nlead, &rcond, z__, &ierr);
printf("\n %5i %5i %5i %5i \n", nlead[0], nlead[1], nlead[2], nlead[3]);
for (i__ = 1; i__ <= 4; ++i__) {
printf("\n %16.7e %16.7e \n %16.7e %16.7e \n",
a_ref(i__, 1), a_ref(i__, 2), a_ref(i__, 3), a_ref(i__, 4));
}
printf("\n %5i %16.7e \n", ierr, rcond);
return 0;
} /* main */
Для данной матрицы известны собственные значения:
λ1 = 2.32274
λ2 = 0.24226
λ3 = 0.63828
λ4 = 0.79670
Спектральное число обусловленности λ1 / λ2 ≈ 9.588
Результаты:
| 1.0 0.42 0.54 0.66 |
| -0.42 0.82360 0.09320 0.16280 |
a_ref = | -0.54 -0.11316 0.69785 -0.15482 |
| -0.66 -0.19767 0.22186 0.49787 |
nlead = (1, 2, 3, 4)
rcond = 0.09880