NL: Файл lu.h


Функции
double	lu_decomp (double *A, size_t n, size_t p, int *sgn)
	-разложение матрицы.
void	lu_solve (double *LU, size_t n, size_t p, double *b)
	Решение системы линейных уравнений.
void	lu_improve (double *LU, size_t n, size_t p, double b, double x, double *r)
	Один шаг итерационного уточнения решения системы линейных уравнений.
double	lu_det (double **LU, size_t n, int sgn)
	Определитель матрицы.
void	lu_invert (double *LU, size_t n, size_t p, double **B)
	Обращение матрицы.
void	lu_mult_col (double *LU, size_t n, size_t p, double v, double res)
	Умножение матрицы , представленной массивами и , на вектор .

Подробное описание

Файл содержит функции, реализующие и использующие $LU$

-разложение квадратной матрицы.

Произвольная квадратная матрица $A$

порядка

допускает разложение вида $PA=LU$

, где

- матрица перестановок, $L$

- нижнетреугольная матрица с единичной диагональю, $U$

- верхнетреугольная матрица. Такое разложение используют для решения систем линейных уравнений $Ax=b$

, заменяя данную систему на пару систем $Ly=Pb$

Функции

double lu_decomp	(	double **	A,
		size_t	n,
		size_t *	p,
		int *	sgn
	)

$LU$ -разложение матрицы.

Функция находит $LU$ -разложение и оценку числа обусловленности квадратной матрицы $A$ порядка $n$ . На выходе диагональ и верхнетреугольная часть матрицы $A$ содержат элементы матрицы $U$ . Нижнетреугольная часть матрицы $A$ содержит поддиагональные элементы матрицы $L$ . Диагональные элементы матрицы $L$ равны $1$ . Матрица перестановок $P$ представлена вектором перестановок $p$ , в котором $p[j]$ - индекс ведущей строки на $j$ -ой итерации. На выходе $sgn$ - определитель матрицы $P$ . Функция возвращает нижнюю оценку числа обусловленности матрицы $A$ порядка $n$ . Как правило, данная оценка отличается от настоящего числа обусловленности не более чем в $10$ раз.

Для разложения используется алгоритм Гаусса с выбором главного элемента по столбцу. Для оценки числа обусловленности используется алгоритм, описанный в [FMM, KMN].

Трудоемкость: $\frac{2}{3}n^3 + O(n^2)$ операций с плавающей запятой

Примеры:: xlu.c.

double lu_det	(	double **	LU,
		size_t	n,
		int	sgn
	)

Определитель матрицы.

Функция возвращает определитель квадратной матрицы $A$ порядка $n$ . Ииспользуется предварительно найденное $LU$ -разложение матрицы A, представленное матрицей $LU$ и определителем $sgn$ матрицы перестановок $P$ (сам вектор перестановок не нужен), возвращаемыми функцией lu_decomp.

Трудоемкость: $2n + O(1)$ операций с плавающей запятой

void lu_improve	(	double **	LU,
		size_t	n,
		size_t *	p,
		double *	b,
		double *	x,
		double *	r
	)

Один шаг итерационного уточнения решения системы линейных уравнений.

Функция выполняет один шаг итерационного уточнения решения $x$ системы линейных уравнений $Ax=b$ . Используется предварительно найденное $LU$ -разложение матрицы $A$ , представленное матрицей $LU$ и вектором перестановок $p$ , возвращаемыми функцией lu_decomp. Возвращается улучшенное решение $x$ и вектор невязки $r=Ax-b$ .

Используется следующий алгоритм: вначале вычисляется вектор $r'=Ax-b$ , после этого с помощью функции lu_solve решается система $Ar'=b$ и, наконец, определяется улучшенное решение $x=x-r$ и вектор невязки $Ax-b$ .

Трудоемкость: $2n^2 + O(n)$ операций с плавающей запятой

void lu_invert	(	double **	LU,
		size_t	n,
		size_t *	p,
		double **	B
	)

Обращение матрицы.

Функция находит матрицу $B$ , обратную к матрице $A$ . Используется предварительно найденное $LU$ -разложение матрицы $A$ , представленное матрицей $LU$ и вектором перестановок $p$ , возвращаемыми функцией lu_decomp.

Трудоемкость: $2n^2 + O(n)$ операций с плавающей запятой

void lu_mult_col	(	double **	LU,
		size_t	n,
		size_t *	p,
		double *	v,
		double *	res
	)

Умножение матрицы $A$ , представленной массивами $LU$ и $p$ , на вектор $v$ .

Функция находит произведение $A\cdot v$ . Используется предварительно найденное $LU$ -разложение матрицы $A$ , представленное матрицей $LU$ и вектором перестановок $p$ , возвращаемыми функцией lu_decomp. Результат возвращается в векторе $res$ .

Трудоемкость: $2n^2 + O(n)$ операций с плавающей запятой

void lu_solve	(	double **	LU,
		size_t	n,
		size_t *	p,
		double *	b
	)

Решение системы линейных уравнений.

Функция находит решение квадратной системы линейных уравнений $Ax=b$ порядка $n$ , используя предварительно найденное $LU$ -разложение матрицы $A$ , представленное матрицей $LU$ и вектором перестановок $p$ , возвращаемыми функцией lu_decomp. Решение $x$ записывается на месте вектора правой части $b$ . Матрица $L$ U и вектор $p$ не меняются данной функцией и могут быть использованы для решения системы с другой правой частью.

Используется алгоритм прямой и обратной подстановки.

Трудоемкость: $2n^2 + O(n)$ операций с плавающей запятой

Примеры:: xlu.c.