Newton Raphson gives answers divided by 0
4 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
I'm trying to estimate 4 angles and 2 distances using the newton raphson method with phi4 and theta as the input. I have 6 functions that equal to 0 and when i try substituting to calculate dx (the array of all my values) i get a division by 0. I transposed some arrays to columns, but id be surprised if thats the issue. The functions themselves check out and im mostly using matlab functions to make this work so I don't really know what could be wrong
Angle=3;
phi4=pi/180*Angle;
theta=linspace(0,2*pi, 72);
%First guesses
phi1=ones(1, 72);
phi2=ones(1, 72);
phi3=ones(1, 72);
phi5=ones(1, 72);
d=ones(1,72);
s=ones(1, 72);
syms phi1V phi2V phi3V phi5V dV sV;
xV=[sV, dV, phi1V, phi2V, phi3V, phi5V];
ds=0;
dd=0;
dphi1=0;
dphi2=0;
dphi3=0;
dphi4=0;
dx=[ds, dd, dphi1, dphi2, dphi3, phi5]';
lim=0.01;
i=1;
while(i<73)
%redefine functions for changes in theta
f1=25-sV+dV*cos(phi2V)+3.1*cos(phi3V);
f2=-5+dV*sin(phi2V)+3.1*sin(phi3V);
f3=-25+6.3*cos(theta(i))+27*cos(phi1V)-19.2*cos(phi2V);
f4=11+6.3*sin(theta(i))+27*sin(phi1V)-19.2*sin(phi2V);
f5=-40+6.3*cos(theta(i))+(dV-19.2)*cos(phi2V)+29.1*cos(phi5V)+16*cos(phi4);
f6=-24+6.3*sin(theta(i))+(dV-19.2)*sin(phi2V)+29.1*sin(phi5V)+16*sin(phi4);
f=[f1, f2, f3, f4, f5, f6]';
J=jacobian(f, xV);
M=-inv(J)*f;
x=[s(i), d(i), phi1(i), phi2(i), phi3(i), phi5(i)];
firstloop=1;
while(ds > lim | dd > lim | dphi1 > lim | dphi2 > lim | dphi3 > lim | dphi4 > lim | firstloop==1)
dx=subs(M, xV, x);
x=x+dx;
firstloop=0;
end
i=i+1;
end
0 commentaires
Réponse acceptée
Sam Chak
le 10 Mar 2024
Hi @Jean-Paul
During the computation of M, certain elements involve divisions by . In the line 'dx = subs(M, xV, x)', an attempt is made to substitute the variables in 'xV = [sV, dV, phi1V, phi2V, phi3V, phi5V]' with 'x = [1, 1, 1, 1, 1, 1]', which leads to . Consequently, it triggers the 'Division-by-zero' error message.
0 commentaires
Plus de réponses (2)
Torsten
le 10 Mar 2024
Modifié(e) : Torsten
le 10 Mar 2024
Your linear systems in the course of the Newton iterations cannot be solved since the Jacobian is rank-deficient (it has only rank 5 instead of rank 6 as needed). This usually indicates that the equations you try to solve are not independent or even contradictory.
Angle=3;
phi4=pi/180*Angle;
theta=linspace(0,2*pi, 72);
%First guesses
phi1=ones(1, 72);
phi2=ones(1, 72);
phi3=ones(1, 72);
phi5=ones(1, 72);
d=ones(1,72);
s=ones(1, 72);
syms phi1V phi2V phi3V phi5V dV sV;
xV=[sV, dV, phi1V, phi2V, phi3V, phi5V];
ds=0;
dd=0;
dphi1=0;
dphi2=0;
dphi3=0;
dphi4=0;
%dx=[ds, dd, dphi1, dphi2, dphi3, phi5]';
lim=0.01;
i=1;
while(i<73)
%redefine functions for changes in theta
f1=25-sV+dV*cos(phi2V)+3.1*cos(phi3V);
f2=-5+dV*sin(phi2V)+3.1*sin(phi3V);
f3=-25+6.3*cos(theta(i))+27*cos(phi1V)-19.2*cos(phi2V);
f4=11+6.3*sin(theta(i))+27*sin(phi1V)-19.2*sin(phi2V);
f5=-40+6.3*cos(theta(i))+(dV-19.2)*cos(phi2V)+29.1*cos(phi5V)+16*cos(phi4);
f6=-24+6.3*sin(theta(i))+(dV-19.2)*sin(phi2V)+29.1*sin(phi5V)+16*sin(phi4);
f=[f1, f2, f3, f4, f5, f6]';
J=jacobian(f, xV);
%M=-inv(J)*f
x=[s(i), d(i), phi1(i), phi2(i), phi3(i), phi5(i)]
firstloop=1;
while(ds > lim | dd > lim | dphi1 > lim | dphi2 > lim | dphi3 > lim | dphi4 > lim | firstloop==1)
dx=-subs(J,xV,x)\subs(f,xV,x);%subs(M, xV, x)
rank(subs(J,xV,x))
rank([subs(J,xV,x),subs(f,xV,x)])
x=x+dx;
firstloop=0;
end
i=i+1;
end
0 commentaires
Matt J
le 11 Mar 2024
Modifié(e) : Matt J
le 11 Mar 2024
Raw Newton-Raphson is a pretty terrible method, but I assume you'll tell us you have no choice...
[f,M]=preanalysis(); %Symbolic pre-analysis
%%%%Numeric stuff -- Newton-Raphson
Angle=3;
phi4=pi/180*Angle;
theta=linspace(0,2*pi, 72);
lim=1e-6;
i=1;
clear solutions
while(i<73)
F=@(x) f(x, theta(i),phi4);
Update=@(x) M(x, theta(i),phi4);
x=F(zeros(6,1))/1000; %Initial guess - needs some experimentation
dx=inf;
firstloop=1;
while( all(abs(dx)>lim) | firstloop==1)
dx=Update(x);
if any(~isfinite(dx)),
warning 'Bad update. Aborting...'
break;
end
x=x+dx;
x(1:4)=atan2( sin(x(1:4)) , cos(x(1:4))); %normalize the 2*pi periodicity
firstloop=0;
end
solutions{i}=x;
i=i+1;
end
solutions=cell2mat(solutions);
solutions(:,1:20)
function [F,M]=preanalysis()
%%%Symbolic stuff
syms phi1 phi2 phi3 phi5 d s theta phi4 real
f1=25-s+d*cos(phi2)+3.1*cos(phi3);
f2=-5+d*sin(phi2)+3.1*sin(phi3);
f3=-25+6.3*cos(theta)+27*cos(phi1)-19.2*cos(phi2);
f4=11+6.3*sin(theta)+27*sin(phi1)-19.2*sin(phi2);
f5=-40+6.3*cos(theta)+(d-19.2)*cos(phi2)+29.1*cos(phi5)+16*cos(phi4);
f6=-24+6.3*sin(theta)+(d-19.2)*sin(phi2)+29.1*sin(phi5)+16*sin(phi4);
f=[f1, f2, f3, f4, f5, f6]';
x=[phi1 phi2 phi3 phi5 d s]';
J=jacobian(f, x);
M=matlabFunction(-J\f,'Vars',{x,theta,phi4});
F=matlabFunction(f,'Vars',{x,theta,phi4});
end
0 commentaires
Voir également
Catégories
En savoir plus sur Symbolic Math Toolbox 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!