Метод локальной линеаризации

В численном анализе метод локальной линеаризации (ЛЛ) представляет собой общую стратегию разработки числовых интеграторов для дифференциальных уравнений, основанную на локальной (кусочной) линеаризации данного уравнения на последовательных интервалах времени. Затем числовые интеграторы итеративно определяются как решение полученного кусочно-линейного уравнения в конце каждого последовательного интервала. Метод ЛЛ был разработан для различных уравнений, таких как обыкновенные , запаздывающие , случайные и стохастические дифференциальные уравнения. Интеграторы LL являются ключевым компонентом в реализации методов вывода для оценки неизвестных параметров и ненаблюдаемых переменных дифференциальных уравнений с учетом временных рядов (потенциально зашумленных) наблюдений. Схемы LL идеально подходят для работы со сложными моделями в различных областях, таких как нейронаука , финансы , управление лесным хозяйством , техника управления , математическая статистика и т. д.

Предыстория

Дифференциальные уравнения стали важным математическим инструментом для описания эволюции во времени ряда явлений, например, вращения планет вокруг Солнца, динамики цен активов на рынке, пожара нейронов, распространения эпидемий и т. д. Однако поскольку точные решения этих уравнений обычно неизвестны, необходимы численные аппроксимации к ним, полученные с помощью численных интеграторов. В настоящее время многие приложения в инженерных и прикладных науках, ориентированные на динамические исследования, требуют разработки эффективных численных интеграторов, максимально сохраняющих динамику этих уравнений. С этой основной мотивацией были разработаны интеграторы локальной линеаризации.

линеаризации высокого

Метод локальной линеаризации высокого порядка (HOLL) представляет собой обобщение метода локальной линеаризации, ориентированное на получение интеграторов высокого порядка для дифференциальных уравнений, сохраняющих устойчивость и динамику линейных уравнений. Интеграторы получаются путем разделения на последовательных интервалах времени решения x исходного уравнения на две части: решение z локально линеаризованного уравнения плюс аппроксимация невязки высокого порядка. $\mathbf {r} =\mathbf {x} -\mathbf {z}$.

Схема локальной линеаризации

Схема локальной линеаризации (LL) — это окончательный рекурсивный алгоритм , который позволяет численно реализовать дискретизацию , полученную на основе метода LL или HOLL, для класса дифференциальных уравнений.

Методы LL для ODE

Рассмотрим d -мерное обыкновенное дифференциальное уравнение (ОДУ)

${\frac {d\mathbf {x} \left(t\right)}{dt}}=\mathbf {f} \left(t,\mathbf {x} \left(t\right)\right),\qquad t\in \left[t_{0},T\right],\qquad \qquad \qquad \qquad (4.1)$

с начальным состоянием $\mathbf {x} (t_{0})=\mathbf {x} _{0}$, где $\mathbf {f}$ является дифференцируемой функцией.

Позволять $\left(t\right)_{h}=\{t_{n}:n=0,..,N\}$ быть временной дискретизацией временного интервала $[t_{0},T]$ с максимальным размером шага h таким, что $t_{n}<t_{n+1}$ и $h_{n}=t_{n+1}-t_{n}\leq h$. После локальной линеаризации уравнения (4.1) на шаге по времени $t_{n}$ формула вариации констант дает

$\mathbf {x} (t_{n}+h)=\mathbf {x} (t_{n})+\mathbf {\phi } (t_{n},\mathbf {x} (t_{n});h)+\mathbf {r} (t_{n},\mathbf {x} (t_{n});h),$

где

$\mathbf {\phi } (t_{n},\mathbf {z} _{n};h)=\int \limits _{0}^{h}e^{\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})(h-s)}(\mathbf {f} (t_{n},\mathbf {z} _{n})+\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})s)\,ds\qquad$

результат линейного приближения, и

$\mathbf {r} (t_{n},\mathbf {z} _{n};h)=\int \limits _{0}^{h}e^{\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})(h-s)}\mathbf {g} _{n}(s,\mathbf {x} (t_{n}+s))\,ds,\qquad \qquad \qquad (4.2)$

– остаток линейного приближения. Здесь, $\mathbf {f} _{\mathbf {x} }$ и $\mathbf {f} _{t}$ обозначают частные производные f по переменным x и t соответственно, и $\mathbf {g} _{n}(s,\mathbf {u} )=\mathbf {f} (s,\mathbf {u} )-\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})\mathbf {u} -\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})(s-t_{n})-\mathbf {f} (t_{n},\mathbf {z} _{n})+\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})\mathbf {z} _{n}.$

Локальная

Для временной дискретизации $\left(t\right)_{h}$, Локальная линейная дискретизация ОДУ (4.1) в каждой точке $t_{n+1}\in \left(t\right)_{h}$ определяется рекурсивным выражением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h_{n}),\qquad {\text{ with }}\quad \mathbf {z} _{0}=\mathbf {x} _{0}.\qquad \qquad \qquad \qquad (4.3)$

Локальная линейная дискретизация (4.3) сходится с порядком 2 к решению нелинейных ОДУ, но соответствует решению линейных ОДУ. Рекурсия (4.3) также известна как экспоненциальная дискретизация Эйлера.

высокого

Для временной дискретизации $(t)_{h},$ точке локальная линейная (ХОЛЛ) дискретизация ОДУ (4.1) высокого порядка в каждой $t_{n+1}\in (t)_{h}$ определяется рекурсивным выражением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h_{n})+{\widetilde {\mathbf {r} }}(t_{n},\mathbf {z} _{n};h_{n}),\qquad {\text{ with }}\quad \mathbf {z} _{0}=\mathbf {x} _{0},\qquad \qquad \qquad (4.4)$

где ${\tilde {r}}$ это заказ $\alpha$ (> 2 ) аппроксимация невязки r $(i.e.,\left\vert \mathbf {r} (t_{n},\mathbf {z} _{n};h)-{\widetilde {\mathbf {r} }}(t_{n},\mathbf {z} _{n};h)\right\vert \propto h^{\alpha +1}).$ Дискретизация ХОЛЛА (4.4) сходится с порядком $\alpha$ к решению нелинейных ОДУ, но соответствует решению линейных ОДУ.

Дискретизация ХОЛЛА может быть получена двумя способами: 1) (квадратурный) путем аппроксимации интегрального представления (4.2) функции r ; и 2) (на основе интегратора) с использованием числового интегратора для дифференциального представления r, определенного формулой

${\frac {d\mathbf {r} (t)}{dt}}=\mathbf {q} (t_{n},\mathbf {z} _{n};t\mathbf {,\mathbf {r} } (t)\mathbf {)} ,\qquad {\text{ with }}\qquad \mathbf {r} (t_{n})=\mathbf {0,} \qquad \qquad \qquad (4.5)$

для всех $t\in \lbrack t_{k},t_{k+1}]$, где

$\mathbf {q} (t_{n},\mathbf {z} _{n};s\mathbf {,\xi } )=\mathbf {f} (s,\mathbf {z} _{n}+\mathbf {\phi } \left(t_{n},\mathbf {z} _{n};s-t_{n}\right)+\mathbf {\xi } )-\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})\mathbf {\phi } (t_{n},\mathbf {z} _{n};s-t_{n})-\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})(s-t_{n})-\mathbf {f} (t_{n},\mathbf {z} _{n}).$

Дискретизация ХОЛЛА, например, следующая:

• Локально линеаризованная дискретизация Рунге-Кутты

$\qquad \mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h_{n})+h_{n}\sum _{j=1}^{s}b_{j}\mathbf {k} _{j},\quad {\text{ with }}\quad \mathbf {k} _{i}=\mathbf {q} (t_{n},\mathbf {z} _{n};{\text{ }}t_{n}+c_{i}h_{n}\mathbf {,} \mathbf {} h_{n}\sum _{j=1}^{i-1}a_{ij}\mathbf {k} _{j}),$

