Колебательный контур под внешним воздействием периодического сигнала и белого шума 

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

In[1]:=

X CapacitorChargeModel = ItoProcess[{ \[DifferentialD]\[ScriptCapitalQ][t] == \[ScriptCapitalI][ t] \[DifferentialD]t, \[ScriptCapitalL] \[DifferentialD]\[ScriptCapitalI][ t] + \[ScriptCapitalR] \[ScriptCapitalI][ t] \[DifferentialD]t + 1/\[ScriptCapitalC] \[ScriptCapitalQ][t] \[DifferentialD]t == Sin[\[Omega] t] \[DifferentialD]t + \[Sigma] \[DifferentialD]w[ t]}, \[ScriptCapitalQ][ t], {{\[ScriptCapitalI], \[ScriptCapitalQ]}, {0, 0}}, t, w \[Distributed] WienerProcess[]]

Out[1]=

Зададим значения параметров модели.

In[2]:=

X params = {\[ScriptCapitalL] -> 2, \[ScriptCapitalC] -> 2, \[ScriptCapitalR] -> 1/2, \[Omega] -> 3/2, \[Sigma] -> 1/4};

Найдём функцию среднего, то есть среднее значение заряда конденсатора как функции от времени.

In[3]:=

X mqf = Simplify[Mean[CapacitorChargeModel[t] /. params]]

Out[3]=

Решим детерминистическую задачу в отсутствии шума, и сравним решение с функцией среднего.

In[4]:=

X detq = (Q[t] /. First[DSolve[\[ScriptCapitalL] Q''[t] + \[ScriptCapitalR] Q'[t] + Q[t]/\[ScriptCapitalC] == Sin[\[Omega] t] \[And] Q[0] == 0 && Q'[0] == 0 /. params, Q[t], t]]);

In[5]:=

X Simplify[detq == mqf]

Out[5]=

Найдём функцию дисперсии.

In[6]:=

X vf = Simplify[Variance[CapacitorChargeModel[t] /. params]]

Out[6]=

Построим случайные реализации модели, используя схему Клодена-Платена-Шурца, с порядком сильной сходимости .

In[7]:=

X path = RandomFunction[ CapacitorChargeModel /. params, {0., 22., 0.002}, 12, Method -> "KloedenPlatenSchurz"]

Out[7]=

Визуализация случайных траекторий.

In[8]:=

X pls = ListLinePlot[path, PlotRange -> All, ImageSize -> 400]

Out[8]=

Визуализация 95% доверительного интервала относительно функции среднего.

In[9]:=

X epl = With[{q = Quantile[NormalDistribution[], (1 + .95)/2]}, Plot[{mqf - q Sqrt[vf], mqf, mqf + q Sqrt[vf]} /. params // Evaluate, {t, 0, 22}, PlotRange -> All, PlotStyle -> Directive[Thick, Black], FillingStyle -> Directive[Lighter[Blue, .5], Opacity[.15]], Filling -> {1 -> {3}}, ImageSize -> 400]]

Out[9]=

Сочетание этих графиков.

In[10]:=

X Show[{epl, pls}, PlotRange -> All]

Out[10]=