![](/profiles/wavemetrics/themes/wavemetrics/logo.png)
Stochastic ODE: Euler-Maruyama illustrated
![](/sites/default/files/styles/thumbnail/public/default_images/Artboard%201_1.png?itok=jeHOCIXy)
bech
If others are interested in stochastic methods, I encourage you to post more examples. This will perhaps persuade the Igor folks to add explicit support for stochastic equations. (While these methods are easy to program in Igor as it stands, they are much slower than what a C-coded version could do. And these calculations tend to require LOTS of iterations.....)
Note that this was developed in 6.10b07, but I see no reason that it won't run in Igor 6.0x.
//________________________________________________________________ function em(M,N,D,Tf) // vector implementation of Euler-Maruyama variable M,N,D,Tf // M = paths, N = time points // D = noise strength, Tf = final sim. time // try M=1, N=1e4, D=0.01, Tf=1000 // SetRandomSeed 0.5 // useful for debugging variable dt = Tf/N, Xzero = 0 // dt = time step; Xzero = initial condition wave pars // declare parameter wave make /o/n=(M,N) dW dW = sqrt(2*D*dt)*gnoise(1) // M realizations of noise terms dW(t) Integrate /meth=2 /dim=1 dW /D=W // integral is just sum of dW's Duplicate /o W, Xem Xem[][0] = Xzero // initialize all paths Xem[][1,N] = Xem[p][q-1] + dt*myfunc(pars, (q-1)*dt, Xem[p][q-1]) + dW[p][q-1] // this is the crucial step for the EM method! MatrixOp /o Xem0 = row(Xem,0)^t // get a single trajectory to graph SetScale/P x 0,dt,"", Xem0 killwaves dW,W,Xem // get rid of big waves end //________________________________________________________________ Function myfunc(pw,t,x) wave pw; variable t,x // t for time-dependent forces; not needed here variable a=pw[0] // barrier height (try a=0.3) return a*x - x^3 // force = -dU/dx, U = -a*x^2/2 + x^4/4 End //________________________________________________________________ Window path_graphs() : Graph PauseUpdate; Silent 1 // building window... Display /W=(533,50,1004,389) Xem0 ModifyGraph width=360,height=216 ModifyGraph marker=19 ModifyGraph msize=1 ModifyGraph nticks(bottom)=2 ModifyGraph fSize=18 ModifyGraph manTick(left)={0,1,0,0},manMinor(left)={1,0} ModifyGraph manTick(bottom)={0,500,0,0},manMinor(bottom)={4,0} Label left "Particle position (x)" Label bottom "Time (t)" SetAxis left -1.5,1.5 ControlBar 42 SetVariable setvar0,pos={67,8},size={214,27},title="barrier height a ",fSize=18 SetVariable setvar0,limits={0,1,0.1},value= pars[0] EndMacro
![Sample output showing hops across a potential barrier](/sites/default/files/styles/content_body/public/images-imported/path_graphs.preview.jpg?itok=nwVjyNUz)
![](/sites/default/files/forum.png)
Forum
![](/sites/default/files/support.png)
Support
![](/sites/default/files/gallery.png)
Gallery
Igor Pro 9
Learn More
Igor XOP Toolkit
Learn More
Igor NIDAQ Tools MX
Learn More