Прямые методы решения систем линейных алгебраических уравнений

Постановка задачи

Дана СЛАУ:

${\begin{cases}a_{11}x_{1}+a_{12}x_{2}+...+a_{1m}x_{m}=b_{1}\\\dots \\a_{m1}x_{1}+a_{m2}x_{2}+...+a_{mm}x_{m}=b_{m}\end{cases}}$

$A={\begin{pmatrix}a_{11}&\dots &a_{1m}\\\vdots &\ddots &\vdots \\a_{m1}&\dots &a_{mm}\end{pmatrix}}\quad x={\begin{pmatrix}x_{1}\\\vdots \\x_{m}\end{pmatrix}}\quad b={\begin{pmatrix}b_{1}\\\vdots \\b_{m}\end{pmatrix}}$

$Ax=b$ Предполагаем, что $A$ - невырожденная ($\Rightarrow$ задача корректна). $x^{*}=(x_{1}^{*}\dots x_{m}^{*})^{T}$ - приближенное решение задачи. Погрешность $e=x-x^{*}$. Невязка $r=b-Ax^{*}\Rightarrow r=Ax-Ax^{*}=A(x-x^{*})\Rightarrow e=A^{-1}r$

Определение. $\|x\|$ называется нормой вектора $x$, если обладает свойствами:

1. $\|x\|\geq 0$, $\|x\|=0\Leftrightarrow x=0$
2. $\|\alpha x\|=|\alpha |\|x\|\forall x\forall \alpha$
3. $\|x+y\|\leq \|x\|+\|y\|\forall x,y$

Три варианта задания нормы:

$\|x\|_{1}=\sum _{i=1}^{m}|x_{i}|,\quad \|x\|_{2}=(\sum _{i=1}^{m}|x_{i}|^{2})^{\frac {1}{2}},\quad \|x\|_{\infty }=\max {1\leq i\leq m}|x_{i}|$

Скалярное произведение: $(x,y)=x_{1}y_{1}+...+x_{m}y_{m}$. Погрешности: $\delta (x^{*})={\frac {\|x-x^{*}\|}{\|x\|}}$, $\Delta (x^{*})=\|x-x^{*}\|$

Определение. Сходимость по норме: $\{x^{(n)}\}_{n=1}^{\infty }$ - последовательность векторов $x^{(n)}\to x$ при $n\to \infty$, если $\Delta (x^{(n)})=\|x^{(n)}-x\|\to 0$ при $n\to \infty$.

Норма матрицы $\|A\|=\max {x\neq 0}{\frac {\|Ax\|}{\|x\|}}$. Свойства:

1. $\|A\|\geq 0\quad \|A\|=0\Leftrightarrow A=0$
2. $\|\alpha A\|=|\alpha |\|A\|\forall A,\alpha$
3. $\|A+B\|\leq \|A\|+\|B\|$
4. $\|A\cdot B\|\leq \|A\|\cdot \|B\|$
5. $\|Ax\|\leq \|A\|\cdot \|x\|\quad (\|x\|\neq 0\Rightarrow {\frac {\|Ax\|}{\|x\|}}\leq \|A\|-{\text{из определения}})$

$\|A\|_{1}=\max {1\leq j\leq m}\sum _{i=1}^{m}|a_{ij}|,\quad \|A\|_{2}=\max {1\leq j\leq m}{\sqrt {\lambda _{j}(A^{T}A)}},\quad \|A\|_{\infty }=\max {1\leq j\leq m}\sum _{i=1}^{m}|a_{j}i|$

Это подчиненные нормы $\|A\|_{E}={\sqrt {\sum _{i,j=1}^{m}\|a_{i}j\|^{2}}}$

Обусловленность задачи решения СЛАУ

Лемма. $\Delta (x^{*})\leq \|A^{-1}\|\cdot \|r\|$, где $r=b-Ax^{*}$

Доказательство. $x-x^{*}=e=A^{-1}r\Rightarrow \Delta (x^{*})=\|x-x^{*}\|=\|A^{-1}r\|\leq \|A^{-1}\|\cdot \|r\|$

