1.3.1 to 1.3.2 Euler's Method for Solving ODEs Numerically, Euler's method using Matlab
#-----------------------------------[br]# Define ODE y' = f[br]#-----------------------------------[br]f(t,y) = sin(1-t)*(y-1)*(y+1)[br][br]#-----------------------------------[br]# Define increment[br]#-----------------------------------[br]Δt = Slider(0.01, 1, 0.01, 1, 400, false, true, false, false)[br]SetValue(Δt, 0.05)[br][br]#-----------------------------------[br]# Euler's method[br]#-----------------------------------[br]n = 50[br]Execute( Join( {"T0 = 0"}, Zip( "T"+ (k+1) +" = T"+ k +" + Δt", k, 0..n ) ) ) [br][br]Execute( Join( {"Y0 = 0.4"}, Zip( "Y"+ (k+1) +" = Y"+ k +"+ f(T"+ k +", Y"+ k +") * Δt", k, 0..n ) ) ) [br][br]Execute( Join( {"P0 = (T0, Y0)"}, Zip( "P"+ (k+1) +" = (T"+ (k+1) +", Y"+ (k+1) +")", k, 0..n ) ) ) [br][br]# 4th Degree Runge-Kutta solution with step = 0.001[br]sol = SolveODE(f, T0, Y0, 3, 0.001)[br][br]#-----------------------------------[br]# Settings[br]#-----------------------------------[br]Execute( Join( {"X0 = (T0, 0)"}, Zip( "X"+ (k+1) +" = (T"+ (k+1) +", 0)", k, 0..n ) ) )[br][br]# Texts[br]T = Zip("y_{"+k"}", k, 0..n)[br]U = Zip("t_{"+k"}", k, 0..n)[br]Execute(Zip("SetCaption(P"+ k +", Element(T, "+ (k+1) +"))", k, 0..n))[br]Execute(Zip("SetCaption(X"+ k +", Element(U, "+ (k+1) +"))", k, 0..n))[br][br]b = "Blue"[br]r = "Red"[br]Execute(Zip("SetColor(X"+ k +", b)", k, 0..n))[br]Execute(Zip("SetColor(P"+ k +", r)", k, 0..n))[br][br]Execute(Zip("SetPointSize(X"+ k +", 5)", k, 0..n))[br]Execute(Zip("SetPointSize(P"+ k +", 5)", k, 0..n))[br][br]Execute(Zip("SetPointStyle(X"+ k +", 10)", k, 0..n))[br]Execute(Zip("SetPointStyle(P"+ k +", 10)", k, 0..n))[br][br]Execute(Zip("S"+ k +" = Segment(P"+ k +", X"+ k +")", k, 0..n))[br]Execute(Zip("SetLineStyle(S"+ k +", 2)", k, 0..n))[br]Execute(Zip("SetLineThickness(S"+ k +", 2)", k, 0..n))[br][br]SL = Join({P0}, CellRange(P1, P50))[br][br]appx = Polyline(SL)[br][br]#-----------------------------------------------[br]Finally add in Update this script to Δt [br]#-----------------------------------------------[br]If( Δt <= 0.04, Execute(Zip("ShowLabel(X"+ k +", false)", k, 1..n)), Execute(Zip("ShowLabel(X"+ k +", true)", k, 1..n)))[br][br]If( Δt <= 0.04, Execute(Zip("ShowLabel(P"+ k +", false)", k, 1..n)), Execute(Zip("ShowLabel(P"+ k +", true)", k, 1..n)))