которое получается в результате решения (4.5) по s-этапной явной схеме Рунге–Кутты (РК) с коэффициентами $\mathbf {c} =\left[c_{i}\right],\mathbf {A} =\left[a_{ij}\right]\quad and\quad \mathbf {b} =\left[b_{j}\right]$.

• Локальная линейная дискретизация Тейлора

$$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h_{n})+\int _{0}^{h_{n}}e^{(h_{n}-s)\mathbf {f} _{\mathbf {x} }\left(t_{n},\mathbf {z} _{n}\right)}\sum _{j=2}^{p}{\frac {\mathbf {c} _{n,j}}{j!}}s^{j}\,ds,{\text{ with }}\mathbf {c} _{n,j}=\left({\frac {d^{j+1}\mathbf {x} (t)}{dt^{j+1}}}-\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n}){\frac {d^{j}\mathbf {x} (t)}{dt^{j}}}\right)\mid _{t=\mathbf {z} _{n}},$$

что является результатом аппроксимации $\mathbf {g} _{n}$в (4.2) его порядка p усеченным разложением Тейлора .

• Дискретизация экспоненциального распространения многошагового типа

${\textstyle \mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h)+h\sum _{j=0}^{p-1}\gamma _{j}\nabla ^{j}\mathbf {g} _{n}(t_{n},\mathbf {z} _{n}),\quad with\quad \gamma _{j}=(-1)^{j}\int \limits _{0}^{1}e^{(1-\theta )h\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})}\left({\begin{array}{c}-\theta \\j\end{array}}\right)d\theta ,}$

что является результатом интерполяции $\mathbf {g} _{n}$в (4.2) полиномом степени p на $t_{n},\ldots ,t_{n-p+1}$, где $\nabla ^{j}\mathbf {g} _{n}(t_{m},\mathbf {z} _{m})$ обозначает j -ю обратную разность $\mathbf {g} _{n}(t_{m},\mathbf {z} _{m})$.

• Дискретизация экспоненциального распространения типа Рунге Кутты

${\textstyle \mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h)+h\sum _{j=0}^{p-1}\gamma _{j,p}\nabla ^{j}\mathbf {g} _{n}(t_{n},\mathbf {z} _{n}),\quad {\text{ with }}\quad \gamma _{j,p}=\int \limits _{0}^{1}e^{(1-\theta )h\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})}\left({\begin{array}{c}\theta p\\j\end{array}}\right)d\theta ,}$

что является результатом интерполяции $\mathbf {g} _{n}$в (4.2) полиномом степени p на $t_{n},\ldots ,t_{n}+(p-1)h/p$,

• Линеализованная экспоненциальная дискретизация Адамса

${\textstyle \mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h)+h\sum _{j=1}^{p-1}\sum _{l=1}^{j}{\frac {\gamma _{j+1}}{l}}\nabla ^{l}\mathbf {g} _{n}(t_{n},\mathbf {z} _{n}),\quad {\text{ with }}\quad \gamma _{j+1}=(-1)^{j+1}\int \limits _{0}^{1}e^{(1-\theta )h\mathbf {f} _{\mathbf {x} }\left(t_{n},\mathbf {z} _{n}\right)}\theta \left({\begin{array}{c}-\theta \\j\end{array}}\right)d\theta ,}$

что является результатом интерполяции $\mathbf {g} _{n}$в (4.2) полиномом Эрмита степени p на $t_{n},\ldots ,t_{n-p+1}$.

Схемы локальной линеаризации

Вся численная реализация $\mathbf {y} _{n}$ дискретизации LL (или HOLL) $\mathbf {z} _{n}$ включает приближения ${\widetilde {\phi }}_{j}$ к интегралам $\phi _{j}$ формы

$\phi _{j}(\mathbf {A} ,h)=\int \limits _{0}^{h}e^{(h-s)\mathbf {A} }s^{j-1}\,ds,\qquad j=1,2\ldots ,$

где A — матрица размера d × d . Каждая численная реализация $\mathbf {y} _{n}$ LL (или HOLL) $\mathbf {z} _{n}$ любого порядка обычно называется схемой локальной линеаризации .

Вычисление интегралов с

Среди ряда алгоритмов вычисления интегралов $\phi _{j}$предпочтительны те, которые основаны на аппроксимациях рациональных подпространств Паде и Крылова для экспоненциальной матрицы. При этом центральную роль играет выражение

$\sum \nolimits _{i=1}^{l}\phi _{i}(\mathbf {A} ,h)\mathbf {a} _{i}=\mathbf {L} e^{h\mathbf {H} }\mathbf {r} ,$

где $\mathbf {a} _{i}$ являются d -мерными векторами,

$\mathbf {H} ={\begin{bmatrix}\mathbf {A} &\mathbf {v} _{l}&\mathbf {v} _{l-1}&\cdots &\mathbf {v} _{1}\\\mathbf {0} &\mathbf {0} &1&\cdots &0\\\mathbf {0} &\mathbf {0} &0&\ddots &0\\\vdots &\vdots &\vdots &\ddots &1\\\mathbf {0} &\mathbf {0} &0&\cdots &0\end{bmatrix}}\in \mathbb {R} ^{(d+l)\times (d+l)},$

$\mathbf {L} =[\mathbf {I} \quad \mathbf {0} _{d\times l}]$, $\mathbf {r} =[\mathbf {0} _{1\times (d+l-1)}\quad 1]^{\intercal },$ $\mathbf {v} _{i}=\mathbf {a} _{i}(i-1)!$, существование $\mathbf {I}$ -мерная единичная матрица d .

Если $\mathbf {P} _{p,q}(2^{-k}\mathbf {H} h)$ обозначает ( p ; q ) Паде -аппроксимацию $e^{2^{-k}\mathbf {H} h}$ и k — наименьшее натуральное число такое, что $|2^{-k}\mathbf {H} h|\leq {\frac {1}{2}},then$

$\left\vert \sum \nolimits _{i=1}^{l}\phi _{i}(\mathbf {A} ,h)\mathbf {a} _{i}-\mathbf {L} \left(\mathbf {\mathbf {P} } _{p,q}(2^{-k}\mathbf {H} h)\right)^{2^{k}}\mathbf {r} \right\vert \varpropto h^{p+q+1}.$

Если $\mathbf {\mathbf {k} } _{m,k}^{p,q}(h,\mathbf {H} ,\mathbf {r} )$ обозначает (m; p; q; k) Крылова–Паде аппроксимацию $e^{h\mathbf {H} }\mathbf {r}$, затем

$\left\vert \sum \nolimits _{i=1}^{l}\phi _{i}(\mathbf {A} ,h)\mathbf {a} _{i}-\mathbf {L\mathbf {k} } _{m,k}^{p,q}(h,\mathbf {H} ,\mathbf {r} )\right\vert \varpropto h^{\min({m,p+q+1})},$

где $m\leq d$ — размерность подпространства Крылова.

Схемы ЛЛ Ордена-2

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r,}$ $\qquad \qquad (4.6)$

где матрицы $\mathbf {M} _{n}$, L и r определяются как

$\mathbf {M} _{n}={\begin{bmatrix}\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n})&\mathbf {f} _{t}(t_{n},\mathbf {y} _{n})&\mathbf {f} (t_{n},\mathbf {y} _{n})\\0&0&1\\0&0&0\end{bmatrix}}\in \mathbb {R} ^{(d+2)\times (d+2)},$

$\mathbf {L} =\left[{\begin{array}{ll}\mathbf {I} &\mathbf {0} _{d\times 2}\end{array}}\right]$ и $\mathbf {r} ^{\intercal }=\left[{\begin{array}{ll}\mathbf {0} _{1\times (d+1)}&1\end{array}}\right]$ с $p+q>1$ . Для больших систем ОДУ

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {M} _{n},\mathbf {r} )\mathbf {,} \qquad {\text{ with }}\qquad m_{n}>2.$

Схемы Ордена-3 LL-Тейлора

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L} _{1}(\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {T} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r} _{1}\mathbf {,}$ $\qquad \qquad (4.7)$

где для автономных ОДУ матрицы $\mathbf {T} _{n},\mathbf {L} _{1}$ и $\mathbf {r} _{1}$ определяются как

$\mathbf {T} _{n}=\left[{\begin{array}{cccc}\mathbf {f} _{\mathbf {x} }(\mathbf {y} _{n})&(\mathbf {I} \otimes \mathbf {f} ^{\intercal }(\mathbf {y} _{n}))\mathbf {f} _{\mathbf {xx} }(\mathbf {y} _{n})\mathbf {f} (\mathbf {y} _{n})&\mathbf {0} &\mathbf {f} (\mathbf {y} _{n})\\0&0&0&0\\0&0&0&1\\0&0&0&0\end{array}}\right]\in \mathbb {R} ^{(d+3)\times (d+3)},$

$\mathbf {L} _{1}=\left[{\begin{array}{ll}\mathbf {I} &\mathbf {0} _{d\times 3}\end{array}}\right]\quad and\quad \mathbf {r} _{1}^{\intercal }=\left[{\begin{array}{ll}\mathbf {0} _{1\times (d+2)}&1\end{array}}\right]$. Здесь, $\mathbf {f} _{\mathbf {xx} }$ обозначает вторую производную f по x и p + q > 2 . Для больших систем ОДУ

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {T} _{n},\mathbf {r} )\mathbf {,} \qquad {\text{ with }}\qquad m_{n}>3.$

Схемы ЛЛ-РК Приказ-4

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {u} _{4}+{\frac {h_{n}}{6}}(2\mathbf {k} _{2}+2\mathbf {k} _{3}+\mathbf {k} _{4}),\quad$ $\qquad \qquad (4.8)$

где

$\mathbf {u} _{j}=\mathbf {L} (\mathbf {P} _{p,q}(2^{-\kappa _{j}}\mathbf {M} _{n}c_{j}h_{n}))^{2^{\kappa _{j}}}\mathbf {r}$

и

$\mathbf {k} _{j}=\mathbf {f} \left(t_{n}+c_{j}h_{n},\mathbf {y} _{n}+\mathbf {u} _{j}+c_{j}h_{n}\mathbf {k} _{j-1}\right)-\mathbf {f} \left(t_{n},\mathbf {y} _{n}\right)-\mathbf {f} _{\mathbf {x} }\left(t_{n},\mathbf {y} _{n}\right)\mathbf {u} _{j}\ -\mathbf {f} _{t}\left(t_{n},\mathbf {y} _{n}\right)c_{j}h_{n},$

с $\mathbf {k} _{1}\equiv \mathbf {0} ,c=\left[{\begin{array}{cccc}0&{\frac {1}{2}}&{\frac {1}{2}}&1\end{array}}\right],$ и р + q > 3 . Для больших систем ОДУ вектор $\mathbf {u} _{j}$ в приведенной выше схеме заменяется на $\mathbf {u} _{j}=\mathbf {L\mathbf {k} } _{m_{j},k_{j}}^{p,q}(c_{j}h_{n},\mathbf {M} _{n},\mathbf {r} )$ с $m_{j}>4.$

Кутты Дормана и

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {u} _{s}+h_{n}\sum _{j=1}^{s}b_{j}\mathbf {k} _{j}\qquad {\text{ and }}\qquad {\widehat {\mathbf {y} }}_{n+1}=\mathbf {y} _{n}+\mathbf {u} _{s}+h_{n}\sum _{j=1}^{s}{\widehat {b}}_{j}\mathbf {k} _{j},\quad$$\qquad \qquad (4.9)$

где s = 7 – количество ступеней,

$\mathbf {k} _{j}=\mathbf {f(} t_{n}+c_{j}h_{n},\mathbf {y} _{n}+\mathbf {u} _{j}+h_{n}\sum _{i=1}^{s-1}a_{j,i}\mathbf {k} _{i})-\mathbf {f} \left(t_{n},\mathbf {y} _{n}\right)-\mathbf {f} _{\mathbf {x} }\left(t_{n},\mathbf {y} _{n}\right)\mathbf {u} _{j}\ -\mathbf {f} _{t}\left(t_{n},\mathbf {y} _{n}\right)c_{j}h_{n},$

с $\mathbf {k} _{1}\equiv \mathbf {0}$, и $a_{j,i},b_{j},{\widehat {b}}_{j}\quad and\quad c_{j}$ – коэффициенты Рунге–Кутты Дормана и Принса , p + q > 4. Вектор $\mathbf {u} _{j}$ в приведенной выше схеме вычисляется с помощью аппроксимации Паде или Крайлора – Паде для малых или больших систем ОДУ соответственно.

Стабильность и динамика

Рис. 1 Фазовый портрет (пунктирная линия) и приближенный фазовый портрет (сплошная линия) нелинейного ОДУ (4.10)-(4.11), рассчитанного по схеме LL 2-го порядка (4.2), классической схеме Рюгена-Кутты RK 4-го порядка 4 и схемы LLRK 4 порядка 4 (4.8) с размером шага h=1/2 и p=q=6.

По построению дискретизации LL и HOLL наследуют стабильность и динамику линейных ОДУ, но это не относится к схемам LL в целом. С $p\leq q\leq p+2$, ЛЛ-схемы (4.6)-(4.9) A -стабильны . При q = p + 1 или q = p + 2 ЛЛ-схемы (4.6)–(4.9) также являются L -стабильными . Для линейных ОДУ схемы ЛЛ (4.6)–(4.9) сходятся с порядком p + q . Кроме того, при p = q = 6 и $m_{n}$ = d , все вышеописанные схемы LL поддаются «точному вычислению» (с точностью до арифметики с плавающей запятой ) линейных ОДУ на современных персональных компьютерах. Сюда входят жесткие и сильно колебательные линейные уравнения. Более того, ЛЛ-схемы (4.6)-(4.9) регулярны для линейных ОДУ и наследуют симплектическую структуру гамильтоновых гармонических осцилляторов . Эти схемы LL также сохраняют линеаризацию и демонстрируют лучшее воспроизведение устойчивых и неустойчивых многообразий вокруг гиперболических точек равновесия и периодических орбит , чем другие численные схемы с тем же размером шага. Например, на рис. 1 показан фазовый портрет ОДУ

${\begin{aligned}&{\frac {dx_{1}}{dt}}=-2x_{1}+x_{2}+1-\mu f(x_{1},\lambda )\qquad \qquad (4.10)\\[6pt]&{\frac {dx_{2}}{dt}}=x_{1}-2x_{2}+1-\mu f(x_{2},\lambda )\qquad \qquad \quad (4.11)\end{aligned}}$

с $f(u,\lambda )=u(1+u+\lambda u^{2})^{-1}$, $\mu =15$ и $\lambda =57$и его аппроксимация различными схемами. Эта система имеет две устойчивые стационарные точки и одну неустойчивую стационарную точку в области $0\leq x_{1},x_{2}\leq 1$.

Методы LL для DDE

Рассмотрим d -мерное дифференциальное уравнение с запаздыванием (DDE).

${\frac {d\mathbf {x} (t)}{dt}}=\mathbf {f} (t,\mathbf {x} (t),\mathbf {x} _{t}(-\tau _{1}),\ldots ,\mathbf {x} _{t}(-\tau _{m})),\qquad t\in [t_{0},T],\qquad \qquad (5.1)$

с m постоянными задержками $\tau _{i}>0$ и начальное состояние $\mathbf {x} _{t_{0}}(s)=\mathbf {\varphi } (s)$ для всех $s\in [-\tau ,0],$ где f — дифференцируемая функция, $\mathbf {x} _{t}:[-\tau ,0]\longrightarrow \mathbb {R} ^{d}$ — функция сегмента, определенная как

$\mathbf {x} _{t}(s):=\mathbf {x} (t+s),{\text{ }}s\in [-\tau ,0],$

для всех $t\in [t_{0},T],\mathbf {\varphi } :[-\tau ,0]\longrightarrow \mathbb {R} ^{d}$ — заданная функция, и $\tau =\max \left\{\tau _{1},\ldots ,\tau _{m}\right\}.$

Локальная

Для временной дискретизации $(t)_{h}$ , Локальная линейная дискретизация ДДУ (5.1) в каждой точке $t_{n+1}\in (t)_{h}$ определяется рекурсивным выражением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\Phi (t_{n},\mathbf {z} _{n},h_{n};{\widetilde {\mathbf {z} }}_{t_{n}}^{1},\ldots ,{\widetilde {\mathbf {z} }}_{t_{n}}^{m}),\qquad \qquad (5.2)$

где

$\Phi (t_{n},\mathbf {z} _{n},h_{n};{\widetilde {\mathbf {z} }}_{t_{n}}^{1},\ldots ,{\widetilde {\mathbf {z} }}_{t_{n}}^{m})=\int \limits _{0}^{h_{n}}e^{\mathbf {A} _{n}(h_{n}-u)}\left[\sum \limits _{i=1}^{m}\mathbf {B} _{n}^{i}({\widetilde {\mathbf {z} }}_{t_{n}}^{i}(u-\tau _{i})-{\widetilde {\mathbf {z} }}_{t_{n}}^{i}(-\tau _{i}))+\mathbf {d} _{n}\right]\,du+\int \limits _{0}^{h_{n}}\int \limits _{0}^{u}e^{\mathbf {A} _{n}(h_{n}-u)}\mathbf {c} _{n}\,dr\,du$

${\widetilde {\mathbf {z} }}_{t_{n}}^{i}:\left[-\tau _{i},0\right]\longrightarrow \mathbb {R} ^{d}$ — функция сегмента, определенная как

${\widetilde {\mathbf {z} }}_{t_{n}}^{i}(s):={\widetilde {\mathbf {z} }}^{i}(t_{n}+s),{\text{ }}s\in [-\tau _{i},0],$

и ${\widetilde {\mathbf {z} }}^{i}:\left[t_{n}-\tau _{i},t_{n}\right]\longrightarrow \mathbb {R} ^{d}$является подходящим приближением к $\mathbf {x} (t)$ для всех $t\in \lbrack t_{n}-\tau _{i},t_{n}]$ такой, что ${\widetilde {\mathbf {z} }}^{i}(t_{n})=\mathbf {z} _{n}.$ Здесь,

$\mathbf {A} _{n}=\mathbf {f} _{x}(t_{n},\mathbf {z} _{n},{\widetilde {\mathbf {z} }}_{t_{n}}^{1}(-\tau _{1}),\ldots ,{\widetilde {\mathbf {z} }}_{t_{n}}^{m}(-\tau _{d})),{\text{ }}\mathbf {B} _{n}^{i}=\mathbf {f} _{x_{t}(-\tau _{i})}(t_{n},\mathbf {z} _{n},{\widetilde {\mathbf {z} }}_{t_{n}}^{1}(-\tau _{1}),\ldots ,{\widetilde {\mathbf {z} }}_{t_{n}}^{m}(-\tau _{d}))$

являются постоянными матрицами и

$\mathbf {c} _{n}=\mathbf {f} _{t}(t_{n},\mathbf {z} _{n},{\widetilde {\mathbf {z} }}_{t_{n}}^{1}(-\tau _{1}),\ldots ,{\widetilde {\mathbf {z} }}_{t_{n}}^{m}(-\tau _{d})){\text{ and }}\mathbf {d} _{n}=\mathbf {f(} t_{n},\mathbf {z} _{n},{\widetilde {\mathbf {z} }}_{t_{n}}^{1}(-\tau _{1}),\ldots ,{\widetilde {\mathbf {z} }}_{t_{n}}^{m}(-\tau _{d}))$

являются постоянными векторами. $\mathbf {f} _{t},\mathbf {f} _{x}\quad and\quad \mathbf {f} _{x_{t}(-\tau _{i})}$ обозначаем соответственно частные производные f по переменным t и x , а $\mathbf {x} _{t}(-\tau _{i})$. Локальная линейная дискретизация (5.2) сходится к решению (5.1) с порядком $\alpha =\min\{2,r\},$ если ${\widetilde {\mathbf {z} }}_{t_{n}}^{i}$ приближает $\mathbf {z} _{t_{n}}^{i}$ с заказом $r\quad (i.e.,\left\vert \mathbf {z} _{t_{n}}^{i}\mathbf {(} u-\tau _{i}\mathbf {)} -{\widetilde {\mathbf {z} }}_{t_{n}}^{i}\mathbf {(} u-\tau _{i}\mathbf {)} \right\vert \propto h_{n}^{r}$ для всех $u\in \lbrack 0,h_{n}])$.

Схемы локальной линеаризации

Рис. 2. Примерные пути Марчука и др. (1991) противовирусная иммунная модель, описываемая жесткой системой десятимерных нелинейных DDE с пятью временными задержками: вверху — Рунге–Кутты (2,3) непрерывная схема ; внизу — схема ЛЛ (5.3). Размер шага h = 0,01 фиксированный, а p = q = 6.

В зависимости от приближений ${\widetilde {\mathbf {z} }}_{t_{n}}^{i}$ и об алгоритме вычисления $\mathbf {\phi }$ могут быть определены различные схемы локальной линеаризации. Каждая численная реализация $\mathbf {y} _{n}$ локальной линейной дискретизации $\mathbf {z} _{n}$ обычно называется схемой локальной линеаризации .

Схемы LL с полиномом 2-го порядка

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r,} \quad$ $\qquad (5.3)$

где матрицы $\mathbf {M} _{n},\mathbf {L}$ и $\mathbf {r}$ определяются как

$\mathbf {M} _{n}={\begin{bmatrix}\mathbf {A} _{n}&\mathbf {c} _{n}+\sum \limits _{i=1}^{m}\mathbf {B} _{n}^{i}\mathbf {\alpha } _{n}^{i}&\mathbf {d} _{n}\\0&0&1\\0&0&0\end{bmatrix}}\in \mathbb {R} ^{(d+2)\times (d+2)},$

$\mathbf {L} =\left[{\begin{array}{ll}\mathbf {I} &\mathbf {0} _{d\times 2}\end{array}}\right]$ и $\mathbf {r} ^{\intercal }=\left[{\begin{array}{ll}\mathbf {0} _{1\times (d+1)}&1\end{array}}\right],h_{n}\leq \tau$, и $p+q>1$. Здесь матрицы $\mathbf {A} _{n}$, $\mathbf {B} _{n}^{i}$, $\mathbf {c} _{n}$ и $\mathbf {d} _{n}$ определяются как в (5.2), но заменяя $\mathbf {z}$ к $\mathbf {y}$ и $\mathbf {\alpha } _{n}^{i}=(\mathbf {y} (t_{n+1}-\tau _{i})-\mathbf {y} (t_{n}-\tau _{i}))/h_{n},$ где

$\mathbf {y} \left(t\right)=\mathbf {y} _{n_{t}}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n_{t}}(t-t_{n_{t}})))^{2^{k_{n}}}\mathbf {r} ,$

с $n_{t}=\max\{n=0,1,2,...,:t_{n}\leq t{\text{ and }}t_{n}\in \left(t\right)_{h}\}$, — локальная линейная аппроксимация решения (5.1), определенная с помощью схемы ЛЛ (5.3) для всех $t\in \lbrack t_{0},t_{n}]$ и по $\mathbf {y} \left(t\right)=\mathbf {\varphi } \left(t\right)$ для $t\in \left[t_{0}-\tau ,t_{0}\right]$. Для больших систем DDE

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {M} _{n},\mathbf {r} )\quad and\quad \mathbf {y} \left(t\right)=\mathbf {y} _{n_{t}}+\mathbf {L\mathbf {k} } _{m_{n_{t}},k_{n_{t}}}^{p,q}(t-t_{n_{t}},\mathbf {M} _{n_{t}},\mathbf {r} ),$

с $p+q>1$ и $m_{n}>2$. Рис. 2 иллюстрирует устойчивость схемы ЛЛ (5.3) и явной схемы аналогичного порядка при интегрировании жесткой системы ДДУ.

Методы LL для RDE

Рассмотрим d -мерное случайное дифференциальное уравнение (СДУ)

${\frac {d\mathbf {x} \left(t\right)}{dt}}=\mathbf {f} (\mathbf {x} (t),\mathbf {\xi } (t)),\quad t\in \left[t_{0},T\right],\qquad \qquad \qquad (6.1)$

с начальным состоянием $\mathbf {x} (t_{0})=\mathbf {x} _{0},$ где $\mathbf {\xi }$ — k -мерный сепарабельный конечный непрерывный случайный процесс , а f — дифференцируемая функция. Предположим, что реализация (путь) $\mathbf {\xi }$ дано.

Локальная

Для временной дискретизации $\left(t\right)_{h}$, Локальная линейная дискретизация РДУ (6.1) в каждой точке $t_{n+1}\in \left(t\right)_{h}$ определяется рекурсивным выражением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } (t_{n},\mathbf {z} _{n};h_{n}),\qquad {\text{ with }}\qquad \mathbf {z} _{0}=\mathbf {x} _{0},$

где

$\mathbf {\phi } (t_{n},\mathbf {z} _{n};h_{n})=\int \limits _{0}^{h_{n}}e^{\mathbf {f} _{\mathbf {x} }(\mathbf {z} _{n},\mathbf {\xi } (t_{n}))(h_{n}-u)}(\mathbf {f(z} _{n},\mathbf {\xi } (t_{n}))+\mathbf {f} _{\mathbf {\xi } }(\mathbf {z} _{n},\mathbf {\xi } (t_{n}))({\widetilde {\mathbf {\xi } }}(t_{n}+u)-{\widetilde {\mathbf {\xi } }}(t_{n})))\,du$

и ${\widetilde {\mathbf {\xi } }}$ является приближением к процессу $\mathbf {\xi }$ для всех $t\in \left[t_{0},T\right].$ Здесь, $\mathbf {f} _{x}$ и $\mathbf {f} _{\xi }$ обозначаем частные производные $\mathbf {f}$ относительно $\mathbf {x}$ и $\xi$, соответственно.

Схемы локальной линеаризации

Рис. 3. Фазовый портрет траекторий схем Эйлера и ЛЛ при интегрировании нелинейного СДУ (6.2)–(6.3) с размером шага h = 1/32 и p = q = 6.

В зависимости от приближений ${\widetilde {\mathbf {\xi } }}$ к процессу $\mathbf {\xi }$ и алгоритма вычисления $\mathbf {\phi }$, могут быть определены различные схемы локальной линеаризации. Каждая численная реализация $\mathbf {y} _{n}$ локальной линейной дискретизации $\mathbf {z} _{n}$ обычно называется схемой локальной линеаризации.

Схемы LL

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r,} \quad$

где матрицы $\mathbf {M} _{n},\quad \mathbf {L} \quad and\quad \mathbf {r}$ определяются как

$\mathbf {M} _{n}=\left[{\begin{array}{ccc}\mathbf {f} _{\mathbf {x} }\left(\mathbf {y} _{n},\mathbf {\xi } (t_{n})\right)&\mathbf {f} _{\mathbf {\xi } }(\mathbf {y} _{n},\mathbf {\xi } (t_{n})(\mathbf {\xi } (t_{n+1})-\mathbf {\xi } (t_{n}))/h_{n}&\mathbf {f} \left(\mathbf {y} _{n},\mathbf {\xi } (t_{n})\right)\\0&0&1\\0&0&0\end{array}}\right]$

$\mathbf {L} =\left[{\begin{array}{ll}\mathbf {I} &\mathbf {0} _{d\times 2}\end{array}}\right]$, $\mathbf {r} ^{\intercal }=\left[{\begin{array}{ll}\mathbf {0} _{1\times (d+1)}&1\end{array}}\right]$, и p+q>1 . Для больших систем RDE

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {M} _{n},\mathbf {r} ),\quad p+q>1\quad and\quad m_{n}>2.$

Скорость сходимости обеих схем равна $min\{2,2\gamma \}$, где находится $\gamma$ показатель условия Холдера $\mathbf {\xi }$.

На рис. 3 представлен фазовый портрет РДЭ.

${\frac {dx_{1}}{dt}}=-x_{2}+\left(1-x_{1}^{2}-x_{2}^{2}\right)x_{1}\sin(w^{H}(t))^{2},\quad \qquad x_{1}(0)=0.8\qquad (6.2)$

${\frac {dx_{2}}{dt}}=x_{1}+(1-x_{1}^{2}-x_{2}^{2})x_{2}\sin(w^{H}(t))^{2},\qquad \qquad x_{2}(0)=0.1,\qquad (6.3)$

и ее аппроксимация двумя численными схемами, где $w^{H}$ обозначает дробный броуновский процесс с показателем Херста H=0,45 .

методы LL

Рассмотрим d -мерное стохастическое дифференциальное уравнение (СДУ)

$d\mathbf {x} (t)=\mathbf {f} (t,\mathbf {x} (t))dt+\sum \limits _{i=1}^{m}\mathbf {g} _{i}(t)d\mathbf {w} ^{i}(t),\quad t\in \left[t_{0},T\right],\qquad \qquad \qquad (7.1)$

с начальным состоянием $\mathbf {x} (t_{0})=\mathbf {x} _{0}$, где коэффициент дрейфа $\mathbf {f}$ и коэффициент диффузии $\mathbf {g} _{i}$ являются дифференцируемыми функциями, а $\mathbf {w=(\mathbf {w} } ^{1},\ldots ,\mathbf {w} ^{m}\mathbf {)}$ является m -мерным стандартным винеровским процессом .

Локальная

Для временной дискретизации $\left(t\right)_{h}$ , заказ- $\mathbb {\gamma }$ (=1,1.5) Сильная локальная линейная дискретизация решения СДУ (7.1) определяется рекуррентным соотношением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } _{\mathbb {\gamma } }(t_{n},\mathbf {z} _{n};h_{n})+\mathbf {\xi } (t_{n},\mathbf {z} _{n};h_{n}),\quad with\quad \mathbf {z} _{0}=\mathbf {x} _{0},$

где

$\mathbf {\phi } _{\mathbb {\gamma } }(t_{n},\mathbf {z} _{n};\delta )=\int _{0}^{\delta }e^{\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n})(\delta -u)}(\mathbf {f(} t_{n},\mathbf {z} _{n})+\mathbf {a} ^{\mathbb {\gamma } }(t_{n},\mathbf {z} _{n})u)du$

и

$\mathbf {\xi } \left(t_{n},\mathbf {z} _{n};\delta \right)=\sum \limits _{i=1}^{m}\int \nolimits _{t_{n}}^{t_{n}+\delta }e^{\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})(t_{n}+\delta -u)}\mathbf {g} _{i}(u)\,d\mathbf {w} ^{i}(u).$

Здесь,

$\mathbf {a} ^{\mathbb {\gamma } }(t_{n},\mathbf {z} _{n})=\left\{{\begin{array}{cl}\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})&{\text{for }}\qquad \mathbb {\gamma } =1\\\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})+{\frac {1}{2}}\sum \limits _{j=1}^{m}(\mathbf {I} \otimes \mathbf {g} _{j}^{\intercal }(t_{n}))\mathbf {f} _{\mathbf {xx} }(t_{n},\mathbf {z} _{n})\mathbf {g} _{j}(t_{n})&{\text{for }}\quad \mathbb {\gamma } =1.5,\end{array}}\right.$

$\mathbf {f} _{\mathbf {x} },\mathbf {f} _{t}$ обозначаем частные производные $\mathbf {f}$ относительно переменных $\mathbf {x}$ и t соответственно и $\mathbf {f} _{\mathbf {xx} }$ матрица Гессе $\mathbf {f}$ относительно $\mathbf {x}$. Сильная локальная линейная дискретизация $\mathbf {z} _{n+1}$ сходится с порядком $\mathbb {\gamma }$ (= 1, 1,5) к решению (7.1).

высокого

После локальной линеаризации дрейфового члена (7.1) при $(t_{n},\mathbf {z} _{n})$, уравнение для невязки $\mathbf {r}$ дается

$d\mathbf {r} (t)=\mathbf {q} _{\gamma }(t_{n},\mathbf {z} _{n};t\mathbf {,\mathbf {r} } (t))\,dt+\sum \limits _{i=1}^{m}\mathbf {g} _{i}(t)\,d\mathbf {w} ^{i}(t)\mathbf {,} \qquad \mathbf {r} (t_{n})=\mathbf {0}$

для всех $t\in \lbrack t_{n},t_{n+1}]$, где

$\mathbf {q} _{\gamma }(t_{n},\mathbf {z} _{n};s\mathbf {,\xi } )=\mathbf {f} (s,\mathbf {z} _{n}+\mathbf {\phi } _{\gamma }(t_{n},\mathbf {z} _{n};s-t_{n})+\mathbf {\xi } )-\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})\mathbf {\phi } _{\gamma }(t_{n},\mathbf {z} _{n};s-t_{n})-\mathbf {a} ^{\gamma }(t_{n},\mathbf {z} _{n})(s-t_{n})-\mathbf {f} (t_{n},\mathbf {z} _{n}).$

Локальная линейная дискретизация высокого порядка СДУ (7.1) в каждой точке $t_{n+1}\in (t)_{h}$ затем определяется рекурсивным выражением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } _{\gamma }(t_{n},\mathbf {z} _{n};h_{n})+{\widetilde {\mathbf {r} }}(t_{n},\mathbf {z} _{n};h_{n}),\qquad {\text{ with }}\qquad \mathbf {z} _{0}=\mathbf {x} _{0},$

где ${\widetilde {\mathbf {r} }}$ является сильным приближением к остатку $\mathbf {r}$ порядка $\alpha$ выше 1,5 . Сильная дискретизация ХОЛЛА. $\mathbf {z} _{n+1}$ сходится с порядком $\alpha$ к решению (7.1).

Схемы локальной линеаризации

В зависимости от способа расчета $\mathbf {\phi } _{\mathbb {\gamma } }$ , $\mathbf {\xi }$ и ${\widetilde {\mathbf {r} }}$ могут быть получены различные численные схемы. Каждая численная реализация $\mathbf {y} _{n}$ сильной локальной линейной дискретизации $\mathbf {z} _{n}$ любого порядка обычно называется схемой сильной локальной линеаризации (SLL) .

Схемы SLL заказа 1

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r+} \sum \limits _{i=1}^{m}\mathbf {g} _{i}(t_{n})\Delta \mathbf {w} _{n}^{i},\quad$ $\qquad \qquad (7.2)$

где матрицы $\mathbf {M} _{n}$, $\mathbf {L}$ и $\mathbf {r}$ определяются как в (4.6), $\Delta \mathbf {w} _{n}^{i}$ представляет собой iid нулевое среднее гауссовской случайной величины с дисперсией $h_{n}$, и p + q > 1. Для больших систем СДУ в приведенной выше схеме $(\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r}$ заменяется на $\mathbf {\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {M} _{n},\mathbf {r} )$.

Схемы SLL заказа 1.5

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r} +\sum \limits _{i=1}^{m}\left(\mathbf {g} _{i}(t_{n})\Delta \mathbf {w} _{n}^{i}\mathbf {f} _{\mathbf {x} }(t_{n},{\widetilde {\mathbf {y} }}_{n})\mathbf {g} _{i}(t_{n})\Delta \mathbf {z} _{n}^{i}+{\frac {d\mathbf {g} _{i}(t_{n})}{dt}}(\Delta \mathbf {w} _{n}^{i}h_{n}-\Delta \mathbf {z} _{n}^{i})\right),\qquad \qquad (7.3)$

где матрицы $\mathbf {M} _{n}$, $\mathbf {L}$ и $\mathbf {r}$ определяются как

$\mathbf {M} _{n}={\begin{bmatrix}\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n})&\mathbf {f} _{t}(t_{n},\mathbf {y} _{n})+{\frac {1}{2}}\sum \limits _{j=1}^{m}\left(\mathbf {I} \otimes \mathbf {g} _{j}^{\intercal }(t_{n})\right)\mathbf {f} _{\mathbf {xx} }(t_{n},\mathbf {y} _{n})\mathbf {g} _{j}(t_{n})&\mathbf {f} (t_{n},\mathbf {y} _{n})\\0&0&1\\0&0&0\end{bmatrix}}\in \mathbb {R} ^{(d+2)\times (d+2)},$

$\mathbf {L} =\left[{\begin{array}{ll}\mathbf {I} &\mathbf {0} _{d\times 2}\end{array}}\right],\mathbf {r} ^{\intercal }=\left[{\begin{array}{ll}\mathbf {0} _{1\times (d+1)}&1\end{array}}\right]$, $\Delta \mathbf {z} _{n}^{i}$ является iid нулевым средним гауссовой случайной величиной с дисперсией $E\left((\Delta \mathbf {z} _{n}^{i})^{2}\right)={\frac {1}{3}}h_{n}^{3}$ и ковариация $E(\Delta \mathbf {w} _{n}^{i}\Delta \mathbf {z} _{n}^{i})={\frac {1}{2}}h_{n}^{2}$ и p+q>1 . Для больших систем СДУ в приведенной выше схеме $(\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r}$ заменяется на $\mathbf {\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {M} _{n},\mathbf {r} )$.

Заказать 2 схемы SLL-Тейлора

$\mathbf {y} _{t_{n+1}}=\mathbf {y} _{n}+\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r} +\sum \limits _{j=1}^{m}\mathbf {g} _{j}\left(t_{n}\right)\Delta \mathbf {w} _{n}^{j}+\sum \limits _{j=1}^{m}\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n})\mathbf {g} _{j}\left(t_{n}\right){\widetilde {J}}_{\left(j,0\right)}+\sum \limits _{j=1}^{m}{\frac {d\mathbf {g} _{_{j}}}{dt}}\left(t_{n}\right){\widetilde {J}}_{\left(0,j\right)}$

$\qquad \qquad +\sum \limits _{j_{1},j_{2}=1}^{m}\left(\mathbf {I} \otimes \mathbf {g} _{j_{2}}^{\intercal }\left(t_{n}\right)\right)\mathbf {f} _{\mathbf {xx} }(t_{n},\mathbf {y} _{n})\mathbf {g} _{j_{1}}\left(t_{n}\right){\widetilde {J}}_{\left(j_{1},j_{2},0\right),}\qquad \qquad (7.4)$

где $\mathbf {M} _{n}$, $\mathbf {L}$, $\mathbf {r}$ и $\Delta \mathbf {w} _{n}^{i}$ определяются как в схемах SLL первого порядка, и ${\widetilde {J}}_{\alpha }$ является приближением второго порядка кратного интеграла Стратонова $J_{\alpha }$.

Заказать 2 схемы СЛЛ-РК

Рис. 4, вверху : Эволюция областей на фазовой плоскости гармонического осциллятора (7.6) при ε=0 и ω=σ=1. Изображения исходного единичного круга (зеленый) получены в три момента времени T по точному решению (черный) и по схемам SLL1 (синий) и Неявный Эйлер (красный) с h=0,05 . Внизу : ожидаемое значение энергии (сплошная линия) вдоль решения нелинейного осциллятора (7.6) с ε=1 и ω=100 и его аппроксимация (кружки), вычисленная методом Монте-Карло с помощью 10 000 симуляций схемы SLL1 с h. =1/2 и p=q=6 .

Для СДУ с одним винеровским шумом (m=1 )

