Использование Dfun (якобиана) в одеинте Python - PullRequest
0 голосов
/ 10 марта 2020

Это довольно общий вопрос, часть которого, вероятно, относится к любому численному моделированию связанных ODE, а часть может относиться только к методу odeint в библиотеке scipy.integrate Python.

Во-первых, как odeint использует вручную введенный якобиан (параметр Dfunc) и почему он так сильно ускоряет работу больших систем ODE?

Во-вторых, и более уместно для моего указания c проблема, если функция Якоби слегка неверна, odeint даст неправильное решение или просто замедлит его? На глаз (анимация результата симуляции) я не могу обнаружить разницу; Я надеюсь, что это потому, что якобиан прав, но я не могу быть полностью уверен.

1 Ответ

2 голосов
/ 11 марта 2020

LSODA, метод ODEPACK, вызываемый odeint, использует неявные многошаговые методы для решения y'=f(t,y). Это означает, что на каждом шаге необходимо решать нелинейную систему уравнений, по существу, уравнение с фиксированной точкой

y[i+1] = h*b[0]*f(t[i+1],y[i+1]) + C

, где C - линейная комбинация предыдущих точек y[i-k] и значений f(t[i-k],y[i-k]), k=0,..,q, которая является константой в текущем шаге для вычисления y[i+1].

Теперь в любом уравнении с фиксированной точкой y=g(y), если оно сжимается в интересующей области, вы можете выполнить итерацию это и найти фиксированную точку yfix в качестве предела. Сходимость будет линейной, когда коэффициент сжатия определяется в основном нормой / спектральным радиусом якобиана g'(yfix). Теперь представьте, что разбиение g известно на линейную часть и нелинейный остаток (возможно, с небольшой линейной частью),

g(y) = A*y + r(y) = yfix + A*(y-yfix) + [r(y)-r(yfix)]

, где «маленький» относится к последнему члену в Последнее разложение. Линейная часть в уравнении с фиксированной точкой может быть перенесена на другую сторону, так что получается новое уравнение с фиксированной точкой

y = (I-A)^(-1) * r(y)

. Его коэффициент сходимости теперь зависит от расстояния до фиксированной точки: чем ближе вы подходите, тем быстрее сходимость. Это может иметь ненулевую нижнюю границу и не быть квадратичным c, но это будет быстрее, чем исходная итерация.

...