Newton Raphson gives answers divided by 0
5 views (last 30 days)
Show older comments
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 Comments
Accepted Answer
Sam Chak
on 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 Comments
More Answers (2)
Torsten
on 10 Mar 2024
Edited: Torsten
on 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 Comments
Matt J
on 11 Mar 2024
Edited: Matt J
on 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 Comments
See Also
Categories
Find more on Symbolic Math Toolbox in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!