$\mathbf {y} _{t_{n+1}}=\mathbf {y} _{n}+{\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};h_{n})+{\frac {h_{n}}{2}}\left(\mathbf {k} _{1}+\mathbf {k} _{2}\right)+\mathbf {g} \left(t_{n}\right)\Delta w_{n}+{\frac {\left(\mathbf {g} \left(t_{n+1}\right)-\mathbf {g} \left(t_{n}\right)\right)}{h_{n}}}J_{\left(0,1\right)}\quad (7.5)$

$\quad \quad \quad$

где

$\mathbf {k} _{1}=\mathbf {f} (t_{n}+{\frac {h_{n}}{2}},\mathbf {y} _{n}+{\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};{\frac {h_{n}}{2}})+\gamma _{+})-\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n}){\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};{\frac {h_{n}}{2}})-\mathbf {f} \left(t_{n},\mathbf {y} _{n}\right)-\mathbf {f} _{t}\left(t_{n},\mathbf {y} _{n}\right){\frac {h_{n}}{2}},$

$\mathbf {k} _{2}=\mathbf {f} (t_{n}+{\frac {h_{n}}{2}},\mathbf {y} _{n}+{\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};{\frac {h_{n}}{2}})+\gamma _{-})-\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n}){\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};{\frac {h_{n}}{2}})-\mathbf {f} \left(t_{n},\mathbf {y} _{n}\right)-\mathbf {f} _{t}\left(t_{n},\mathbf {y} _{n}\right){\frac {h_{n}}{2}},$

с $\gamma _{\pm }={\frac {1}{h_{n}}}\mathbf {g} \left(t_{n}\right){\Bigl (}{\widetilde {J}}_{\left(1,0\right)}\pm {\sqrt {2{\widetilde {J}}_{\left(1,1,0\right)}h_{n}-{\widetilde {J}}_{\left(1,0\right)}^{2}}}{\Bigr )}$.

Здесь, ${\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};h_{n})=\mathbf {L} (\mathbf {P} _{p,q}(2^{-k_{n}}\mathbf {M} _{n}h_{n}))^{2^{k_{n}}}\mathbf {r}$ для низкоразмерных СДУ и ${\widetilde {\mathbf {\phi } }}(t_{n},\mathbf {y} _{n};h_{n})=\mathbf {L\mathbf {k} } _{m_{n},k_{n}}^{p,q}(h_{n},\mathbf {M} _{n},\mathbf {r} )$ для больших систем СДУ, где $\mathbf {M} _{n}$, $\mathbf {L}$, $\mathbf {r}$, $\Delta \mathbf {w} _{n}^{i}$ и ${\widetilde {J}}_{\alpha }$ определяются как в схемах SLL-Тейлора порядка 2 , p+q>1 и $m_{n}>2$.

Стабильность и динамика

По построению сильные дискретизации LL и HOLL наследуют устойчивость и динамику линейных СДУ, но это не относится к сильным схемам LL в целом. схемы ЛЛ (7.2)-(7.5) с $p\leq q\leq p+2$ являются A -стабильными, включая жесткие и высококолебательные линейные уравнения. Более того, для линейных СДУ со случайными аттракторами эти схемы также имеют случайный аттрактор, сходящийся по вероятности к точному при уменьшении размера шага и сохраняющий эргодичность этих уравнений для любого размера шага. Эти схемы также воспроизводят существенные динамические свойства простых и связанных гармонических осцилляторов, такие как линейный рост энергии вдоль путей, колебательное поведение около 0, симплектическая структура гамильтоновых осцилляторов и среднее значение путей. Для нелинейных СДУ с малым шумом (т. е. (7.1) с $\mathbf {g} _{i}(t)\approx 0$), пути этих схем SLL по сути представляют собой неслучайные пути схемы LL (4.6) для ОДУ плюс небольшое возмущение, связанное с малым шумом. В этой ситуации для путей схемы SLL становятся актуальными динамические свойства этой детерминированной схемы, такие как сохранение линеаризации и сохранение точной динамики решения вокруг гиперболических точек равновесия и периодических орбит. Например, на рис. 4 показана эволюция доменов на фазовой плоскости и энергия стохастического осциллятора.

${\begin{array}{ll}dx(t)=y(t)dt,&x_{1}(0)=0.01\\dy(t)=-(\omega ^{2}x(t)+\epsilon x^{4}(t))dt+\sigma dw_{t},&x_{1}(0)=0.1,\end{array}}\qquad \qquad (7.6)$

и их аппроксимации двумя численными схемами.

методы LL

Рассмотрим d -мерное стохастическое дифференциальное уравнение

$d\mathbf {x} (t)=\mathbf {f} (t,\mathbf {x} (t))dt+\sum \limits _{i=1}^{m}\mathbf {g} _{i}(t)d\mathbf {w} ^{i}(t),\qquad t\in \left[t_{0},T\right],\qquad \qquad (8.1)$

с начальным состоянием $\mathbf {x} (t_{0})=\mathbf {x} _{0}$, где коэффициент дрейфа $\mathbf {f}$ и коэффициент диффузии $\mathbf {g} _{i}$ являются дифференцируемыми функциями, а $\mathbf {w=(\mathbf {w} } ^{1},\ldots ,\mathbf {w} ^{m}\mathbf {)}$ представляет собой m -мерный стандартный винеровский процесс.

Локальная

Для временной дискретизации $\left(t\right)_{h}$, заказ- $\mathbb {\beta }$ $(=1,2)$ Слабая локальная линейная дискретизация решения СДУ (8.1) определяется рекурсивным соотношением

$\mathbf {z} _{n+1}=\mathbf {z} _{n}+\mathbf {\phi } _{\mathbb {\beta } }(t_{n},\mathbf {z} _{n};h_{n})+\mathbf {\eta } (t_{n},\mathbf {z} _{n};h_{n}),\quad with\quad \mathbf {z} _{0}=\mathbf {x} _{0},$

где

$\mathbf {\phi } _{\mathbb {\beta } }(t_{n},\mathbf {z} _{n};\delta )=\int _{0}^{\delta }e^{\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})(\delta -u)}(\mathbf {f(} t_{n},\mathbf {z} _{n})+\mathbf {b} ^{\mathbb {\beta } }(t_{n},\mathbf {z} _{n})u)du$

с

$\mathbf {b} ^{\mathbb {\beta } }(t_{n},\mathbf {z} _{n})={\begin{cases}\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})&{\text{for }}\mathbb {\beta } =1\\\mathbf {f} _{t}(t_{n},\mathbf {z} _{n})+{\frac {1}{2}}\sum \limits _{j=1}^{m}\left(\mathbf {I} \otimes \mathbf {g} _{j}^{\intercal }\left(t_{n}\right)\right)\mathbf {f} _{\mathbf {xx} }(t_{n},\mathbf {z} _{n})\mathbf {g} _{j}\left(t_{n}\right)&{\text{for }}\mathbb {\beta } =2,\end{cases}}$

и $\mathbf {\eta } (t_{n},\mathbf {z} _{n};\delta )$ представляет собой стохастический процесс с нулевым средним с матрицей дисперсии

$\mathbf {\Sigma } (t_{n},\mathbf {z} _{n};\delta )=\int \limits _{0}^{\delta }e^{\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {z} _{n})(\delta -s)}\mathbf {G} (t_{n}+s)\mathbf {G} ^{\intercal }(t_{n}+s)e^{\mathbf {f} _{\mathbf {x} }^{\intercal }(t_{n},\mathbf {z} _{n})(\delta -s)}ds.$

Здесь, $\mathbf {f} _{\mathbf {x} }$, $\mathbf {f} _{t}$ обозначаем частные производные $\mathbf {f}$ относительно переменных $\mathbf {x}$ и t соответственно $\mathbf {f} _{\mathbf {xx} }$ матрица Гессе $\mathbf {f}$ относительно $\mathbf {x}$, и $\mathbf {G} (t)=[\mathbf {g} _{1}(t),\ldots ,\mathbf {g} _{m}(t)]$. Слабая локальная линейная дискретизация $\mathbf {z} _{n+1}$ сходится с порядком $\mathbb {\beta }$ (=1,2) к решению (8.1).

