Алгоритм решения задачи по численным методам

stm5539978-02

Задача.
Дана система линейных уравнений Ax=b с нижнетреугольной матрицей A. Уравнений много, скажем, больше 1000. Требуется написать алгоритм решения этой системы устойчивый к шумам в b.
Отцы численных методов , посоветуйте алгоритм или книжку какую-нибудь на эту тему!

gygo

ну так если она нижнетреугольная то из первого
уравнения можно сразу получить x1 и тд
зачем тут ЧМЫ?

stm5539978-02

Я даже больше скажу, на диагоналях этой нижнетреугольной матрицы стоят одинаковые числа!

gygo

ты лучше скажи зачем тебе тут
чмы когда здесь в явном виде все выражается?

stm5539978-02

Затем, что размер матрицы очень большой.

gygo

какая разница какой размер? хотя мб я чего то не понял
x1=b1/a11
x2=(b2-a21*x1)/a22
x3=(b3-a31*x1-a32*x2)/a33
и тд

griz_a

Типа ему не хочется производить 499500 вычитаний....

gygo

подозреваю что таким образом будет наименьшее число операций

griz_a

Я тоже
Но я в алгоритмах не больно-то шарю..

kachokslava

Метод нужен точный? Или итерационный сойдёт?
про методы Зейделя, релаксации, простой итерации почитай - они там при некоторых условиях хорошо сходятся. причём условия на матрицу, а не на правую часть.

sfmike

А смысл? Все равно там будет нужно умножать матрицу на вектор... А это ничуть не меньше операций чем было предложено. Тем более итераций может быть много.

kachokslava

Всё-таки Гаусс, даже с треугольной матрицей - это N^3. А итерационный метод сходится со скоростью геометрической прогрессии. пусть начальное приближение отстоит от точного ответа на 1000 (в некоторой норме каждый шаг это расстояние сокращается вдвое - за сколько шагов оно сойдётся до 1e-10 в этой норме? log(1E13)/log(2)~= 43 шага. а размерность матрицы >1000. сравни: 1000^3 и 43*1000^2

sfmike

Какой гаусс? Ты о чем? У него же матрица нижнетреугольная! N^2 операций.

kachokslava

ок, стормозил.. да - это у него обратный ход.

valex1971

разрулили

_NEKTO_

вроде гауссом большие СЛУ не модно решать из-за большой погрешности

a101

Не модно решать, потому что в треугольную с погрешностью переходим. А тут она уже треугольная, просто обратный ход. Он быстрый, и так как на диагонали уже одинаковые числа, то все будет довольно хорошо.

stm5539978-02

Проблема вот в чем.
Фиксируем А и x, вычисляем отклик b.
Затем
1) забываем, что мы знаем x и решаем систему Ax=b, применяя некоторый Алгоритм, получаем ответ x_1
2) зашумляем белым шумом вектор b, и снова решаем систему с помощью того же Алгоритма, получаем ответ x_2
Алгоритм нужен такой, чтобы в обоих случаях (x_1 и x_2) отклонения от x были минимально возможными.
Извините, что я с самого начала не четко поставил задачу.

gygo

смотри здесь
Самарский А.А. - Введение в численные методы стр 83.

Barmaglot

У тебя всегда есть оценка
 
|x2-x1|<|A^-1|*|b2-b1|
и это в отсутствии ошибок округления
Ошибки же решения будут всегда ничтожны при любом методе,
если у тебя конечно не матрица Гильберта
То есть если x2` x1` - их приближенные решения,то

|x2`-x1`|<|x2`-x2+x2-(x1`-x1+x1)|<2*eps+|A^-1|*|b2-b1|

stm5539978-02

Из полезного для себя в этой книжке нашел только постановку задачи, совпадающую с моей на странице 90 под номером 2.

stm5539978-02

Ты что-то доказываешь? Сформулируй что, так будет проще понять, что ты имеешь ввиду.

sfmike

