Стохастический осциллятор

Теория броуновского движения <<	Оглавление	>> Дрожание земной оси

---

$\textstyle \bullet$ Множество механических, электромагнитных, биологических и социальных систем описываются осцилляторными уравнениями. Для определённости мы рассмотрим одномерный механический осциллятор массой $\textstyle m$, подверженный трению и внешним стохастическим воздействиям. Определение импульса и закон Ньютона в этом случае имеют вид:

[[File:oscillator.png

где сила состоит из трёх компонент:

Невозможно разобрать выражение (синтаксическая ошибка): {\displaystyle F = - \underbrace{(k + {\rm Noise}_1)\cdot x}_{сила\;упругости} \;-\; \underbrace{(2\lambda + {\rm Noise}_2)\cdot p}_{сила\;трения} \;+\; \underbrace{{\rm Noise}_3.}_{внешняя\;сила}}

Сила упругости пропорциональна величине отклонения от положения равновесия $\textstyle x$. Мы будем считать, что коэффициент упругости $\textstyle k$ испытывает стохастические изменения, которые символически обозначены членом Noise$\textstyle _{1}$. Знак минус перед упругой силой означает, что она стремится вернуть частицу назад, к положению равновесия. Сила сопротивления тем больше, чем больше скорость (импульс) частицы. Так происходит при движении в среде (воздух, вода). Сопротивление стремится остановить движение. Будем также предполагать, что коэффициент сопротивления подвержен стохастическим воздействиям Noise$\textstyle _{2}$. Наконец, третья составляющая силы — это шум Noise$\textstyle _{3}$, который может быть, например, внешними случайными толчками.

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

Будем работать в системе единиц, для которой $\textstyle m=1$, $\textstyle k=1$ ($\textstyle \lessdot$ C). Стохастические уравнения движения в этом случае имеют вид:

$\left\{{\begin{array}{l}dx=p\,dt\\dp=-x\,dt-2\lambda \,p\,dt+\sigma _{1}\,x\,\delta W_{1}+\sigma _{2}\,p\,\delta W_{2}+\sigma _{3}\,\delta W_{3},\end{array}}\right.$

где $\textstyle \sigma _{1}$ — волатильность коэффициента упругости, $\textstyle \sigma _{2}$ — силы трения, а $\textstyle \sigma _{3}$ — внешнего шума. Винеровские переменные $\textstyle \delta W_{1}$, $\textstyle \delta W_{2}$ и $\textstyle \delta W_{3}$ представляют собой изменения трёх независимых процессов.

$\textstyle \bullet$ Рассмотрим сначала общий случай, записав систему:

$d\mathbf {x} =\mathbf {a} \,dt+\mathbf {b} \cdot \delta \mathbf {W} ,$

со следующими векторами и матрицами:

$\mathbf {x} ={\begin{pmatrix}x\\p\\\end{pmatrix}},\;\;\;\;\mathbf {a} ={\begin{pmatrix}p\\-x-2\lambda p\\\end{pmatrix}},\;\;\;\;\mathbf {b} ={\begin{pmatrix}0&0&0\\\sigma _{1}x&\sigma _{2}p&\sigma _{3}\\\end{pmatrix}},\;\;\;\;\delta \mathbf {W} ={\begin{pmatrix}\delta W_{1}\\\delta W_{2}\\\delta W_{3}\\\end{pmatrix}}.$

Для функции $\textstyle F(\mathbf {x} )=F(x,p)$ координат и импульсов воспользуемся динамическим уравнением для средних (стр. \pageref{aver_F_n_dim}):

${\frac {d}{dt}}{\left\langle F(\mathbf {x} )\right\rangle }=\left\langle {\frac {\partial F}{\partial \mathbf {x} }}\cdot \,\mathbf {a} +{\frac {1}{2}}\mathrm {Tr} \,\left[\mathbf {b} ^{T}\cdot {\frac {\partial ^{2}F}{\partial \mathbf {x} ^{2}}}\cdot \mathbf {b} \right]\right\rangle ,\;\;\;\;\;\;\;\;\;{\frac {\partial ^{2}F}{\partial \mathbf {x} ^{2}}}={\begin{pmatrix}F_{xx}&F_{xp}\\F_{px}&F_{pp}\end{pmatrix}},$

где $\textstyle F_{xx}$ — вторая производная по $\textstyle x$, $\textstyle F_{xp}$ — производная по $\textstyle x$ и $\textstyle p$, и т.д. Подставляя матрицы и перемножая их, получаем ($\textstyle \lessdot$ H):

Невозможно разобрать выражение (синтаксическая ошибка): {\displaystyle \frac{d}{dt}\, \Bigl<F(x,p)\Bigr> =\Bigl< p F_x - (x+2\lambda p)F_p\Bigr> + \frac{1}{2}\Bigl<(\sigma^2_1 x^2 + \sigma^2_2 p^2 + \sigma^2_3) F_{pp}\Bigr>. }

Выбор $\textstyle F=x$ и $\textstyle F=p$ приводит к системе уравнений, совпадающих с детерминированными (снос линеен!):

$\left\{{\begin{array}{l}{\dot {\left\langle x\right\rangle }}=\;\,\left\langle p\right\rangle \\{\dot {\left\langle p\right\rangle }}=-\left\langle x\right\rangle -2\lambda \,\left\langle p\right\rangle .\\\end{array}}\right.$

Её решение с начальными условиями $\textstyle x_{0}=x(0)$, $\textstyle p_{0}=p(0)$ имеет вид:

$\left\{{\begin{array}{l}\left\langle x\right\rangle =\left(x_{0}\,\cos \omega t+{\frac {p_{0}+\lambda x_{0}}{\omega }}\,\sin \omega t\right)\,e^{-\lambda t}\\\left\langle p\right\rangle =\left(p_{0}\,\cos \omega t-{\frac {x_{0}+\lambda p_{0}}{\omega }}\,\sin \omega t\right)\,e^{-\lambda t},\end{array}}\right.$

(EQN)

где $\textstyle \omega ={\sqrt {1-\lambda ^{2}}}$ (мы считаем, что трение мало и $\textstyle \lambda <1$). При выводе () можно воспользоваться алгоритмом на стр. \pageref{sec_line_n_dim_models} или привести систему к одному дифференциальному уравнению второго порядка ($\textstyle \lessdot$ H).

Выбор $\textstyle F=x^{2},\;p^{2},\;xp$ приводит к системе уравнений для моментов:

$\left\{{\begin{array}{l}{\dot {\left\langle x^{2}\right\rangle }}=2\left\langle xp\right\rangle \\{\dot {\left\langle xp\right\rangle }}=\left\langle p^{2}\right\rangle -\left\langle x^{2}\right\rangle -2\lambda \,\left\langle xp\right\rangle \\{\dot {\left\langle p^{2}\right\rangle }}=-2\left\langle xp\right\rangle +\sigma _{1}^{2}\left\langle x^{2}\right\rangle +(\sigma _{2}^{2}-4\lambda )\,\left\langle p^{2}\right\rangle +\sigma _{3}^{2}.\\\end{array}}\right.$

(EQN)

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

$\textstyle \bullet$ Если $\textstyle 4\lambda >\sigma _{1}^{2}+\sigma _{2}^{2}$, система имеет стационарный режим при $\textstyle t\to \infty$, в котором:

$\left\langle x^{2}\right\rangle =\left\langle p^{2}\right\rangle ={\frac {\sigma _{3}^{2}}{4\lambda -\sigma _{1}^{2}-\sigma _{2}^{2}}},\;\;\;\;\;\;\;\;\;\;\left\langle xp\right\rangle =0.$

(EQN)

При $\textstyle \lambda >0$ средние стремятся к нулю, и матрица дисперсии оказывается диагональной:

$\mathbf {D} ={\frac {\sigma _{3}^{2}}{4\lambda -\sigma _{1}^{2}-\sigma _{2}^{2}}}{\begin{pmatrix}1&0\\0&1\\\end{pmatrix}}.$

Каждая разновидность шума увеличивает дисперсии, но по-разному. Трение $\textstyle \lambda$ играет стабилизирующую роль, уменьшая $\textstyle \mathbf {D}$.

Заметим, что динамика при $\textstyle t\to \infty$ продолжается только, если существует внешний шум ($\textstyle \sigma _{3}\neq 0$). Если $\textstyle \sigma _{3}=0$, стационарный режим также существует, но он вырождается в полное затухание колебаний, и $\textstyle \left\langle x^{2}\right\rangle =\left\langle p^{2}\right\rangle =0$. Причина подобного поведения та же, что и у логистического уравнения (стр. \pageref{why_logistic_stop}).

$\textstyle \bullet$ Пусть детерминированной составляющей трения нет $\textstyle \lambda =0$, а флуктуации упругости и трения имеют одинаковые амплитуды $\textstyle \sigma _{1}=\sigma _{2}=\sigma$. Введём энергию гармонического осциллятора:

$E={\frac {x^{2}+p^{2}}{2}}.$

Из () следует, что её среднее значение удовлетворяет уравнению:

${\frac {d}{dt}}\left\langle E\right\rangle =\sigma ^{2}\left\langle E\right\rangle +{\frac {\sigma _{3}^{2}}{2}},$

а, следовательно, возрастает со временем:

$\left\langle E\right\rangle =\left(E_{0}+{\frac {\sigma _{3}^{2}}{2\sigma ^{2}}}\right)e^{\sigma ^{2}t}-{\frac {\sigma _{3}^{2}}{2\sigma ^{2}}},\;\;\;\;\;\;\;\;\;\;\;E_{0}={\frac {x_{0}^{2}+p_{2}^{2}}{2}}.$

Если стохастическое воздействие обусловлено только внешними толчками ($\textstyle \sigma _{1}=\sigma _{2}=0$), то рост не такой быстрый и аналогичен винеровскому увеличению неопределённости $\textstyle \left\langle E\right\rangle =E_{0}+\sigma _{3}^{2}\,t/2$. Так же, как и броуновская частица под внешним воздействием в среднем удаляется от начального положения, так и квадрат амплитуды осциллятора при $\textstyle \lambda =0$ в среднем увеличивается.

$\textstyle \bullet$ Если существуют только внешние толчки ($\textstyle \sigma _{1}=\sigma _{2}=0$), то стохастика имеет постоянную волатильность $\textstyle \sigma _{3}=\sigma$:

$\left\{{\begin{array}{l}dx=p\,dt\\dp=-x\,dt-2\lambda \,p\,dt+\sigma \,\delta W.\end{array}}\right.$

Подобную систему мы рассматривали в шестой главе (стр. \pageref{stochastic_oscillator}). Она обладает точным решением, которое выражается через две независимые гауссовы переменные. Воспользуемся общим алгоритмом решения системы линейных уравнений (см. стр. \pageref{sec_line_n_dim_models}) с матрицами:

$\mathbf {A} ={\begin{pmatrix}0&1\\-1&-2\lambda \\\end{pmatrix}},\;\;\;\;\;\;\mathbf {B} ={\begin{pmatrix}0&0\\0&\sigma \\\end{pmatrix}}.$

Чтобы найти $\textstyle e^{\mathbf {A} t}$, продифференцируем () по $\textstyle x_{0}$ и $\textstyle p_{0}$:

$e^{\mathbf {A} t}={\begin{pmatrix}\omega \cos \omega t+\lambda \sin \omega t&\sin \omega t\\-\sin \omega t&\omega \cos \omega t-\lambda \sin \omega t\\\end{pmatrix}}\cdot {\frac {e^{-\lambda t}}{\omega }}.$

При помощи этой матрицы, интегрируя (), стр. \pageref{n_sys_line_disp}, можно найти дисперсию координаты и импульса:

$\left\{{\begin{matrix}D_{xx}(t)\\D_{pp}(t)\\\end{matrix}}\right\}={\frac {\sigma ^{2}}{4\lambda }}-{\frac {\sigma ^{2}}{4\lambda \omega ^{2}}}\,{\bigl [}1-\lambda ^{2}\cos(2\omega t)\pm \lambda \omega \sin(2\omega t){\bigr ]}\,e^{-2\lambda t}.$

Верхний знак соответствует дисперсии для $\textstyle x$: $\textstyle D_{xx}=\left\langle x^{2}\right\rangle -\left\langle x\right\rangle ^{2}$, а нижний — для $\textstyle p$: $\textstyle D_{pp}=\left\langle p^{2}\right\rangle -\left\langle p\right\rangle ^{2}$. Дисперсия произведения динамических переменных $\textstyle D_{xp}(t)=\left\langle xp\right\rangle -\left\langle x\right\rangle \left\langle p\right\rangle$ имеет вид:

$D_{xp}(t)={\frac {\sigma ^{2}}{2\omega ^{2}}}\,\sin ^{2}(\omega t)\,e^{-2\lambda t}$

и стремится к нулю при $\textstyle t\to \infty$ и $\textstyle \lambda >0$. В результате, в стационарном режиме ($\textstyle t\to \infty$) матрица дисперсий диагональна (), поэтому автоковариационная матрица $\textstyle \mathrm {cov} \,\tau )$ равна $\textstyle e^{\mathbf {A} ^{T}|\tau |}$ с множителем $\textstyle \sigma ^{2}/4\lambda$.

При отсутствии трения $\textstyle \lambda =0$, $\textstyle \omega =1$:

$\mathbf {D} (t)={\frac {\sigma ^{2}}{2}}{\begin{pmatrix}t-\sin t\cos t&\sin ^{2}t\\\sin ^{2}t&t+\sin t\cos t\\\end{pmatrix}},\;\;\;\;\;\;\;e^{\mathbf {A} t}={\begin{pmatrix}\cos t&\sin t\\-\sin t&\cos t\\\end{pmatrix}}$

и, как мы видели выше, стационарного режима нет. Дисперсии по $\textstyle x$ и $\textstyle p$ растут во времени, совершая периодические колебания. Автоковариационная матрица $\textstyle \mathrm {cov} \,t,t+\tau )$ получается перемножением $\textstyle \mathbf {D} (t)$ и $\textstyle e^{\mathbf {A} ^{T}|\tau |}$.

---

Теория броуновского движения <<	Оглавление	>> Дрожание земной оси

---

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