Теорема. $x^{*}$ - точное решение системы $Ax^{*}=b^{*}$, где $b^{*}$ - приближенное к $b$.

Тогда верно: $\Delta (x^{*})\leq \nu _{\Delta }\Delta (b^{*})\quad \delta (x^{*})\leq \nu _{\delta }\delta (b*)$, где $\nu _{\Delta }=\|A^{-1}\|,\nu _{\delta }=\|A^{-1}\|\cdot \|b\|/\|x\|$

Доказательство. Здесь $r=b-Ax^{*}=b-b^{*}$. Из леммы: $\Delta (x^{*})\leq \|A^{-1}\|\cdot \Delta (b^{*})$. Разделим на $\|x\|$: ${\frac {\Delta (x^{*})}{\|x\|}}\leq {\frac {\|A^{-1}\|\cdot \|b\|}{\|x\|}}\cdot {\frac {\Delta (b^{*})}{\|b\|}}$

$\nu _{\Delta }=\|A^{-1}\|$ - абсолютное число обусловленности для $Ax=b$. $\nu _{\delta }=\nu {\delta }(x)={\frac {\|A^{-1}\|\cdot \|b\|}{\|x\|}}$ - естественное число обусловленности (относительное). $\max {x\neq 0}\nu _{\delta }(x)=\|A^{-1}\|\cdot \|A\|=\nu (A)=cond(A)$ - число обусловленности.

Теорема. $\delta (x^{*})\leq cond(A)\cdot \delta (b^{*})$. Факты:

1. $cond(E)=1\quad (E^{-1}=E,\|E\|=1)$
2. $cond(A)\geq 1\quad (E=A^{-1}A,\|E\|\leq \|A^{-1}\|\cdot \|A\|)$
3. $\forall \alpha \neq 0:cond(\alpha A)=cond(A)\quad (\|\alpha A\|\cdot \|\alpha A^{-1}\|=|\alpha |\|A\|\cdot |\alpha |^{-1}\|A^{-1}\|)$

Теорема. $x^{*}$ - точное решение. $A_{*}x^{*}=b$ - приближено заданная матрица $A_{*}$. Тогда $\delta ^{*}(x^{*})\leq cond(A)\cdot \delta (A_{*})$, где $\delta ^{*}(x^{*})={\frac {\|x-x^{*}\|}{\|x^{*}\|}},\delta (A_{*})={\frac {\|A-A_{*}\|}{\|A\|}}$

Доказательство. $r=b-Ax^{*}=A_{*}x^{*}-Ax^{*}=(A_{*}-A)x^{*}$

$\delta ^{*}(x^{*})={\frac {\|x-x^{*}\|}{\|x^{*}\|}}\leq \|A^{-1}\|\cdot \|(A_{*}-A)x^{*}\|/\|x^{*}\|\leq {\frac {\|a^{-1}\|\cdot \|A_{*}-A\|\cdot \|x^{*}\|}{\|x^{*}\|}}=cond(A)\cdot \delta (A_{*})$

Методы решения

Метод Гаусса

Прямой ход:

1-й шаг - пусть $a_{11}\neq 0$. Найдем $\mu _{i1}={\frac {a_{i1}}{a{11}}},i={\overline {2,m}}$. Из $2..m$ уравнений вычитаем первое, домноженное на $\mu _{i1}$. Получаем:

${\begin{matrix}a_{11}x_{1}+a_{12}x_{2}+\dots +a_{1m}x_{m}=b_{1}\\a_{22}^{(1)}x_{2}+\dots +a_{2m}^{(1)}x_{m}=b_{2}^{(1)}\\\dots \\a_{m2}^{(1)}x_{2}+\dots +a_{mm}^{(1)}x_{m}=b_{m}^{(1)}\end{matrix}}$

$a_{ij}^{(1)}=a_{ij}-\mu _{i1}a_{ij},\quad b_{i}^{(1)}=b_{i}-\mu _{i1}b_{1}$