Если у тебя хорошая матрица(а у тебя такая то маленькие погрешности в исходных данных приведут к маленьким погрешностям ответа.

stm5539978-02

т.е предлагается не изголяться и в качестве Алгортима использовать очевидный прямой метод?

sfmike

Да.

stm5539978-02

У меня сомнения по этому поводу.
Задача состоит в том (см выше чтобы в обоих случаях (зашумленного и незашумленного вектора b) решения x_1 и x_2 менее всего отклонялись от правильного ответа x. Если использовать прямой метод, то x_1=x, т.е r(x_1,x)=0 и это круто, но r(x_2,x) не будет нулем.
Почему не существует Алгоритма, который давал бы не точные ответы y_1 и y_2, но такие, что r(y_1,x) + r(y_2,x) < r(x,x_2) ?

stm5539978-02

Не, ребята! Тут все н так просто
И главной особенностью является огромный размер матрицы.
Когда я начинаю делать обратный ход (или как он там называется, то, что вы советуете то k-ая компонета x_1 все более отклоняется от k-ой компонеты x при возрастании k (c экспоненциальной скоростью). Т.е даже без всяких зашумлений решать ситему Ax=b простым подсчетом следующих компонент через предыдущие при больших размерах матрицы А нельзя! Надо что-то придумывать...

Barmaglot

А то, что оценка в моем последнем неравенстве фактически достижима. И как бы ты не хотел малости нормы |x2`-x1`|, она зависит только от свойств оператора A и не зависит от метода решения. То есть и при хорошей матрице эта величина может быть не малой.
Если же матрица плохо обусловлена, например 1 – на главной диагонали и a>1 – на побочной, то возмущение в одной лишь компоненте правой части

b2_n = b1_n+eps
приведет к такому возмущению решения при использовании метода прогонки

x2_n = x1_n+a^(n-1)*eps
теперь возьми a=2, n=1000 и посчитай
Это означает, что:
а) прямые методы дадут решение сколь угодно далекое от точного.
б) итерационные методы никогда не сойдутся

stm5539978-02

И что же делать?
Про матрицу А известно, что на диагоналях стоят одинаковые числа, те $a_{i,j}=z_{i-j}$ и эти z_i нормально распределены N(0,1).

gygo

Это означает, что:
а) прямые методы дадут решение сколь угодно далекое от точного.
б) итерационные методы никогда не сойдутся
вообще то число обусловленности влияет на скорость сходимости,
а не на саму сходимость(если я не ошибаюсь)

Barmaglot

Жаль, что ты сразу не выписал численные значения элементов стоящих на диагоналях. В моем примере с двумя диагоналями, если a<1, то как раз матрица очень хорошая.

stm5539978-02

понятно, но матрица, на самом деле плохая.... (большая!)

Barmaglot

возьми самый лучший показатель сходимости,

(M^(1/2)-m^(1/2/(M^(1/2)+m^(1/2
да посчитай сколько нужно итераций, если хотя бы [m,M]=[1,10^9]
На практике сходимости не будет по многим причинам:
1) Ни один компьютер не сможет это обсчитать
2) Сходимость быстро прекратится из-за ошибок округления, относительная невязка достигнет своего наименьшего значения и процесс начнет крутиться на одном месте.
3) Предположим мы совершили трудовой подвиг и дожали относительную невязку до 10^-16, но тогда ошибка может быть все равно сколь угодно большой, в силу известного неравенства связывающего ошибку, невязку и число обусловленности.

sfmike

Может тебе длинная арифметика поможет? gsl к примеру.

stm5539978-02

Что это такое? В смысле где посмотреть.

Barmaglot

Да?
Плохой мы называем матрицу с большим числом обусловленности,
а не высокой размерности

sfmike

Про gsl это я соврал. В ней этого нет. Смотри например тут: http://www.ginac.de/CLN/

stm5539978-02

Хоть горшком назови... не в этом суть.
Если матрица 20*20 или 30*30, то решая систему первым приходящим на ум методом, мы получим ответ с достаточной точностью (проверял в матлабе). Если же порядок системы 1000, то будет накапливаться машинная ошибка (или как там "вы" ее называете). Вопрос: как с этим можно бороться?

sfmike

Я тебе написал как с этим можно бороться

Barmaglot

Ёлы, да не должно там ничего копиться, если матрица хорошая, а если копится, то как раз и означает что плохая.
К тому же (см. мои посты выше даже решив ее длинной арифметикой (у меня юзалась gmp-4.1) или чем угодно, нужные тебе величины не будут маленькими.
Оставить комментарий
Имя или ник:
Комментарий: