Hello,
I am trying to implement ode45 to integrate my numerical function by following the ODE with Time-Dependent Terms in the ode45 documentation. I have a function with the form:
Where and are 1xn vectors that change at each time step and is a constant 1xn vector with each element being 9.81 . As suggested in the documentation, I created the vector for in the main script as follows: gt = linspace(0,0.1,25);
g = -tan(phi).*(zdd + gravity);
I am not sure about gt. I chose it arbitrarily.
Then, I defined my ode function as follows:
function dydt = myode(t,gt,g)
g = interp1(gt,g,t);
dydt = g;
Then, in the main script again, I tried to integrate using ode45 as follows:
tspan = [0 1.99];
ic = 0;
opts = odeset('RelTol',1e-2,'AbsTol',1e-4);
[t,y] = ode45(@(t) myode(t,gt,g), tspan, ic, opts);
However, if I try to run the script, I get an error saying that there are too many input arguments. Can someone tell me what I am doing wrong?
I understand that I can use cumtrapz() in my case to solve to do the numerical integration. However, my professor asked me to use ode45() since it offers a very accurate integration.
Moreover, I would like to point out that my main goal is to obtain . So, how can I use ode45() twice in this case?
If you are interested to find out how phi, zdd and gravity are defined, please refer to the code below:
I am solving an optimization problem using fmincon() to find the 9 parameters
z_hover1, z_start, z_end, z_hover2, phi_start, phi_end, t1, t2 ,t3
where t1,t2 and t3 are the times required to perform the trajectory in each of the three different phases of the flipping trajectory.
After finding the 9 parameters above, I can use them to create a flipping trajectory for a drone using a separate script. It should also be noted that I used the following function sevaral times in the script to come in order to find the trajectory of the position, velocity, acceleration, jerk and snap:
function traj = trajectory(start,goal,time)
global step;
[A,b] = build_matrix(start,goal,time);
coeff_p = fliplr((A\b)');
t = step:step:time;
coeff_v = polyder(coeff_p);
coeff_a = polyder(coeff_v);
coeff_j = polyder(coeff_a);
coeff_s = polyder(coeff_j);
traj = [polyval(coeff_p,t);
polyval(coeff_v,t);
polyval(coeff_a,t);
polyval(coeff_j,t);
polyval(coeff_s,t)];
end
So, the script which is then used to find the trajetory along y, z and phi is:
z1_start = [z_hover1 0 0 0 0];
z1_end = [z_start ((z_start-z_end)/t2 + g*t2/2) -g 0 0];
z1 = trajectory(z1_start,z1_end,t1);
phi1_start = [0 0 0 0 0];
phi1_end = [phi_start (phi_end-phi_start)/t2 0 0 0];
phi1 = trajectory(phi1_start,phi1_end,t1);
coeff_z2 = [-g/2 ((z_end-z_start)/t2+g*t2/2) z_start];
coeff_zd2 = polyder(coeff_z2);
coeff_zdd2 = polyder(coeff_zd2);
z2 = polyval(coeff_z2,step:step:t2);
zd2 = polyval(coeff_zd2,step:step:t2);
zdd2 = polyval(coeff_zdd2,step:step:t2);
coeff_phi2 = [(phi_end-phi_start)/t2 phi_start];
coeff_phid2 = polyder(coeff_phi2);
coeff_phidd2 = polyder(coeff_phid2);
phi2 = polyval(coeff_phi2,step:step:t2);
phid2 = polyval(coeff_phid2,step:step:t2);
phidd2 = polyval(coeff_phidd2,step:step:t2);
z3_start = [z_end ((z_start-z_end)/t2-g*t2/2) -g 0 0];
z3_end = [z_hover2 0 0 0 0];
z3 = trajectory(z3_start,z3_end,t3);
phi3_start = [phi_end (phi_end-phi_start)/t2 0 0 0];
phi3_end = [2*flips*pi 0 0 0 0];
phi3 = trajectory(phi3_start,phi3_end,t3);
z = [z1(1,:) z2 z3(1,:)];
zd = [z1(2,:) zd2 z3(2,:)];
zdd = [z1(3,:) zdd2 z3(3,:)];
phi = [phi1(1,:) phi2 phi3(1,:)];
phid = [phi1(2,:) phid2 phi3(2,:)];
phidd = [phi1(3,:) phidd2 phi3(3,:)];
gravity = g*ones(size(z));
ydd = -tan(phi).*(zdd + gravity);
As you can see, all that is left is to integrate ydd twice in order to the trajectory along y. And, this is where I am facing my problem.
Thank you for taking the time to read my question.