2-й шаг - пусть $a_{22}^{(1)}\neq 0,\quad \mu _{i2}={\frac {a_{i2}^{(1)}}{a_{22}^{(1)}}}$

k-й шаг: $a_{kk}^{(k-1)}\neq 0,\quad \nu _{ik}={\frac {a_{ik}^{(k-1)}}{a_{kk}^{(k-1)}}},i={\overline {k+1,m}}$

За (m-1) шаг получаем:

${\begin{matrix}a_{11}x_{1}+a_{12}x_{2}+a_{33}x_{3}+\dots +a_{1m}x_{m}=b_{1}\\a_{22}^{(1)}x_{2}+a_{33}^{(1)}x_{3}+\dots +a_{2m}^{(1)}x_{m}=b_{2}^{(1)}\\a_{33}^{(2)}x_{3}+\dots +a_{2m}^{(2)}x_{m}=b_{2}^{(2)}\\\dots \\a_{mm}^{(m-1)}x_{m}=b_{m}^{(m-1)}\end{matrix}}$

Обратный ход:

Из последнего уравнения находим $x_{m}$, подставляем в предпоследнее, находим $x_{m-1}$. $x_{m}=b_{m}^{(m-1)}/a_{mm}^{(m-1)},x_{k}=(b_{k}^{(k-1)}-a_{k,k+1}^{(k-1)}x_{k+1}-...-a_{km}^{(k-1)}x_{m})/a_{kk}^{k-1},k={\overline {m-1,1}}$

Трудоемкость:

1. m-1 деление, (m-1)m умножений, (m-1)m вычитаний $\Rightarrow Q_{1}=2(m-1)^{2}+3(m-1)$

k. $Q_{k}=2(m-k)^{2}+3(m-k)$

$Q=\sum _{k=1}^{n-1}Q_{k}=2\sum _{k=1}^{m-1}(m-k)^{2}+3\sum _{k=1}^{m-1}(m-k)=2\sum _{k=1}^{m-1}k^{2}+3\sum _{k=1}^{m-1}k={\frac {2(m-1)m(2m-1)}{6}}+{\frac {3(m-1)m}{2}}\approx {\frac {2}{3}}m^{3}$

Обратный ход: $m^{2}$ операций. Необходимость выбора главных элементов: главный элемент не должен быть равен нулю. Если он близок к нулю, растет погрешность (это для ЭВМ).

Матричная форма метода Гаусса. (LU - разложение)

После первого шага получаем $A^{(1)}x=b^{(1)}$. Введем:

$M_{1}={\begin{pmatrix}1&0&0&\dots &0\\-\mu _{21}&1&0&\dots &0\\-\mu _{31}&0&1&\dots &0\\\dots \\-\mu _{m1}&0&0&\dots &1\end{pmatrix}}M_{2}={\begin{pmatrix}1&0&0&\dots &0\\0&1&0&\dots &0\\0&-\mu _{32}&1&\dots &0\\\dots \\0&-\mu _{m2}&0&\dots &1\end{pmatrix}}$

$M_{m-1}={\begin{pmatrix}1&0&\dots &0&0\\0&1&\dots &0&0\\0&0&\dots &0&0\\\dots \\0&0&\dots &-\mu _{m,m-1}&1\end{pmatrix}}$

$A^{(1)}=M_{1}A,b^{(1)}=M_{1}b,\quad A^{(m-1)}=M_{m-1}A^{(m-2)},b^{(m-1)}=M_{m-1}b^{(m-2)}$

Начальная система: $M_{m-1}\cdot ...\cdot M_{2}\cdot M_{1}\cdot Ax=M_{m-1}\cdot ...\cdot M_{2}\cdot M_{1}\cdot b$.

$M_{m-1}\cdot ...\cdot M_{2}\cdot M_{1}\cdot A$ - верхняя треугольная матрица $=A^{(m-1)}$.

