fmincon optimization not running properly, stops beacuase of StepTolerance
1 view (last 30 days)
Show older comments
I'm designing a strain wave gear for my bsc diploma, but I cannot run this optimization, and I ran out of people I can ask.
This is my code so far:
%Input paramaters:
beta = 62.31*pi/180;
r_0 = 24.519; %eq 3
e = 0.44;
n1 = 1;
n2 = 2;
n3 = 3;
%Input equalites and inequalites from the research:
ineq1_1 = (4*n1^2-1)*cos(2*n1*beta); %eq 8
ineq1_2 = (4*n2^2-1)*cos(2*n2*beta); %eq 8
ineq1_3 = (4*n3^2-1)*cos(2*n3*beta); %eq 8
ineq2_1 = 4*n1^2*(4*n1^2-1)*cos(n1*pi); %eq 9
ineq2_2 = 4*n2^2*(4*n2^2-1)*cos(n2*pi); %eq 9
ineq2_3 = 4*n3^2*(4*n3^2-1)*cos(n3*pi); %eq 9
lambda1 = (4*n1^2-1)*cos(2*n1*beta);
lambda2 = (4*n2^2-1)*cos(2*n2*beta);
lambda3 = (4*n3^2-1)*cos(2*n3*beta);
coef1 = (1/r_0^2)*lambda1; %eq 6
coef2 = (1/r_0^2)*lambda2; %eq 6
coef3 = (1/r_0^2)*lambda3; %eq 6
%Constarints:
x02 = [10,5,1];
lb = [];
ub = [];
nonlcon = [];
Aeq2 = [1,1,1];
beq2 = e;
A2 = [ineq1_1,ineq1_2,ineq1_3;ineq2_1,ineq2_2,ineq2_3];
b2 = [r_0,0];
fun2 = @(a)coef1*a(1)+coef2*a(2)+coef3*a(3);
options = optimoptions('fmincon','EnableFeasibilityMode',true,'Algorithm','sqp','StepTolerance',1e-2);
x2 = fmincon(fun2,x02,A2,b2,Aeq2,beq2,lb,ub,nonlcon,options);
It stops with this message:
Local minimum possible. Constraints satisfied.
fmincon stopped because the size of the current step is less than
the value of the step size tolerance and constraints are
satisfied to within the value of the constraint tolerance.
<stopping criteria details>
Changing the StepTolerance didn't fix it. It should give back 0.450; -0,005; -0,005.
Before this code, I wrote one with only two unknowns, where I calculated the ineqs, lambdas and coefs by hand, and I even forgot to convert degrees to radian and my calculator was in degress instead of radian, but this gave me the good results (0,0463; -0,023):
x02 = [100,100];
Aeq2 = [1,1];
beq2 = 0.44;
A2 = [-1.704,-5.317;11.98,238.55];
b2 = [24.519,0];
fun2 = @(a)-2.835*10^-3*a(1)-8.844*10^-3*a(2);
x2 = fmincon(fun2,x02,A2,b2,Aeq2,beq2);
But this one expanded to have 3 unknowns didn't work either.
I attached the research where the whole process is described and I worked from this.
0 Comments
Answers (3)
Walter Roberson
on 18 Apr 2023
fmincon is a general minimizer. It does not stop when the function value is zero: it is fine with going pretty negative.
As such, the only natural way for fmincon to be sure it has reached the best minimum possible would be if the function returned negative infinity.
In all other cases, fmincon has been coded with a number of different artificial ways to know when to give up. For example it knows to give up if the number of function evaluations reaches a configured maximum.
One of the artificial stopping criteria is a combination of circumstances: every time it changes any of the model parameters a little, the function value gets worse; and in trying to pinpoint the best possible spot near where it is, the size of the steps it is taking has reached a configured minimum.
In other words, as far as fmincon can tell, you are in a minimum (the gradient curls upwards in all directions) and maybe you could get a better position by adjusting the position half an iota of a smidgen, but you are darn close to the minimum near the current position.
This does not establish that you have reached a global minima. You cannot tell whether you have reached a global minima without a mathematical analysis of the function, which is not something that the technology can potentially do (you would need symbolic analysis to have a hope of that)
If you have not reached the model parameters that you expected, and those values are more than a trifle different, then adjusting the minimum step size is not going to help. You might have reached a local minima that is too deep for the algorithm to escape. Or you might have programmed the equations incorrectly.
Zaikun Zhang
on 18 Apr 2023
Edited: Zaikun Zhang
on 19 Apr 2023
I conducted a quick experiment using the function and data in the question. The solvers are fmincon of MATLAB and PRIMA, which is available at https://github.com/libprima/prima and https://www.mathworks.com/matlabcentral/fileexchange/127983-prima .
Here is the code I used (the defition of the problem is essnetially the same as that provided by @Levente Varga).
%Input paramaters:
beta = 62.31*pi/180;
r_0 = 24.519; %eq 3
e = 0.44;
n1 = 1;
n2 = 2;
n3 = 3;
%Input equalites and inequalites from the research:
ineq1_1 = (4*n1^2-1)*cos(2*n1*beta); %eq 8
ineq1_2 = (4*n2^2-1)*cos(2*n2*beta); %eq 8
ineq1_3 = (4*n3^2-1)*cos(2*n3*beta); %eq 8
ineq2_1 = 4*n1^2*(4*n1^2-1)*cos(n1*pi); %eq 9
ineq2_2 = 4*n2^2*(4*n2^2-1)*cos(n2*pi); %eq 9
ineq2_3 = 4*n3^2*(4*n3^2-1)*cos(n3*pi); %eq 9
lambda1 = (4*n1^2-1)*cos(2*n1*beta);
lambda2 = (4*n2^2-1)*cos(2*n2*beta);
lambda3 = (4*n3^2-1)*cos(2*n3*beta);
coef1 = (1/r_0^2)*lambda1; %eq 6
coef2 = (1/r_0^2)*lambda2; %eq 6
coef3 = (1/r_0^2)*lambda3; %eq 6
%Constarints:
x02 = [10,5,1]';
lb = [];
ub = [];
nonlcon = [];
Aeq2 = [1,1,1];
beq2 = e;
A2 = [ineq1_1,ineq1_2,ineq1_3;ineq2_1,ineq2_2,ineq2_3];
b2 = [r_0,0]'; % N.B.: Mathematically speaking, b2 should be a column.
fun2 = @(a)coef1*a(1)+coef2*a(2)+coef3*a(3);
fprintf('\nTesting fmincon:\n')
options = optimoptions('fmincon','EnableFeasibilityMode',true,'Algorithm','sqp','StepTolerance',1e-2);
x_fmincon = fmincon(fun2,x02,A2,b2,Aeq2,beq2,lb,ub,nonlcon,options);
fprintf('Solution provided by fmincon: x = [%g; %g; %g]\n', x_fmincon(1), x_fmincon(2), x_fmincon(3));
fprintf('f(x) = %.16f\n', fun2(x_fmincon));
fprintf('Constraint violation: %.16f\n\n', max([0; A2*x_fmincon-b2; abs(Aeq2*x_fmincon-beq2)]));
fprintf('Testing PRIMA:\n')
warning('off', 'prima:ReviseX0');
x_prima = prima(fun2,x02,A2,b2,Aeq2,beq2,lb,ub,nonlcon);
fprintf('Solution provided by PRIMA: x = [%g; %g; %g]\n', x_prima(1), x_prima(2), x_prima(3));
fprintf('f(x) = %.16f\n', fun2(x_prima));
fprintf('Constraint violation: %.16f\n\n', max([0; A2*x_prima-b2; abs(Aeq2*x_prima-beq2)]));
x = [0.450; -0.005; -0.005];
fprintf('Reference solution: x = [0.450; -0.005; -0.005]\n');
fprintf('f(x) = %.16f\n', fun2(x));
fprintf('Constraint violation: %.16f\n\n', max([0; A2*x-b2; abs(Aeq2*x-beq2)]));
Here is the result I obtained.
Testing fmincon:
Local minimum possible. Constraints satisfied.
fmincon stopped because the size of the current step is less than
the value of the step size tolerance and constraints are
satisfied to within the value of the constraint tolerance.
Solution provided by fmincon: x = [3.48821; -2.53259; -0.515619]
f(x) = -0.0166358927551636
Constraint violation: 0.0000000000000001
Testing PRIMA:
Solution provided by PRIMA: x = [8.81918; -6.96796; -1.41122]
f(x) = -0.0431461784076583
Constraint violation: 0.0000000000000005
Reference solution: x = [0.450; -0.005; -0.005]
f(x) = -0.0015141812204616
Constraint violation: 0.0000000000000000
Conclusions:
- Both fmincon and PRIMA achieve feasiblity up to the machine precision.
- Both fmincon and PRIMA provide feasible solutions much better than the reference point [0.450; -0.005; -0.005] in terms of the objective function value.
- PRIMA achives a much smaller function value than fmincon.
I hope this is helpful.
1 Comment
Walter Roberson
on 18 Apr 2023
format long g
%Input paramaters:
beta = 62.31*pi/180;
r_0 = 24.519; %eq 3
e = 0.44;
n1 = 1;
n2 = 2;
n3 = 3;
%Input equalites and inequalites from the research:
ineq1_1 = (4*n1^2-1)*cos(2*n1*beta); %eq 8
ineq1_2 = (4*n2^2-1)*cos(2*n2*beta); %eq 8
ineq1_3 = (4*n3^2-1)*cos(2*n3*beta); %eq 8
ineq2_1 = 4*n1^2*(4*n1^2-1)*cos(n1*pi); %eq 9
ineq2_2 = 4*n2^2*(4*n2^2-1)*cos(n2*pi); %eq 9
ineq2_3 = 4*n3^2*(4*n3^2-1)*cos(n3*pi); %eq 9
lambda1 = (4*n1^2-1)*cos(2*n1*beta);
lambda2 = (4*n2^2-1)*cos(2*n2*beta);
lambda3 = (4*n3^2-1)*cos(2*n3*beta);
coef1 = (1/r_0^2)*lambda1; %eq 6
coef2 = (1/r_0^2)*lambda2; %eq 6
coef3 = (1/r_0^2)*lambda3; %eq 6
%Constarints:
x02 = [10,5,1]';
lb = [];
ub = [];
nonlcon = [];
Aeq2 = [1,1,1];
beq2 = e;
A2 = [ineq1_1,ineq1_2,ineq1_3;ineq2_1,ineq2_2,ineq2_3];
b2 = [r_0,0]'; % N.B.: Mathematically speaking, b2 should be a column.
fun2 = @(a)coef1*a(1)+coef2*a(2)+coef3*a(3);
options = optimoptions('ga', 'MaxGenerations', 1e5);
[x_ga, f_ga] = ga(fun2, length(x02), A2, b2, Aeq2, beq2, lb, ub, nonlcon, options)
%double check
fun2(x_ga)
Unfortunately, most of the time this takes more than 55 second run-time limit for the Answers facility.
With coordiantes of [14144.5938655468 -11767.9316586277 -2376.22120927818] the result is -70.3386331774086
Matt J
on 18 Apr 2023
Edited: Matt J
on 18 Apr 2023
It should give back 0.450; -0,005; -0,005.
If you have that foreknowledge, you should apply lb and ub bounds as I have done below. Also, StepTolerance=1e-2 is way too large.
%Constarints:
x02 = [10,5,1];
lb = [-1,-1,-1];
ub = [1,1,1];
nonlcon = [];
Aeq2 = [1,1,1];
beq2 = e;
A2 = [ineq1_1,ineq1_2,ineq1_3;ineq2_1,ineq2_2,ineq2_3];
b2 = [r_0,0];
fun2 = @(a)coef1*a(1)+coef2*a(2)+coef3*a(3);
options = optimoptions('fmincon','EnableFeasibilityMode',true,...
'Algorithm','sqp','StepTolerance',1e-10);
[x2,fval,flag] = fmincon(fun2,x02,A2,b2,Aeq2,beq2,lb,ub,nonlcon,options)
These are not the x2 values you expect, but if we evaluate the objective fun2() at the expected point, we get a worse value than what the code above is giving you. So, either your expectations are wrong, or your problem definition is.
fun2([0.450; -0.005; -0.005])
ans =
-0.0015
2 Comments
Matt J
on 21 Apr 2023
The paper isn't going to be an effective way to communicate the optimization problem you are trying to solve to an unspecialized audience. It is advisable for you to write the problem in self-contained mathematical notation using the rich text capabilities of the Answers editor, e.g.,
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1362733/image.png)
See Also
Categories
Find more on Optimization 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!