Why results of the 'solve' are wrong?

6 views (last 30 days)
p_max=715.685424949238;
p_min=415.685424949238;
z_max=128.051406871193;
z_min=43.1985931288072;
D1_max=1128.05140687119;
D1_min=1043.19859312881;
angle=45;
y0=0:1:200;
p0=y0/sind(angle);
x0=p0*cosd(angle);
syms a b
q=sqrt((a-p_min)^2+(b-z_min)^2)-D1_min-sqrt((a-p_max)^2+(b-z_max)^2)+D1_max==0; %equation to solve
z_solve=solve(q,b);
z_value=subs(z_solve,a,p0);
z0=eval(z_value);
q1=sqrt((p0-p_min).^2+(z0-z_min).^2)-D1_min-sqrt((p0-p_max).^2+(z0-z_max).^2)+D1_max; %check the reults
  2 Comments
Steven Lord
Steven Lord on 26 Jan 2023
You have not defined a variable named angle so when I run this code it tries to call the angle function with no input arguments and that throws an error. I recommend choosing a different name for your variable and editing your question to include the value of that variable.
ruiheng wu
ruiheng wu on 26 Jan 2023
sorry about that,'angle' is 45.

Sign in to comment.

Accepted Answer

Walter Roberson
Walter Roberson on 26 Jan 2023
p_max=715.685424949238;
(and a bunch of other numbers).
What does that input mean? Does it mean that p_max is exactly ? Does it mean that p_max is some number not precisely defined that is between 715.6854249492375 (inclusive) and 715.6854249492385 (exclusive) ? Or is it intentional that p_max is exactly
fprintf('%.999g\n', p_max)
715.6854249492380404262803494930267333984375
which is how MATLAB will interpret it?
q=sqrt((a-p_min)^2+(b-z_min)^2)-D1_min-sqrt((a-p_max)^2+(b-z_max)^2)+D1_max==0; %equation to solve
You have sqrt() -- the Euclidean distance between (a,b) and (p_min, z_min) . most of the time, sqrt() is going to lead to irrational numbers. So your equation is going to try to exactly balance irrational numbers
z_solve=solve(q,b);
Remember that the purpose of solve() is to find exact closed form algebraic solutions . So the solution it comes up with is not necesarily going to be remotely valid if p_max is exactly or if it is intended to be between 715.6854249492375 (inclusive) and 715.6854249492385 (exclusive). If you did not intend p_max to be exactly that number ending in "375" then you either formulated the problem incorrectly or else you made a mistake in using solve()
solve() is very much a "Garbage In, Garbage Out" solver. It is happy to give you the exact solution to the problem you actually asked it to solve, even though that might not be the problem you intended to solve. For any non-linear problem, do not assume that just because your inputs are "close" to the real inputs, that solve() is going to give you a solution that is "close" to the "real" solution: the exact solution can be a very fine balance between components of the expression.
  6 Comments
ruiheng wu
ruiheng wu on 29 Jan 2023
I made a mistake,'f' should have been 'F'.
Walter Roberson
Walter Roberson on 29 Jan 2023
Q = @(v) sym(v);
angle = Q(45); %The results are correct when the angles are 30 degrees and 90 degrees.
F = Q(1000);
y_max = Q(550);
y_min = Q(250);
R=(y_max-y_min)/2;
y_avg=(y_max+y_min)/2;
x_avg=y_avg*cosd(angle)/sind(angle);
y0 = Q(0:1:200);
d=sqrt(x_avg^2+y_avg^2);
p_max=d+R;
p_min=d-R;
z_max=p_max^2/(4*F);
z_min=p_min^2/(4*F);
D1_min=sqrt(p_min.^2+(F-z_min).^2);
D1_max=sqrt(p_max.^2+(F-z_max).^2);
p0=y0/sind(angle);
x0=p0*cosd(angle);
syms a b
q=sqrt((a-p_min)^2+(b-z_min)^2)-D1_min-sqrt((a-p_max)^2+(b-z_max)^2)+D1_max==0; %equation to solve
z_solve = solve(q,b)
z_solve = 
z_value=subs(z_solve,a,p0);
z0 = subs(z_value);
q1=sqrt((p0-p_min).^2+(z0-z_min).^2)-D1_min-sqrt((p0-p_max).^2+(z0-z_max).^2)+D1_max; %check the reults
simplify(q1(1:5))
ans = 

Sign in to comment.

More Answers (1)

Paul
Paul on 26 Jan 2023
Hi ruiheng,
I would stay out of the double world to first make sure everything is working as expected.
p_max=715.685424949238;
p_min=415.685424949238;
z_max=128.051406871193;
z_min=43.1985931288072;
D1_max=1128.05140687119;
D1_min=1043.19859312881;
angle=45;
y0=0:1:200;
p0=y0/sind(angle);
x0=p0*cosd(angle);
syms a b
q = sqrt((a-p_min)^2+(b-z_min)^2)-D1_min-sqrt((a-p_max)^2+(b-z_max)^2)+D1_max == 0; %equation to solve
Add ReturnConditions
z_solve = solve(q,b,'ReturnConditions',true);
vpa(z_solve.b,5)
ans = 
vpa(z_solve.conditions,5)
ans = 
%z_value = subs(z_solve,a,p0);
%z0 = eval(z_value)
%z0 = z_value(:,5);
%q1 = sqrt((p0-p_min).^2+(z0-z_min).^2)-D1_min-sqrt((p0-p_max).^2+(z0-z_max).^2)+D1_max %check the reults
%simplify(q1)
Only checking a few solutions, timed-out when trying to check all
for ii = 8:12, %1:numel(p0)
% verify that the value of a meets the condition for the solution
cond(ii) = simplify(subs(z_solve.conditions,a,p0(ii)));
% check the value of the solution for the value of a
check(ii) = simplify(subs(lhs(q),[a b],{p0(ii) subs(z_solve.b,a,p0(ii))}));
end
cond(8:12)
ans = 
check(8:12)
ans = 
  1 Comment
ruiheng wu
ruiheng wu on 27 Jan 2023
Thank you for checking the code.This should be the cause of the problem.

Sign in to comment.

Products


Release

R2020a

Community Treasure Hunt

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

Start Hunting!