$L=M_{1}^{-1}M_{2}^{-1}...M_{m-1}^{-1}={\begin{pmatrix}1&0&0&...&0\\\mu _{21}&1&0&...&0\\\mu _{31}&\mu _{32}&1&...&0\\\dots \\\mu _{m1}&\mu _{m2}&\mu _{m3}&...&1\end{pmatrix}}$

Свойства элементарных матриц:

1. $\det M_{i}=1$
2. $M_{i}$ - нижние треугольные матрицы
3. $M_{i}^{-1}$ - тоже, только у $\mu$ вместо "-" стоит "+"

$A=LU$

Теорема. Если все главные миноры матрицы $A\neq 0$, то $\exists !$ нижняя треугольная матрица $L$ вида $(1)$ и верхняя треугольная матрица $U$ такие, что $A=LU$.

Использование:

1. преобразуем $b$ в $b^{(m-1)}\quad b^{(m-1)}=M_{m-1}...M_{2}M_{1}b$
2. Решаем систему $Ux=b^{(m-1)}$

Количество действий: ${\frac {2}{3}}m^{3}+2m^{2}$

Модификации Гаусса

1. на k-м шаге в качестве главного элемента выбираем максимальный по модулю коэффициент $a_{i_{k},k}$ и переставляем строки.
2. на каждом шаге выбираем максимальный по модулю коэффициент из строк с номерами $i=k..m$ и переставляем строки (столбцы)

Масштабирование - делим каждое уравнение на наибольший по модулю коэффициент.

Метод Холецкого

$Ax=b$, $A$ - симметричная ($A=A^{T}$), положительно определенная ($\forall x\neq 0:(Ax,x)>0$ - скалярное произведение $Ax$ на $x$) Проверка:

1. Все главные миноры положительны
2. $|a_{ii}|>\sum _{j\neq i}|a_{ij}|$ - диагональное преобладание
3. $|a_{jj}|>\sum _{i\neq j}|a_{ij}|$ - столбцовое преобладание

Специальное $LU$ - разложение: $A=LL^{T}$.

$L={\begin{pmatrix}l_{11}&0&\dots &0\\l_{21}&l_{22}&\dots &0\\\dots \\l_{m1}&l_{m2}&\dots &l_{mm}\end{pmatrix}}$

${\begin{matrix}&l_{11}^{2}=a_{11}&\\&l_{i1}l_{11}=a_{i1}&i={\overline {2,m}}\\&l_{21}^{2}+l_{22}^{2}=a_{i2}&\\&l_{i1}l_{21}+l_{i2}l_{22}=a_{i2}&i={\overline {3,m}}\\&\dots &\\&l_{k1}^{2}+l_{k2}^{2}+...+l_{kk}^{2}=a_{kk}&\\&l_{i1}l_{k1}+...+l_{ik}l_{kk}=a_{ik}&i={\overline {k+1,m}}\\&\dots &\\&l_{m1}^{2}+...+l_{mm}^{2}=a_{mm}&\end{matrix}}\Rightarrow$

$\Rightarrow {\begin{matrix}&l_{11}={\sqrt {a_{11}}}&\\&l_{i1}=a_{i1}/l_{11}&i={\overline {2,m}}\\&l_{22}={\sqrt {a_{i2}-l_{21}^{2}}}&\\&l_{i2}=(a_{i2}-l_{i1}l_{21})/l_{22}&i={\overline {3,m}}\\&\dots &\\&l_{kk}={\sqrt {a_{kk}-l_{k1}^{2}-...-l_{k,k-1}^{2}}}&\\&l_{ik}=(a_{ik}-l_{i1}l_{k1}-...-l_{i,k-1}l_{k,k-1})/l_{kk}&i={\overline {k+1,m}}\\&\dots &\\&l_{mm}={\sqrt {a_{mm}-l_{m1}^{2}-...-l_{m,m-1}^{2}}}&\end{matrix}}$

После получения этого разложения решаем $Ly=b,L^{T}x=y$. Это решение требует $2m^{2}$ действий

Метод Прогонки

