Need more accurate answer from fmincon/fminsearch

Hi!

I am using fmincon function to find the reflection point on a spherical reflector. I set a source point s at position (6,6,0), field point f at (-6,6,0) and the spherical reflector to get the reflection point r at (0,0,0). fmincon return the value of r that is very small, but not equal to zero. As I increase the radius of the reflector to infinitity, the inaccuracy becomes more obvious and will affect following calculations. I have try to use fminsearch, but still encounter the same problem.

I found out that it happened because fmincon and fminsearch do not return the answer of phi and theta as accurate as I need (pi/2). Is there any way to make the answer of fminsearch more accurate?

rad = 100;      %set the radius of the spherical reflector
q = [0,-100,0];  %set the center of the spherical reflector
s = [6,6,0];  %set the signal source position
f = [-6,6,0]; %set the field point position

%calculate the minimum distance, fval between s and f, through point of reflection, r
L =@(w) sqrt((q(1) + rad.*sin(w(1)).*cos(w(2)) - s(1))^2 + (q(2) + rad.*sin(w(1)).*sin(w(2)) - s(2)).^2 + (q(3) + rad.*cos(w(1)) - s(3)).^2) + sqrt((f(1) - rad.*sin(w(1)).*cos(w(2)) - q(1)).^2 + (f(2) - rad.*sin(w(1)).*sin(w(2)) - q(2)).^2 + (f(3) - rad.*cos(w(1)) - q(1)).^2);
w0= [0 0];
[w,fval] = fmincon(L,w0,[],[],[],[],[0 0],[pi 2*pi]) %w1=theta, w2=phi

%determine point of reflection, r
x = q(1) + rad*sin(w(1))*cos(w(2));
y = q(2) + rad*sin(w(1))*sin(w(2));
z = q(3) + rad*cos(w(1));
r = [x,y,z]

Answer 1

[Walter Roberson]

If you proceed numerically and you push for higher accuracy, then there is a risk that you might head towards 1.5707963263235436013331991489394567906856536865234375 1.5707963261704962487641523694037459790706634521484375 which gives a slightly smaller fval than your L calculates for [pi/2, pi/2] . The theoretical fval at [pi/2, pi/2] is 12*2^(1/2) but because of the finite precision of numeric calculations, it is possible for rounding differences to lower the fval to less than the theoretical.

Slightly improved code:

rad = 100;      %set the radius of the spherical reflector
q = [0,-100,0];  %set the center of the spherical reflector
s = [6,6,0];  %set the signal source position
f = [-6,6,0]; %set the field point position

%calculate the minimum distance, fval between s and f, through point of reflection, r
L =@(w) sqrt((q(1) + rad.*sin(w(1)).*cos(w(2)) - s(1))^2 + (q(2) + rad.*sin(w(1)).*sin(w(2)) - s(2)).^2 + (q(3) + rad.*cos(w(1)) - s(3)).^2) + sqrt((f(1) - rad.*sin(w(1)).*cos(w(2)) - q(1)).^2 + (f(2) - rad.*sin(w(1)).*sin(w(2)) - q(2)).^2 + (f(3) - rad.*cos(w(1)) - q(1)).^2);
%if you use [0 0] as the bound then interior point will notice it is on the
%lower bound and will move it to [0.99 0.99]
w0= [2 2];

%result gets worse but unevenly if you reduce FiniteDifferenceStepSize
opt = optimoptions('fmincon', 'StepTolerance', 1e-14, 'Display', 'iter', 'FiniteDifferenceStepSize', 1e-9);
A = []; b = [];
Aeq = []; beq = [];
lb = [0 0]; ub = [pi 2*pi];
nonlcon = [];
[w, fval, exitflag, output] = fmincon(L,w0, A, b, Aeq, beq, lb, ub, nonlcon, opt); %w1=theta, w2=phi

%determine point of reflection, r
x = q(1) + rad*sin(w(1))*cos(w(2));
y = q(2) + rad*sin(w(1))*sin(w(2));
z = q(3) + rad*cos(w(1));
r = [x,y,z]

[Ashraf Ramzan]

Sorry for the late reply. I've tried the improved code. It does give a slightly smaller fval. Thank you very much.