Схемы локальной линеаризации

В зависимости от способа расчета $\mathbf {\phi } _{\mathbb {\beta } }$ и $\mathbf {\Sigma }$ могут быть получены различные численные схемы. Каждая численная реализация $\mathbf {y} _{n}$ слабой локальной линейной дискретизации $\mathbf {z} _{n}$ обычно называется схемой слабой локальной линеаризации (WLL) .

Схема WLL заказа 1

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {B} _{14}+(\mathbf {B} _{12}\mathbf {B} _{11}^{\intercal })^{1/2}\mathbf {\xi } _{n}$

где для СДУ с автономными коэффициентами диффузии $\mathbf {B} _{11}$, $\mathbf {B} _{12}$ и $\mathbf {B} _{14}$ являются подматрицами, определяемыми разделенной матрицей $\mathbf {B} =\mathbf {P} _{p,q}(2^{-k_{n}}{\mathcal {M}}_{n}h_{n}))^{2^{k_{n}}}$, с

${\mathcal {M}}_{n}=\left[{\begin{array}{cccc}\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n})&\mathbf {GG} ^{\intercal }&\mathbf {f} _{t}(t_{n},\mathbf {y} _{n})&\mathbf {f} (t_{n},\mathbf {y} _{n})\\\mathbf {0} &-\mathbf {f} _{\mathbf {x} }^{\intercal }(t_{n},\mathbf {y} _{n})&\mathbf {0} &\mathbf {0} \\\mathbf {0} &\mathbf {0} &0&1\\\mathbf {0} &\mathbf {0} &0&0\end{array}}\right]\in \mathbb {R} ^{(2d+2)\times (2d+2)},$

и $\{\mathbf {\xi } _{n}\}$ представляет собой последовательность d -мерных независимых двухточечных распределенных случайных векторов, удовлетворяющих $P(\xi _{n}^{k}=\pm 1)={\frac {1}{2}}$.

Схема WLL заказа 2

$\mathbf {y} _{n+1}=\mathbf {y} _{n}+\mathbf {B} _{16}+(\mathbf {B} _{14}\mathbf {B} _{11}^{\intercal })^{1/2}\mathbf {\xi } _{n},$

где $\mathbf {B} _{11}$, $\mathbf {B} _{14}$ и $\mathbf {B} _{16}$ являются подматрицами, определяемыми разделенной матрицей $\mathbf {B} =\mathbf {P} _{p,q}(2^{-k_{n}}{\mathcal {M}}_{n}h_{n}))^{2^{k_{n}}}$ с

${\mathcal {M}}_{n}=\left[{\begin{array}{cccccc}\mathbf {J} &\mathbf {H} _{2}&\mathbf {H} _{1}&\mathbf {H} _{0}&\mathbf {a} _{2}&\mathbf {a} _{1}\\\mathbf {0} &-\mathbf {J} ^{\intercal }&\mathbf {I} &\mathbf {0} &\mathbf {0} &\mathbf {0} \\\mathbf {0} &\mathbf {0} &-\mathbf {J} ^{\intercal }&\mathbf {I} &\mathbf {0} &\mathbf {0} \\\mathbf {0} &\mathbf {0} &\mathbf {0} &-\mathbf {J} ^{\intercal }&\mathbf {0} &\mathbf {0} \\\mathbf {0} &\mathbf {0} &\mathbf {0} &\mathbf {0} &0&1\\\mathbf {0} &\mathbf {0} &\mathbf {0} &\mathbf {0} &0&0\end{array}}\right]\in \mathbb {R} ^{(4d+2)\times (4d+2)},$

$\mathbf {J} =\mathbf {f} _{\mathbf {x} }(t_{n},\mathbf {y} _{n})\qquad \mathbf {a} _{1}=\mathbf {f} (t_{n},\mathbf {y} _{n})\qquad \mathbf {a} _{2}=\mathbf {f} _{t}(t_{n},\mathbf {y} _{n})+{\frac {1}{2}}\sum \limits _{i=1}^{m}(\mathbf {I} \otimes (\mathbf {g} ^{i}(t_{n}))^{\intercal })\mathbf {f} _{\mathbf {xx} }(t_{n},\mathbf {y} _{n})\mathbf {g} ^{i}(t_{n})$

и

$\mathbf {H} _{0}=\mathbf {G} (t_{n})\mathbf {G} ^{\intercal }(t_{n})\qquad \mathbf {H} _{1}=\mathbf {G} (t_{n}){\frac {d\mathbf {G} ^{\intercal }(t_{n})}{dt}}+{\frac {d\mathbf {G} (t_{n})}{dt}}\mathbf {G} ^{\intercal }(t_{n})\qquad \mathbf {H} _{2}={\frac {d\mathbf {G} (t_{n})}{dt}}{\frac {d\mathbf {G} ^{\intercal }(t_{n})}{dt}}{\text{.}}$

Стабильность и динамика

Рис. 5. Приблизительное среднее значение СДУ (8.2), вычисленное методом Монте-Карло при 100 симуляциях различных схем с h=1/16 и p=q=6 .

По построению слабые LL-дискретизации наследуют устойчивость и динамику линейных СДУ, но это не относится к слабым LL-схемам в целом. Схемы WLL, с $p\leq q\leq p+2,$ сохраняет первые два момента линейных СДУ и наследует среднеквадратическую устойчивость или нестабильность, которую может иметь такое решение. Сюда относятся, например, уравнения связанных гармонических осцилляторов, управляемых случайной силой, а также большие системы жестких линейных СДУ, возникающие на основе метода линий для линейных стохастических уравнений в частных производных. Более того, эти схемы WLL сохраняют эргодичность линейных уравнений и являются геометрически эргодическими для некоторых классов нелинейных СДУ. Для нелинейных СДУ с малым шумом (т. е. (8.1) с $\mathbf {g} _{i}(t)\approx 0$), решения этих схем WLL в основном представляют собой неслучайные пути схемы LL (4.6) для ОДУ плюс небольшое возмущение, связанное с малым шумом. В этой ситуации динамические свойства этой детерминированной схемы, такие как сохранение линеаризации и сохранение точной динамики решения вокруг гиперболических точек равновесия и периодических орбит, становятся актуальными для среднего значения схемы WLL. Например, на рис. 5 показано приблизительное среднее значение СДУ.

$dx=-t^{2}x{\text{ }}dt+{\frac {3}{2(t+1)}}e^{-t^{3}/3}{\text{ }}dw_{t},\qquad \qquad x(0)=1,\qquad \quad (8.2)$

рассчитывается по различным схемам.

Исторические заметки

Ниже представлена ​​временная шкала основных разработок метода локальной линеаризации (LL).

• Поуп Д.А. (1963) представляет дискретизацию LL для ОДУ и схему LL, основанную на расширении Тейлора.
• Озаки Т. (1985) представляет метод LL для интегрирования и оценки СДУ. Термин «локальная линеаризация» используется впервые.
• Бискай Р. и др. (1996) переформулировали сильный метод ЛЛ для СДУ.
• Сёдзи И. и Озаки Т. (1997) переформулируют слабый метод LL для СДУ.
• Хохбрук М. и др. (1998) представили схему LL для ОДУ, основанную на аппроксимации подпространств Крылова.
• Хименес Х.К. (2002) представляет схему LL для ОДУ и СДУ, основанную на рациональной аппроксимации Паде.
• Карбонелл FM и др. (2005) представили метод LL для RDE.
• Хименес Дж.С. и др. (2006) представили метод LL для DDE.
• Де ла Круз Х. и др. (2006, 2007) и Токман М. (2006) представляют два класса интеграторов HOLL для ОДУ: основанные на интеграторах и квадратурный.
• Де ла Круз Х. и др. (2010) представили сильный метод ХОЛЛА для СДУ.