Quote:
Originally Posted by Ether
It's a cost/benefit thing. The cost of converting Euler to trapezoidal for this problem is very small. Going from trap to Simpson's would require more rework, make the code less readable, and probably even slow it down for equivalent accuracy.
|
True (except for the last part -- I don't think Simpson's integration would slow it down for equivalent accuracy, but I haven't tried it, so I don't know).
Quote:
|
It would be interesting to do a comparison. My guess is that for this problem, trapezoidal in compiled C would be just as fast for equivalent accuracy as RK in MatLab
|
Having used MATLAB a bit, I think compiled C++ w/ trapezoidal integration would be way faster for the same accuracy than MATLAB code using an RK-based solver.
Here's code that solves it using one of MATLAB's adaptive ODE solvers (4/5th-order Runge-Kutta (I don't recall which formulation it is)):
Spoiler for more code:
Code:
% Set up various useful functions that don't require integration to calculate
angle = @(t) 5*pi/12 - 3.75 + 3.75 * cos(t/2.5); % Obtained analytically (by hand)
forward = @(t) 5*sin(t/2);
strafe_right = @(t) 4*sin(t/2.2);
xvel = @(t) cos(angle(t)) * forward(t) + cos(angle(t) - pi/2) * strafe_right(t);
yvel = @(t) sin(angle(t)) * forward(t) + sin(angle(t) - pi/2) * strafe_right(t);
spd = @(t) sqrt(forward(t)^2 + strafe_right(t)^2);
% Set up for the ODE solution
tspan = [0 30];
x0 = [0; 0; 0]; % xpos, ypos, distance
% The ODE function
dynamics = @(t,x) [xvel(t); yvel(t); spd(t)];
% Call ode45
diffsoln = ode45(dynamics, tspan, x0);
% Function for obtaining the state at a given time (for plotting)
statefun = @(t) deval(diffsoln, t);
% Time values for plot (using a lot of points so it looks good).
times = linspace(tspan(1), tspan(2), 10000);
% X and Y positions
xpos = [1 0 0] * statefun(times);
ypos = [0 1 0] * statefun(times);
% Do the plotting
plot(xpos, ypos, 'LineWidth', 5)