Error in odearguments (line 90) f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0. Error in ode45 (line 115) odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin); Error in pendulum_sim (line 30) [t,f] = ode
35 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
I keep getting this error message whenever I run my code. I have attached my code for reference.
Appreciate some help on this :)
%% Pendulum - Muscle Model
% Create your own pendulum_muscle_equation function. It should be solvable
% with ode45 and the solution you return both theta and theta dot.
% PART 1: Study how co-contraction influences the stability and position/orientation of the forearm. To do this use three values of co-contraction u=0, u=0.5 and u=1.
% Plot the evolution of theta and thetadot for two different set of initial conditions, for which the system:
% 1) Does not fall to the ground for both small (u = 0.5) and large (u = 1) activation level;
% 2) Does not fall to the ground for large (u = 1) activation level and falls for smaller ones (u = 0.5);
figure()
ax1 = subplot(2,1,1);
ax2 = subplot(2,1,2);
%activation
u_set = [0,0;0.5,0.5;1,1];
tspan = [0 5];
n = size(u_set);
% Set the initial conditions
syms theta(t);
initial_cond = theta(0) == 0 ; %% INSERT YOUR CODE %%
options = odeset('Events',@ fallEventsFcn);
for i = 1:n(1)
u_set_i = u_set(i,:);
u1 = u_set_i(1);
u2 = u_set_i(2);
% Compute theta and thetadot using your own pendulum_muscle_equation
[t,f] = ode45(@(t,theta) pendulum_muscle_equation(t,theta,u1,u2),tspan,initial_cond,options);
subplot(ax1)
plot(t(:),f(:,1)*180/pi,'DisplayName',['u1 =' num2str(u_set_i(1)) '; u2 = ' num2str(u_set_i(2))]);
hold on;
subplot(ax2)
plot(t(:),f(:,2),'DisplayName',['u1 =' num2str(u_set_i(1)) '; u2 = ' num2str(u_set_i(2))]);
hold on;
end
ylabel(ax1,'Theta [°]');
ylabel(ax2,'Thetadot [rad/s]');
xlabel(ax2,'Time [s]');
title(ax1,'Angle');
title(ax2,'Velocity');
legend(ax1);
legend(ax2);
% Function to be called to solve ODE% Function to be called to solve ODE
function theta_dot = pendulum_muscle_equation(t, theta, u1,u2)
% Parameters
g = 9.81;
r = 0.3;
m = 1;
a1 = 0.0436;
a2 = 0.09;
alpha = 0.0218;
k_0 = 810.8;
k = 1621.6;
damping_ratio = 0.26;
delta = 2*damping_ratio*sqrt(k*m);
delta_0 = 2*damping_ratio*sqrt(k_0*m);
% use cosine rule to find l1 and l2
l1 = sqrt(a1^2+a2^2+2*a1*a2*sin(theta(1)));
l2 = sqrt(a1^2+a2^2-2*a1*a2*sin(theta(1)));
l_0 = sqrt(a1^2 + a2^2); % muscle inital length
% Find l_dot for both muscles
l1_dot = diff(sqrt(a1^2+a2^2+2*a1*a2*sin(theta(1))),t);
l2_dot = diff(sqrt(a1^2+a2^2-2*a1*a2*sin(theta(1))),t);
% use sine rule to find h1 and h2
h1 = a1*a2*cos(theta(1))/l1;
h2 = a1*a2*cos(theta(1))/l2;
% Find K_tot and D_tot for both muscles
K_tot1 = k_0 + k*u1;
D_tot1 = delta_0 + delta * sqrt(u1);
K_tot2 = k_0 + k*u2;
D_tot2 = delta_0 + delta * sqrt(u2);
% Muscle tensions
if l1 < l_0 - alpha * u1
mu1 = 0; % if muscle length < rest length -> tension = 0
else
mu1 = K_tot1 * (l_0 - alpha*u1 -l1) - D_tot1 * l1_dot;
% if muscle length > rest length -> tension follows muscle model eq
end
if l2 < l_0 - alpha * u2
mu2 = 0; % if muscle length < rest length -> tension = 0
else
mu2 = K_tot2 * (l_0 - alpha*u2 -l2) - D_tot2 * l2_dot;
% if muscle length > rest length -> tension follows muscle model eq
end
% Pendulum angular velocity and acceleration, to be integrated with ode45
theta_dot = zeros(2,1);
theta_dot(1) = theta(2);
theta_dot(2) = (m*g*r*sin(theta(1)) + mu1*h1 - mu2*h2)/(m*r^2);
% total torque about the elbow (load+muscle1+muscle2)
end
0 commentaires
Réponses (1)
Kiran Felix Robert
le 8 Fév 2021
Hi Ai,
The cause for this error is due to an incorrect call to the sine function at line 65.
The argument passed to the sine function is a symbolic equation and not a valid symbolic expression (refer sin function). This error can be resolved if the symbolic sine function is called with a valid argument.
0 commentaires
Voir également
Catégories
En savoir plus sur Ordinary Differential Equations dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!