Comments 11
Система уравнений же жесткая, почему к ней применяется метод Ньютона? Тот же неявный метод Адамса будет более предпочтителен
А можно источник, где написано, что система Лоренца жесткая? Неустойчивая — да, а вот жесткая?
Система Лоренца при классических параметрах системы является нежесткой. Определим собственные числа линейной части системы в пакете Maxima:
Результат выполнения:
Второй список — это кратность для соответствующего собственного числа. Далее вычислим коэффициент жесткости системы
Результат
A: matrix([-10,10,0], [28,-1,0], [0,0,-8/3]);
float(eigenvalues(%));
Результат выполнения:
[[-22.82772345116346,11.82772345116346,-2.666666666666667],[1.0,1.0,1.0]]Второй список — это кратность для соответствующего собственного числа. Далее вычислим коэффициент жесткости системы
float(22.82772345116346/2.666666666666667);Результат
8.560396294186297Ну только хорошо бы не линейную часть смотреть, а на линеаризованную (матрицу Якоби) в окрестности самого цикла. Но результат должен быть примерно таким же.
Получилось так (для начальной точки цикла, которую я определил):
Если взять другую точку этого цикла (которую Вишванат нашел):
A: matrix([-10,10,0], [28-27,-1,-(-2.147375914135)], [2.078143036771,-2.147375914135,-8/3]);
realpart(float(eigenvalues(%)));
[[-3.198123425968013,-10.4239475369641,-0.044595703734554],[1.0,1.0,1.0]]
float(10.4239475369641/0.044595703734554);
233.7433130108301Если взять другую точку этого цикла (которую Вишванат нашел):
A: matrix([-10,10,0], [28-27,-1,-(-13.763610682134)], [-19.578751942452,-13.763610682134,-8/3]);
realpart(float(eigenvalues(%)));
[[1.59085774810119,1.59085774810119,-16.84838216286904],[1.0,1.0,1.0]]
float(16.84838216286904/1.59085774810119);
10.59075343661548
Простите за глупый вопрос — почему просто не метод Рунге Кутта?
Потому что жесткие системы дифференциальных уравнений он «не возьмет». Методы Р-К, в основном, умеют хорошо решать пусть и сложные, многомерные, но нежесткие системы, где малая погрешность в начальных данных не вызовет сильного изменения решения.
Понимаете, у РК-методов фиксированный порядок точности. Есть, конечно, РК высших порядков, — я видел формулы до 8-ого, — но они там очень громоздкие. Чтобы увеличить точность получаемого численного решения, можно уменьшить шаг, но представление вещественных чисел не даст его сделать сколь угодно малым. В методе степенных рядов для систем с квадратичной правой частью по фазовым координатам можно получить простые рекуррентные закономерности вычисления коэффициентов ряда, и проблем нет, чтобы нарастить полином, аппроксимирующий ряд, на несколько членов.
А при каком числе гармоник решения вам перестает хватать стандартной двойной точности для метода Ньютона? В свое время при решении полиномиальных систем уравнений приходилось привлекать четверную и восьмерную точность для систем из ~30 уравнений, у вас же несравнимо больше. Или вы все же производите вычисления в длинной арифметике?
Я сейчас проверил — для 20-ти гармоник (в системе 124 уравнения) хватает стандартной двойной точности (хотя я там настраивал 30 знаков после точки). Вычисления можно производить с разным представлением вещественных чисел (ordinary floating point numbers или arithmetic on bigfloat numbers — так написано в официальной документации пакета) и разной точностью для метода Ньютона. Это делается так:
(fpprec:30,newtonepsilon:bfloat(10^(-fpprec+5)))$Sign up to leave a comment.
Ищем циклы на аттракторе Лоренца в пакете Maxima