Hi everyone,
I'm trying to make a position analysis of a mechanism via Newton-Raphson algorithm. When end of analysis ı can see graphics but ı see "Error nr" at command window. which part is wrong in my code?resim.JPG
%CLEAR: Variables and command window in MATLAB
clc,clear
% INPUT: PHYSICAL PARAMETERS of MECHANISM (meter)
a1=81;a2=20;a3=143;a4=67;a5=95;a6=25;a7=16;b=93;th2=18*pi/180
% INPUT: Maximum Iteration Number Nmax
Nmax=100;
% INPUT: INITIAL GUESS VALUES for th3, th4, th6, s respectively
x=[50*pi/180,40*pi/180,253*pi/180,5];
% INPUT: ERROR TOLERANCE
xe=0.000000001*abs(x);
% INPUT: SYSTEM INPUTS (th2,w2,al2)
dth=5*pi/180;
th2=-30*pi/180:dth:180*pi/180;
w2=-18*ones(1,length(th2));
al2=0*ones(1,length(th2));
%----------------------------------------------
xe=transpose(abs(xe));
kerr=1; %If kerr=1, results are not converged
for k=1:1:length(th2);
for n=1:Nmax
%----------------------------------------------
%Assign initial guess to unknowns
th3(k)=x(1);
th4(k)=x(2);
th6(k)=x(3);
s(k)=x(4);
% INPUT: JACOBIAN Matrix
J=zeros(4,4);
J(1,1)=-a3*sin(th3(k));J(1,2)=a4*sin(th4(k));
J(2,1)=a3*cos(th3(k));J(2,2)=-a4*cos(th4(k));
J(3,1)=-a5*sin(th3(k))-a6*sin(th3(k)+pi/2);J(3,3)=-a7*sin(th6(k)-pi/2)+s(k)*sin(th6(k));J(3,4)=-cos(th6(k));
J(4,1)=a5*cos(th3(k))+a6*cos(th3(k)+pi/2);J(4,3)=a7*cos(th6(k)-pi/2)-s(k)*cos(th6(k));J(4,4)=-sin(th6(k));
% INPUT: Function f
f=zeros(4,1);
f(1,1)=-(a2*cos(th2(k))+a3*cos(th3(k))-a4*cos(th4(k))-a1);
f(2,1)=-(a2*sin(th2(k))+a3*sin(th3(k))-a4*sin(th4(k)));
f(3,1)=-(a2*cos(th2(k))+a5*cos(th3(k))+a6*cos(th3(k)+pi/2)+a7*cos(th6(k)-pi/2)-s(k)*cos(th6(k)));
f(4,1)=-(a2*sin(th2(k))+a5*sin(th3(k))+a6*sin(th3(k)+pi/2)+a7*sin(th6(k)-pi/2)-s(k)*sin(th6(k))-b);
%----------------------------------------------
eps=inv(J)*f;x=x+transpose(eps);
if abs(eps)<xe
kerr=0;break
end
end
if kerr==1
'Error nr'
end
th3(k)=x(1);
th4(k)=x(2);
th6(k)=x(3);
s(k)=x(4);
%---velocity---------------------------
fv(1,1)=-w2(k)*a2.*sin(th2(k));
fv(2,1)=w2(k)*a2.*cos(th2(k));
fv(3,1)=-w2(k)*a2.*sin(th2(k));
fv(4,1)=-w2(k)*a2.*cos(th2(k));
vel=inv(J)*fv;
w3(k)=vel(1);
w4(k)=vel(2);
w6(k)=vel(3);
Vd(k)=vel(4);
%---acceleration---------------------------
fa(1,1)=-al2(k)*a2*sin(th2(k))-w2(k)^2*a2*cos(th2(k))-w3(k)^2*a3*cos(th3(k))+w4(k)^2*a4*cos(th4(k));
fa(2,1)=al2(k)*a2*cos(th2(k))-w2(k)^2*a2*sin(th2(k))-w3(k)^2*a3*sin(th3(k))+w4(k)^2*a4*sin(th4(k));
fa(3,1)=-al2(k)*a2*sin(th2(k))-w2(k)^2*a2*cos(th2(k))-a5*w3(k)^2*cos(th3(k))-a6*w3(k)^2*cos(th2(k)+pi/2)+2*Vd(k)*w6(k)*sin(th6(k))+s(k)*w6(k)^2*cos(th6(k));
fa(4,1)=al2(k)*a2*cos(th2(k))-a2*w2(k)^2*sin(th2(k))-a5*w3(k)^2*sin(th3(k))-a6*w3(k)^2*sin(th2(k)+pi/2)+a7*w6(k)^2*sin(th6(k)-pi/2)-2*Vd(k)*w6(k)*cos(th6(k))-s(k)*w6(k)^2*sin(th6(k));
acc=inv(J)*fa;
al3(k)=acc(1);
al4(k)=acc(2);
al6(k)=acc(3);
ad(k)=acc(4);
end
% Angle: radian --> degree
th2d=th2*180/pi;
th3d=th3*180/pi;
th4d=th4*180/pi;
th6d=th6*180/pi;
%--------Plots---------------
figure(1),
subplot(4,3,1),plot(th2d,th3d,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\theta_3(^o)'),xlim([0 180])
subplot(4,3,2),plot(th2d,w3,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\omega_3(r/s)'),xlim([0 180])
subplot(4,3,3),plot(th2d,al3,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\alpha_3(r/s^2)'),xlim([0 180])
subplot(4,3,4),plot(th2d,th4d,'r','linewidth',2),xlabel('\theta_2(^o)'),ylabel('\theta_4 (^o)'),xlim([0 180])
subplot(4,3,5),plot(th2d,w4,'r','linewidth',2),xlabel('\theta_2(^o)'),ylabel('\omega_4 (r/s)'),xlim([0 180])
subplot(4,3,6),plot(th2d,al4,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\alpha_4(r/s^2)'),xlim([0 180])
subplot(4,3,7),plot(th2d,th6,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\theta_6(^o)'),xlim([0 180])
subplot(4,3,8),plot(th2d,w6,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\omega_6(r/s)'),xlim([0 180])
subplot(4,3,9),plot(th2d,al6,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('\alpha_6(r/s^2)'),xlim([0 180])
subplot(4,3,10),plot(th2d,s,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('s(cm)'),xlim([0 180])
subplot(4,3,11),plot(th2d,Vd,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('V_d(cm/s)'),xlim([0 180])
subplot(4,3,12),plot(th2d,ad,'r','linewidth',2),xlabel('\theta_2 (^o)'),ylabel('a_d(cm/s^2)'),xlim([0 180])

Réponses (1)

KALYAN ACHARJYA
KALYAN ACHARJYA le 17 Août 2019
Modifié(e) : KALYAN ACHARJYA le 17 Août 2019

0 votes

Its not a Matlab Error.Just display the string as 'Error nr', when ker = 1
Diaplay from this part of the code
if kerr==1
'Error nr'
end
if you delete that part of the code, still, there will be no issue on results. Here it just used as flag to inform when ker=1;

Catégories

En savoir plus sur Mathematics dans Centre d'aide et File Exchange

Produits

Version

R2016b

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by