首先,你的初始存量和流出量是长度单位,我猜应该是毫升吧?
关于随机游走的问题,1)完全随机游走,经典应用是股票价格了。2)有吸收壁游走,经典应用是破产概率模型了。3)弹射壁游走,这个应用不熟悉。4)就是你这种了,有短暂停留的吸收壁。
对于求解的理论方法,1、2大家都很容易找到资料,3参照Steven E.Shreve给出的反射不等式,4我就不知道了。
我感觉2楼的有点问题,就像上面说的2、4是不一样的。她给的很专业,可能是我猜错了。下面简单说一下4。
比如在第1秒内,无降雨。第2秒内有一个雨滴体积为0.3ml,冲击时刻为1.8秒。那么在第一秒结束时存量为1+0-1=0。在第2秒内的前0.8秒是没有任何流出的,后0.2秒的流出速度是1ml/s,共流出min(1,0.3)*(1-0.8)=0.06ml,剩余0.24ml。如果在1s内有2、3等多个雨点到达,就要计算不同的冲击间隔,每个间隔的初始存量、流出量
我不知道dual ruin problem是不是描述这样的问题。
根据其次泊松过程和指数分布的关系,$\lambda$=1,所以每次冲击时间间隔服从$\lambda'$=1的指数分布。
由于本人的知识有限,不能给出纯理论解,所以使用了随机模拟的方法,matlab程序见下方。
以100模拟次数基数,2倍增长模拟次数,结果为:0.42,0.265,0.3525,0.3675,收敛比率为:-13%,-3%,3%。因此大致估计为零的概率应是在0.36左右,而非你朋友和楼上给出的0概率。如果中间给出的物理过程描述没有问题的话,而我对我的的程序比较有信息,所以比较相信这个结果。
此外,你的问题来源可以说明一下吗?我想把这个问题记录在我的手册中,因为比较有意思。
clear;clc
T=60*10; %60秒×10分钟
SimNum=100*2.^[0,1,2,3]; %每模拟SimNum次计算一次P
RecordP=nan(1,4);
RunVol=1; %流出速度
StartVol=1; %初始存量
for k=1:4
EndVol=nan(1,SimNum(k)); %记录最后时刻存量
for j=1:SimNum(k) %模拟SimNum(k)次
RainVolume=zeros(1,T+1); %总时刻+初始时刻,坐标不允许为零。j表示j-1时刻
RainVolume(1)=StartVol; %零时刻等于初始存量
for t=2:T+1 %模拟后续时刻
RainNum=random('poisson',1); %到达的雨滴个数
if RainNum > 0 %有雨点到达
RainPoint=random('gamma',100,0.02,1,RainNum)/(2*10); %雨点半径,cm
RandVol=RainPoint.^3*4*pi/3; %雨点体积,按照球体计算
TimeInterval=2*ones(1,RainNum); %每个雨点冲击时间间隔
while sum(TimeInterval) > 1 %总间隔之和不能超过1s,条件事件
TimeInterval=random('exp',1,1,RainNum); %齐次泊松过程和指数分布
end
if sum(TimeInterval) < 1 %判断每秒内的最后一个雨点是否为最后时刻到来
TimeInterval=[TimeInterval,1-sum(TimeInterval)];
RainNum=RainNum+1; %非最后时刻到来,需加入到计算流出量的
RandVol=[RandVol,0]; %计算统一
end
else
RandVol=0; %没有雨点到来
end
if RainVolume(t-1) >= RunVol %上一时刻存量大于流出速度,直接扣减
RainVolume(t)=RainVolume(t-1)-RunVol+sum(RandVol);
elseif RainNum > 0 %有雨点到来
RainVolume(t)=RainVolume(t-1);
for s=1:RainNum %在每秒内,按照冲击时刻划分小时间段,计算流出量
RainVolume(t)=RainVolume(t)-TimeInterval(s)*min(RainVolume(t),RunVol);
RainVolume(t)=RainVolume(t)+RandVol(s); %小时间段结束时,下一粒雨点到来。
end
else
RainVolume(t)=0; %上1s结束时存量为0,并且这1s内无雨滴到来
end
end
EndVol(j)=RainVolume(T+1); %T时刻存量
end
RecordP(k)=sum(EndVol==0)/SimNum(k); %T时刻存量为0的总数占比,即概率
end