##################################### # Integrate a first order ODE using # 2nd order Runge-Kutta ##################################### # Define the phase speed f:= t-> cos(t): ## # Set initial condition and step length and number of steps ## y[0]:=1: h:= 0.1: N:= 100: ## # Advance solution using N Runge-Kutta steps ## t[0]:=0: for i from 1 to N do ### # Step from t[i-1] to t[i] using # the provisional step tt inbetween ### # 1. Compute f(t,y) f0:= f(t[i-1]); # 2. Set provisional value yt yt:= y[i-1] + h*f0; # 3. Compute phase speed ft at provisional point ft:= f(t[i-1]+h); # 4. Compute final phase speed ff ff:= (1/2)*f0 + (1/2)*ft; # 5. Set new solution value at t[i] y[i]:= y[i-1] + h*ff; # Update time t[i]:= t[i-1] + h; end do; ## # Build set of data points for plotting ## l:= [[t[n],y[n]] \$n=0..N]: ## # Plot approximate solution against exact solution ## with(plots): plotsetup(x11): #plot(l, style=line,symbol=circle); F:= plot(l, style=point,symbol=circle): G:= plot(1+sin(t),t=0..t[N],colour=blue): display({F,G});