${\begin{matrix}&b_{1}x_{1}+c_{1}x_{2}&=d_{1}\\&a_{2}x_{1}+b_{2}x_{2}+c_{2}x_{3}&=d_{2}\\&\dots \\&a_{i}x_{i-1}+b_{i}x_{i}+c_{i}x_{i+1}&=d_{i}\\&\dots \\&a_{m-1}x_{m-2}+b_{m-1}x_{m-1}+c_{m-1}x_{m}&=d_{m-1}\\&a_{m}x_{m-1}+b_{m}x_{m}&=d_{m}\\\end{matrix}}$

Вывод:

${\begin{matrix}x_{1}=\alpha _{1}x_{2}+\beta _{1}\alpha _{1}=-{\frac {c_{1}}{b_{1}}}\beta _{1}={\frac {d_{1}}{b_{1}}}\\x_{2}=\alpha _{2}x_{3}+\beta _{2}\alpha _{2}=-{\frac {c_{2}}{b_{2}+a_{2}\alpha _{1}}}\beta _{2}={\frac {d_{2}-a_{2}\beta _{1}}{b_{2}+a_{2}\alpha _{1}}}\\\dots \\x_{i}=\alpha _{i}x_{i+1}+\beta _{i}\alpha _{i}=-{\frac {c_{i}}{b_{i}+a_{i}\alpha _{i-1}}}\beta _{i}={\frac {d_{i}-a_{i}\beta _{i-1}}{b_{i}+a_{i}\alpha _{i-1}}}\\\dots \\x_{m-1}=\alpha _{m-1}x_{m}+\beta _{mm-1}\end{matrix}}$

$m:a_{m}(a_{m-1}x_{m}+\beta _{m-1})+b_{m}x_{m}\Rightarrow x_{m}=\beta _{m}={\frac {d_{m}-a_{m}\beta _{m-1}}{b_{m}+a_{m}\alpha _{m-1}}}$

Алгоритм

Прямой ход:

${\begin{aligned}\alpha _{i}&=-{\frac {c_{i}}{\gamma _{i}}}\beta _{i}={\frac {d_{i}}{\gamma _{i}}}\gamma _{i}=b_{i},&i=1;&\\\alpha _{i}&=-{\frac {c_{i}}{\gamma _{i}}}\beta _{i}={\frac {d_{i}-a_{i}\beta _{i-1}}{\gamma _{i}}}\gamma _{i}=b_{i}+a_{i}\alpha _{i-1},&i=2..m-1;&\\\beta _{i}&={\frac {d_{i}-a_{i}\beta _{i-1}}{\gamma _{i}}}\gamma _{i}=b_{i}+a_{i}\alpha _{i-1},&i=m.&\end{aligned}}$

Обратный ход:

${\begin{matrix}i=m&x_{m}=b_{m}\\i=m-1..1&x_{i}=a_{i}x_{i+1}+\beta _{i}\\\end{matrix}}$

Количество операций: $Q=8m$.

Теорема. Пусть $|b_{k}|\geq |a_{k}|+|c_{k}|,|b'_{k}|>|a_{k}|,(1\leq k\leq m)$. Тогда $\gamma _{i}\neq 0$ и $\alpha _{i}\leq 1$ $\forall i=1,2,...,m$

Доказательство. $\gamma _{1}=b_{1}\neq 0,|\alpha _{1}|={\frac {|c_{1}|}{|b_{1}|}}\leq 1$

Пусть $\gamma _{k-1}\neq 0$ и $|\lambda _{k-1}|\leq 1$ для некоторого $k>1$. Тогда $|\gamma _{k}|=|b_{k}+a_{k}\alpha _{k-1}|\geq |b_{k}|-|a_{k}||\alpha _{k-1}|\geq |b_{k}|-|a_{k}|>0$, $|\gamma _{k}|\geq |c_{k}|\Rightarrow \gamma _{k}\neq 0$ и $|\alpha _{k}|={\frac {|c_{k}|}{|\gamma _{k}|}}\leq 1$