#COPYRIGHT NOTICE: Copyright (c) 1997-1998 by Sasha Cyganowski. # 2003-2004 by Lars Gruene #Permission is #granted to anyone to use, modify, or redistribute this software freely, #subject to the following restrictions: 1. The author accepts no #responsibility for any consequences of this software and makes no guarantee #that this software is free of defects. 2. The origin of this software must #not be misrepresented, either by explicit claim or by omission. #3. This notice and the copyright must be included with all copies or #modified versions of this software. 4. This software may not be included or #redistributed as part of any package to be sold for profit unlesss the author #has given explicit written permission to do so. # Sasha Cyganowski, School of Computing and Mathematics, Deakin University, #Waurn Ponds, VIC, Australia. Internet: sash@deakin.edu.au # Telephone: +61 3 5227 1382. # additional routines and updates by Lars Gruene (lars.gruene@uni-bayreuth.de) with(stats): Wiener:=proc(t,p)\ local temp1,temp2,temp3,temp4,temp5,temp6,i; temp6[1]:=0; for i from 1 to p do temp6[i+1]:=temp6[i]+stats[random,normald[0,1]](1)*convert(sqrt(t),float); od; RETURN(convert(temp6,list)); end: stochplot:=proc(s,h,p,z) local i,j,temp1,temp2,Xdata,Ydata; temp1:=Wiener(h,p); temp2:=rhs(s[1]); Ydata[1]:=z; j:=1; for i from 1 to p do Ydata[i+1]:=evalf(subs(Y||j[n]=Ydata[i],Delta[n]=h,DeltaW||j[n]=temp1[i+1]-temp1[i],temp2)); od; Xdata[1]:=0; for i from 1 to p do Xdata[i+1]:=Xdata[i]+h; od; Xdata:=convert(Xdata,list); Ydata:=convert(Ydata,list); plots[display]({statplots[scatterplot](Xdata,Ydata)},axes=FRAME); end: stochplot_time:=proc(s,h,p,z,t0) local i,j,temp1,temp2,Xdata,Ydata; temp1:=Wiener(h,p); temp2:=rhs(s[1]); Ydata[1]:=z; j:=1; for i from 1 to p do Ydata[i+1]:=evalf(subs(Y||j[n]=Ydata[i],Delta[n]=h,DeltaW||j[n]=temp1[i+1]-temp1[i],t=t0+(p-1)*h,temp2)); od; Xdata[1]:=t0; for i from 1 to p do Xdata[i+1]:=Xdata[i]+h; od; Xdata:=convert(Xdata,list); Ydata:=convert(Ydata,list); plots[display]({statplots[scatterplot](Xdata,Ydata)},axes=FRAME); end: