2

我正在尝试在mathematica/matlab中实现一个随机过程的例程。此处编写的任何代码都是为mathematica 编写的,但如果有人可以帮助我在matlab 中对此进行编码(如果他们更熟悉的话),那也可以。但是,如果可能的话,mathematica 是优先考虑的。

我会先说明一下方程式后我想得到什么。

以下是我感兴趣的伊藤随机过程(其中z(t)=[x(t),y(t)]:

在此处输入图像描述

具有以下数量:

在此处输入图像描述

我们还有以下内容(从现在开始我将输入x(t),y(t) as x,y):

在此处输入图像描述

其中,在此处输入图像描述是第一次退出时间, GT(x)x(请参阅最后)的平滑函数

目标:我想得到Q0x在此处输入图像描述只有。

==================================================== ============================= 要找到第一个退出时间在此处输入图像描述,也许可以使用以下方法(由 b.gatessucks 提供)。请注意,以下是mathematica 代码。

退出时间的约束

在此处输入图像描述

const[x_, y_] := And[10^-8 <= y <= 10^-3, 0.9*(Uc) <= a1*x^2/2 - a3*x^4/4 <= 1.1*(Uc)]

x0 = 0.1;   (* starting point for x[t] *) y0 = 0.1; (* starting point for y[t] *)

proc = ItoProcess[  {\[DifferentialD]x[t] == y[t] \[DifferentialD]t,   \[DifferentialD]y[t] == (-G*y[t] - (a1*x[t] - a3*x[t]^3) - eps*b3b*y[t]^3) \[DifferentialD]t + Sqrt[2*eps*G] \[DifferentialD]w[t]},    {t, x[t], y[t], Boole[const[x[t], y[t]]]},   {{x, y}, {x0, y0}}, {t, 0}, w \[Distributed] WienerProcess[]]

退出时间找到:

SeedRandom[3] 
sim = RandomFunction[proc, {0, 1, 0.001}];
First@Select[sim[[2, 1, 1]], #[[4]] == 0 &]

输出{t, x[t], y[t], Boole[const[x[t], y[t]]]}

我需要能够使用上面的代码找到退出时间,然后在Q0. 上面查找退出时间的片段会为正确选择的在此处输入图像描述.

==================================================== ===============================

需要设置的数值任务 find Q0(x)

--> 我们从这个表达式中的积分项开始。首先 find在此处输入图像描述从域内部开始的许多初始条件(现在说 100 - 即 100 次退出时间)(可能使用本文中已经提到的代码段)。现在,可以将 100 个退出时间中的每一个的积分在此处输入图像描述计算为 和 的x函数在此处输入图像描述

--> 现在,积分的期望值Q0是之前评估的所有 100 个积分的样本平均值(根据大数定律)。因此,Q0可以找到并且它应该只是 和 的x函数在此处输入图像描述

==================================================== ===============================

我的问题:

  1. The code above for the exit time only seems to produce non-zero exit times for appropriately chosen initial conditions. 如果有人可以阐明如何选择合适的在此处输入图像描述,以便产生足够的退出时间,那将不胜感激。我真的很想保持上面在 中指定的约束const[x_, y_],但如果似乎没有希望找到易于处理的结果,那么我不介意放松它。

  2. 下面的代码GT(x)导致奇点——我从 DSolve 收到关于不确定表达式的错误消息.... b0[xc] 和 DrhoDy[xc] 都变得不确定,因此 DSolve 在使用初始条件 GT[xc] 时会出现问题它也变得不确定......解决这个问题的任何方式都会很高兴。

  3. 最后,我真的需要有人的帮助来有效地评估数学中的 100 个积分(因为这些项很大)对于之前找到的每个退出时间并采取预期。我不确定如何正确找到退出时间。

==================================================== ============================== 求GT(x):

b1b = 0.9; b3b = .8; a1b = 0.1; a3b = 0.2; G = (1/0.1^2)*
  b1b; a1 = (1/.1^2)*a1b; a3 = (1/.1^2)*a3b; xc = Sqrt[a1/a3]; Uc = 
 a1*xc^2/2 - a3*xc^4/4; L = (G + Sqrt[G^2 + 4*(-(a1 - 3*a3*xc^2))])/2;
y[x_] := (x *a1 - x^3*a3)/(G + L)

b0[x_] := (y[
      x]*(a1*x \[Minus] 
       a3*x ^3)*(1 \[Minus] (a1 \[Minus] 3*a3*x ^2)) \[Minus] 
    G*y[x]*(a1 \[Minus] 3*a3*x ^2)) /(y[
      x] ^2 + (\[Minus]G*y[x] + a1*x \[Minus] a3*x ^3) ^2)


DrhoDy [x_] := 
 Sqrt[y[x]^ 2 /(y[x] ^2 + (\[Minus]G*y[x] + a1*x \[Minus] a3*x ^3) ^2)]


linearequation = 
  y[x]*GT '[x] + b0[x]*GT[x] == G*DrhoDy [x]^2 *GT[x]^ 3;

GT[x] = 
 DSolve[{linearequation, GT[xc] == Sqrt[b0[xc]/( G*DrhoDy [xc]^2 )]}, 
  GT, x]

颂:

在此处输入图像描述

满足初始条件:

在此处输入图像描述

在哪里,

在此处输入图像描述

在此处输入图像描述

4

0 回答 0