Уравнение стохастического осциллятора

Системы стохастических уравнений <<	Оглавление	>> Линейные многомерные модели

---

$\textstyle \bullet$ В качестве примера решения стохастической задачи в двух измерениях $\textstyle n=m=2$ рассмотрим движение по окружности с частотой $\textstyle \omega$ и постепенным уменьшением радиуса. Детерминированная версия подобной спирали может иметь следующую зависимость координат от времени:

Примеры систем с таким поведением мы рассмотрим в следующей главе. Сейчас наша цель — математическое описание стохастической динамики. Для этого найдём решение системы уравнений следующего вида:

$\left\{{\begin{array}{l}dx=(-\lambda \,x-\omega \,y)\,dt\;+\;\sigma \,\delta W_{x}\\dy=(+\omega \,x-\lambda \,y)\,dt\;+\;\sigma \,\delta W_{y}.\\\end{array}}\right.$

Предполагается, что шум $\textstyle \delta W_{x}=\varepsilon _{x}{\sqrt {t}}$, $\textstyle \delta W_{y}=\varepsilon _{y}{\sqrt {t}}$ по каждой координате нескоррелирован. В качестве упражнения ($\textstyle \lessdot$ H) стоит записать эту же систему для скоррелированного шума.

$\textstyle \bullet$ Зависимость среднего значения от времени находим из ():

$\left\{{\begin{array}{l}{\dot {\overline {x}}}=-\lambda \,{\overline {x}}-\omega \,{\overline {y}}\\{\dot {\overline {y}}}=+\omega \,{\overline {x}}-\lambda \,{\overline {y}}.\\\end{array}}\right.$

Умножим второе из уравнений на мнимую единицу $\textstyle \imath$ ($\textstyle \imath ^{2}=-1$) и сложим их. В результате для комплексной величины $\textstyle z=x+\imath y$ и параметра $\textstyle \Lambda =\lambda -\imath \omega$ получим "одномерное" уравнение $\textstyle {\dot {\overline {z}}}=-\Lambda \,{\overline {z}}.$ Оно легко интегрируется:

${\overline {z}}(t)=e^{-\Lambda t}\,z_{0}=e^{-\lambda \,t+\imath \omega \,t}\,z_{0}=e^{-\lambda t}(\cos \omega t+\imath \sin \omega t)z_{0},$

где $\textstyle z_{0}=z(0)=x_{0}+\imath y_{0}$ — начальное условие. Приравняв действительную и мнимую части:

$\left\{{\begin{array}{l}{\overline {x}}(t)=e^{-\lambda t}\cdot (x_{0}\cos \omega t-y_{0}\sin \omega t)\\{\overline {y}}(t)=e^{-\lambda t}\cdot (x_{0}\sin \omega t+y_{0}\cos \omega t),\\\end{array}}\right.$

(EQN)

получаем решение для эволюции средних значений координат в виде спирали. По каждой координате происходят затухающие периодические колебания. Параметр $\textstyle \omega$ является их частотой, $\textstyle \lambda$ — скоростью затухания.

$\textstyle \bullet$ Найдём теперь, как ведёт себя среднее значение квадрата расстояния от начала координат. Воспользуемся матричной записью уравнения () для среднего Невозможно разобрать выражение (синтаксическая ошибка): {\displaystyle \textstyle \dot{\left\langle \mathbf{x}^2\right\rangle }=2\left\langle \mathbf{x} \cdot\mathbf{a} \right\rangle + \mathrm{Tr}\,left\langle \mathbf{b} \cdot\mathbf{b}^{T} \right\rangle } . В нашем случае:

$\mathbf {x} ={\begin{pmatrix}x\\y\end{pmatrix}},\;\;\;\;\;\;\;\;\;\mathbf {a} ={\begin{pmatrix}-\lambda x-\omega y\\+\omega x-\lambda y\\\end{pmatrix}},\;\;\;\;\;\;\;\;\;\mathbf {b} =\sigma {\begin{pmatrix}1&0\\0&1\\\end{pmatrix}}.$

Поэтому получаем:

${\dot {\overline {\mathbf {x} ^{2}}}}=-2\lambda \,{\overline {\mathbf {x} ^{2}}}+2\sigma ^{2}\;\;\;\;\;\;=>\;\;\;\;\;\;\;{\overline {\mathbf {x} ^{2}}}(t)={\frac {\sigma ^{2}}{\lambda }}+\left(x_{0}^{2}+y_{0}^{2}-{\frac {\sigma ^{2}}{\lambda }}\right)\,e^{-2\lambda t}.$

Таким образом, когда колебания "затухнут" ($\textstyle t\to \infty$), останется ненулевая дисперсия квадрата, тем большая, чем медленнее происходило затухание! Этот результат свидетельствует о стабилизирующей роли сильного трения (большие $\textstyle \lambda$) при внешних стохастических толчках. Несмотря на то, что $\textstyle \lambda$ находится в знаменателе, особенности при $\textstyle \lambda =0$ нет. Раскладывая экспоненту, несложно убедиться, что при $\textstyle \lambda =0$ среднее равно $\textstyle {\overline {\mathbf {x} ^{2}}}=x_{0}^{2}+y_{0}^{2}+2\sigma ^{2}t$. Это означает, что, как и в винеровском блуждании, внешние толчки со временем размывают круговую траекторию. Найти асимптотическое решение $\textstyle {\overline {x^{2}}}(t)$, $\textstyle {\overline {y^{2}}}(t)$, $\textstyle {\overline {xy}}(t)$ в ситуации скоррелированного шума предлагается в качестве упражнения ($\textstyle \lessdot$ H).

$\textstyle \bullet$ Выразим решение задачи через гауссовы переменные. В комплексных обозначениях стохастическое уравнение имеет компактный вид:

$dz=-\Lambda \,z\,dt+\sigma \,\delta W,$

где $\textstyle \delta W=\varepsilon \,{\sqrt {dt}}$, $\textstyle \varepsilon =\varepsilon _{x}+\imath \varepsilon _{y}$ — комплексное гауссово число, а $\textstyle z$ и $\textstyle \Lambda$ определены выше. Перейдём, при помощи формулы Ито, к переменной $\textstyle F=ze^{\Lambda t}$. Её динамическое уравнение не будет содержать сноса:

$dF=\sigma e^{\Lambda t}\delta W=S(t)\,\delta W,$

где $\textstyle S(t)=\sigma \,e^{\Lambda t}$. Решим уравнение $\textstyle dF=S(t)\,\delta W$ итерациями ($\textstyle k=1...n$):

$F=F_{0}+\sum S(t_{k-1})\varepsilon (t_{k}){\sqrt {\Delta t}}.$

Как функция $\textstyle S(t)$, так и $\textstyle \varepsilon _{k}$ являются комплексными величинами, поэтому необходима определённая осторожность по сворачиванию этой суммы в одно гауссово число. Распишем действительные и мнимые части:

$\sum [S_{x}+\imath S_{y}][\varepsilon _{x}+\imath \varepsilon _{y}]{\sqrt {\Delta t}}=\sum [(\underbrace {S_{x}\varepsilon _{x}-S_{y}\varepsilon _{y}} _{|S|\,\varepsilon '_{x}})+\imath (\underbrace {S_{x}\varepsilon _{y}+S_{y}\varepsilon _{x}} _{|S|\,\varepsilon '_{y}})]{\sqrt {\Delta t}},$

где $\textstyle |S|={\sqrt {S_{x}^{2}+S_{y}^{2}}}$ — модуль комплексного числа, а $\textstyle \varepsilon '_{x}$, $\textstyle \varepsilon '_{x}$ — новые нескоррелированные $\textstyle \left\langle \varepsilon '_{x}\varepsilon '_{y}\right\rangle =0$ гауссовы числа.

Опуская штрихи и повторяя рассуждения одномерного случая (стр. \pageref{ito_only_t_solution}), мы можем окончательно записать:

$F(t)=F(t_{0})+\left[\int \limits _{t_{0}}^{t}|S(\tau )|^{2}d\tau \right]^{1/2}\cdot \varepsilon ,$

где $\textstyle \varepsilon =\varepsilon _{x}+\imath \varepsilon _{y}$ — по-прежнему комплексная гауссова случайная величина.

Заметим, что действительная и мнимая части выражения $\textstyle \varepsilon '=e^{\imath \alpha }\varepsilon$ являются независимыми гауссовыми числами. Действительно, распишем их в явном виде:

${\begin{array}{l}\varepsilon '_{x}=\varepsilon _{x}\cos \alpha -\varepsilon _{y}\sin \alpha \\\varepsilon '_{y}=\varepsilon _{x}\sin \alpha +\varepsilon _{y}\cos \alpha \end{array}}$

Прямым вычислением проверяем Невозможно разобрать выражение (Ошибка преобразования. Сервер («https://wikimedia.org/api/rest_») сообщил: «Cannot get mml. TeX parse error: Double exponent: use braces to clarify»): {\displaystyle \textstyle \left\langle \varepsilon '_{x}^{2}\right\rangle =\left\langle \varepsilon '_{y}^{2}\right\rangle =1} и $\textstyle \left\langle \varepsilon '_{x}\varepsilon '_{y}\right\rangle =0$. Поэтому множители типа $\textstyle e^{\imath \omega t}$ перед комплексным гауссовым числом можно опускать, так как $\textstyle e^{\imath \alpha }\varepsilon$ статистически эквивалентно просто $\textstyle \varepsilon$

Проводя интегрирование для $\textstyle |S(\tau )|^{2}=\sigma ^{2}\,e^{2\lambda \tau }$ и учитывая, что $\textstyle z=F\cdot e^{-\Lambda t}$, $\textstyle z_{0}=F_{0}$, для $\textstyle t_{0}=0$ получаем:

$z=z_{0}\,e^{-\lambda t+\imath \omega t}+{\frac {\sigma }{\sqrt {2\lambda }}}\cdot {\sqrt {1-e^{-2\lambda t}}}\cdot \varepsilon$

(EQN)

или в явном виде для действительной и мнимой частей: \begin{eqnarray*} x(t) &=& \overline{x}(t) + \frac{\sigma}{\sqrt{2\lambda}} \, \sqrt{1-e^{-2\lambda t}} \cdot \varepsilon_x\\ y(t) &=& \overline{y}(t) + \frac{\sigma}{\sqrt{2\lambda}} \, \sqrt{1-e^{-2\lambda t}} \cdot \varepsilon_y, \end{eqnarray*} где $\textstyle {\overline {x}}(t)$ и $\textstyle {\overline {y}}(t)$ — средние, определяемые выражениями (). В качестве упражнения стоит найти $\textstyle {\overline {x^{2}}}(t)$, $\textstyle {\overline {y^{2}}}(t)$, $\textstyle {\overline {xy}}(t)$ и проверить справедливость уравнений для средних ($\textstyle \lessdot$ H).

Квадрат величины $\textstyle |z|^{2}=x^{2}+y^{2}$ является квадратом радиус-вектора, для которого уже известно среднее значение. Сделаем это ещё раз при помощи ():

$\left\langle |z_{t}|^{2}\right\rangle =|z_{0}|^{2}\,e^{-2\lambda t}+{\frac {\sigma ^{2}}{\lambda }}\cdot \left(1-e^{-2\lambda t}\right).$

Обращаем внимание на то, что $\textstyle \left\langle |\varepsilon |^{2}\right\rangle =\left\langle \varepsilon \varepsilon ^{*}\right\rangle =\left\langle \varepsilon _{x}^{2}+\varepsilon _{y}^{2}\right\rangle =2$, где звёздочка обозначает комплексное сопряжение.

В решении, аналогично одномерному случаю, можно выразить $\textstyle z$ в момент времени $\textstyle t+\tau$ через $\textstyle z$ в момент $\textstyle t$:

$z_{t+\tau }=z_{t}e^{-\lambda \tau +\imath \omega \tau }+{\frac {\sigma }{\sqrt {2\lambda }}}\,{\sqrt {1-e^{-2\lambda \tau }}}\cdot \varepsilon ,$

что легко позволяет вычислить, например, среднее $\textstyle \left\langle z_{t}\,z_{t+\tau }^{*}\right\rangle$ ($\textstyle \lessdot$ H).

$\textstyle \bullet$ При больших временах $\textstyle t\to \infty$ решение забывает начальные условия, и средние стремятся к нулю. Распределение становится стационарным по каждой из координат. Однако это не означает, что периодические свойства системы исчезают. Чтобы в этом убедиться, найдём ковариационную функцию, например, по координате $\textstyle x$. Выражая решение относительно начального момента времени $\textstyle t$, имеем:

$x_{t+\tau }=e^{-\lambda \tau }\cdot (x_{t}\cos \omega \tau -y_{t}\sin \omega \tau )+{\frac {\sigma }{\sqrt {2\lambda }}}\cdot \varepsilon _{x}\,{\sqrt {1-e^{-2\lambda \tau }}}$

Найдем $\textstyle \mathrm {cov} \,{xx}(t,t+\tau )=\left\langle x_{t}x_{t+\tau }\right\rangle -\left\langle x_{t}\right\rangle \left\langle x_{t+\tau }\right\rangle$ в пределе $\textstyle t\to \infty$. Так как в этом случае $\textstyle \left\langle x_{t}\right\rangle =\left\langle y_{t}\right\rangle =\left\langle x_{t}y_{t}\right\rangle =0$, а $\textstyle \left\langle x_{t}^{2}\right\rangle =\sigma ^{2}/2\lambda$, получаем, как и следовало ожидать, стационарную ковариационную функцию, зависящую только от $\textstyle \tau >0$:

$\mathrm {cov} \,{xx}(t,t+\tau )\to \mathrm {cov} \,\tau )={\frac {\sigma ^{2}}{2\lambda }}\,e^{-\lambda \tau }\cos \omega \tau .$

Она оказывается периодической функцией сдвига $\textstyle \tau$. Фурье-образ ковариационной функции характеризует спектральные свойства процесса (стр. \pageref{sect_autocor_fun}):

${\mathcal {S}}(\Omega )={\frac {2}{\pi }}\int \limits _{0}^{\infty }\mathrm {cov} \,\tau )\,\cos(\Omega \tau )\,d\tau ={\frac {\sigma ^{2}}{2\pi }}\,\left[{\frac {1}{\lambda ^{2}+(\Omega +\omega )^{2}}}+{\frac {1}{\lambda ^{2}+(\Omega -\omega )^{2}}}\right].$

Таким образом, спектр имеет максимум в окрестности $\textstyle \Omega =\omega$. Он тем уже, чем меньше параметр затухания $\textstyle \lambda$. Тем не менее, это не строго периодическое движение, так как "типичная" частота "размазана" и сдвинута первым слагаемым в квадратных скобках.

На левом рисунке ниже представлена траектория стохастического осциллятора при достаточно больших временах, когда начальные условия уже "забыты". Справа — колебания по каждой координате:

Системы, обладающие подобным поведением, мы рассмотрим в следующей главе, а сейчас решим многомерное линейное уравнение.

---

Системы стохастических уравнений <<	Оглавление	>> Линейные многомерные модели

---

Стохастический мир - простое введение в стохастические дифференциальные уравнения