Python科学计算：随机微分方程

本来随机微分方程之前，有一延迟微分方程，但是能够主要是在控制理论和数学生物学里面用的多，和我的关系不大，所以就当了解了，这个随机微分方程还是要看一看的。

1. 维纳过程(布朗运动):

这是一种随机行为建模机制，在上的标准标量过程或者运动，是连续依赖于t的随机变量W(t)，并且满足以下准则：

1、

2、如果，那么增量是正态分布的随机变量，其平均值是0，方差是t-s，

3、如果，那么增量和是独立变量

我们现在按着书上的例子，用500个步骤生成和绘制的布朗运动：

import numpy as np
from numpy import random as npr
import matplotlib.pyplot as plt
T=1
N=500
t,dt=np.linspace(0,T,N+1,retstep=True)
dw=npr.normal(0.0,np.sqrt(dt),N+1)
dw[0]=0.0
W=np.cumsum(dw)
fig=plt.figure()
ax=fig.add_subplot(1,1,1)
ax.plot(t,W)
ax.set_xlabel("t")
ax.set_ylabel("W(t)")
plt.show()

注意啊，过程是随机的，所以每次输出都是不一样的：

图1：第一次随机过程

图2：第二次随机过程

哎呀，实在是不好展示，我直接搞个多图吧：

import numpy as np
from numpy import random as npr
import matplotlib.pyplot as plt
T=1
N=500
t,dt=np.linspace(0,T,N+1,retstep=True,dtype=float)
dw1=npr.normal(0.0,np.sqrt(dt),N+1)
dw1[0]=0.0
dw2=npr.normal(0.0,np.sqrt(dt),N+1)
dw2[0]=0.0
dw3=npr.normal(0.0,np.sqrt(dt),N+1)
dw3[0]=0.0
dw4=npr.normal(0.0,np.sqrt(dt),N+1)
dw4[0]=0.0
dw5=npr.normal(0.0,np.sqrt(dt),N+1)
dw5[0]=0.0
dw6=npr.normal(0.0,np.sqrt(dt),N+1)
dw6[0]=0.0
W1=np.cumsum(dw1)
W2=np.cumsum(dw2)
W3=np.cumsum(dw3)
W4=np.cumsum(dw4)
W5=np.cumsum(dw5)
W6=np.cumsum(dw6)
fig=plt.figure()
ax1=fig.add_subplot(2,3,1)
ax1.plot(t,W1)
ax1.set_xlabel("t")
ax1.set_ylabel("W(t)")
ax1.set_title("1")
ax2=fig.add_subplot(2,3,2)
ax2.plot(t,W2)
ax2.set_xlabel("t")
ax2.set_ylabel("W(t)")
ax2.set_title("2")
ax3=fig.add_subplot(2,3,3)
ax3.plot(t,W3)
ax3.set_xlabel("t")
ax3.set_ylabel("W(t)")
ax3.set_title("3")
ax4=fig.add_subplot(2,3,4)
ax4.plot(t,W4)
ax4.set_xlabel("t")
ax4.set_ylabel("W(t)")
ax4.set_title("4")
ax5=fig.add_subplot(2,3,5)
ax5.plot(t,W5)
ax5.set_xlabel("t")
ax5.set_ylabel("W(t)")
ax5.set_title("5")
ax6=fig.add_subplot(2,3,6)
ax6.plot(t,W6)
ax6.set_xlabel("t")
ax6.set_ylabel("W(t)")
ax6.set_title("6")
plt.show()

图3：多次随机过程集中展示

好了，今天还有其他事，就先写到这。