summaryrefslogtreecommitdiff
path: root/macros/fmincon.sci
diff options
context:
space:
mode:
Diffstat (limited to 'macros/fmincon.sci')
-rw-r--r--macros/fmincon.sci404
1 files changed, 228 insertions, 176 deletions
diff --git a/macros/fmincon.sci b/macros/fmincon.sci
index 2393649..9faefc4 100644
--- a/macros/fmincon.sci
+++ b/macros/fmincon.sci
@@ -93,6 +93,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
// <listitem>output.Cpu_Time: The total cpu-time spend during the search</listitem>
// <listitem>output.Objective_Evaluation: The number of Objective Evaluations performed during the search</listitem>
// <listitem>output.Dual_Infeasibility: The Dual Infeasiblity of the final soution</listitem>
+ // <listitem>output.Message: The output message for the problem</listitem>
// </itemizedlist>
//
// The lambda data structure contains the Lagrange multipliers at the end
@@ -260,6 +261,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
// options=list("MaxIter", [1500], "CpuTime", [500], "GradObj", fGrad, "Hessian", lHess,"GradCon", cGrad);
// //Calling Ipopt
// [x,fval,exitflag,output,lambda,grad,hessian] =fmincon(f, x0,A,b,Aeq,beq,lb,ub,nlc,options)
+ // // Press ENTER to continue
// Authors
// R.Vidyadhar , Vignesh Kannan
@@ -268,10 +270,15 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
[lhs , rhs] = argn();
//To check the number of arguments given by the user
- if ( rhs<4 | rhs==5 | rhs==7 | rhs>10 ) then
+ if ( rhs<4 | rhs>10 ) then
errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while it should be 4,6,8,9,10"), "fmincon", rhs);
error(errmsg)
end
+
+ if (rhs==5 | rhs==7) then
+ errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while it should be 4,6,8,9,10s"), "fmincon", rhs);
+ error(errmsg)
+ end
//Storing the Input Parameters
fun = varargin(1);
@@ -299,16 +306,10 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
end
//To check whether the 1st Input argument (fun) is a function or not
- if (type(fun) ~= 13 & type(fun) ~= 11) then
- errmsg = msprintf(gettext("%s: Expected function for Objective (1st Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", fun, "f", 1, "function");
//To check whether the 2nd Input argument (x0) is a vector/scalar
- if (type(x0) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Vector/Scalar for Starting Point (2nd Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", x0, "x0", 2, "constant");
//To check and convert the 2nd Input argument (x0) to a row vector
if((size(x0,1)~=1) & (size(x0,2)~=1)) then
@@ -346,10 +347,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
endfunction
//To check whether the 3rd Input argument (A) is a Matrix/Vector
- if (type(A) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Matrix/Vector for Constraint Matrix A (3rd parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", A, "A", 3, "constant");
//To check for correct size of A(3rd paramter)
if(size(A,2)~=s(2) & size(A,2)~=0) then
@@ -360,10 +358,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
s1=size(A);
//To check whether the 4th Input argument (b) is a vector/scalar
- if (type(b) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Vector/Scalar for b (4th Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", b, "b", 4, "constant");
//To check for the correct size of b (4th paramter) and convert it into a column vector which is required for Ipopt
if(s1(2)==0) then
@@ -389,10 +384,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
end
//To check whether the 5th Input argument (Aeq) is a matrix/vector
- if (type(Aeq) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Matrix/Vector for Equality Constraint Matrix Aeq (5th Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", Aeq, "Aeq", 5, "constant");
//To check for the correct size of Aeq (5th paramter)
if(size(Aeq,2)~=s(2) & size(Aeq,2)~=0) then
@@ -403,10 +395,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
s2=size(Aeq);
//To check whether the 6th Input argument(beq) is a vector/scalar
- if (type(beq) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Vector/Scalar for beq (6th Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", beq, "beq", 6, "constant");
//To check for the correct size of beq(6th paramter) and convert it into a column vector which is required for Ipopt
if(s2(2)==0) then
@@ -433,10 +422,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
//To check whether the 7th Input argument (lb) is a vector/scalar
- if (type(lb) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Vector/Scalar for Lower Bound Vector (7th Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", lb, "lb", 7, "constant");
//To check for the correct size and data of lb (7th paramter) and convert it into a column vector as required by Ipopt
if (size(lb,2)==0) then
@@ -459,10 +445,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
end
//To check whether the 8th Input argument (ub) is a vector/scalar
- if (type(ub) ~= 1) then
- errmsg = msprintf(gettext("%s: Expected Vector/Scalar for Upper Bound Vector (8th Parameter)"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", ub, "ub", 8, "constant");
//To check for the correct size and data of ub (8th paramter) and convert it into a column vector as required by Ipopt
if (size(ub,2)==0) then
@@ -497,7 +480,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
end
if(ub(i)-lb(i)<=1e-6) then
- errmsg = msprintf(gettext("%s: Difference between Upper Bound and Lower bound should be atleast > 10^6 for variable number= %d "), "fmincon", i);
+ errmsg = msprintf(gettext("%s: Difference between Upper Bound and Lower bound should be atleast > 10^-6 for variable number= %d "), "fmincon", i);
error(errmsg)
end
end
@@ -530,7 +513,6 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
no_nlic = size(sample_c,2);
no_nlec = size(sample_ceq,2);
no_nlc = no_nlic + no_nlec;
-
//Constructing a single output variable function for nlc
function y = addnlc(x)
[c,ceq] = nlc(x);
@@ -580,74 +562,13 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
end
//If options has been entered, then check its type for 'list'
- if (type(param) ~= 15) then
- errmsg = msprintf(gettext("%s: Options (10th parameter) should be a list"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", param, "options", 10, "list");
//If options has been entered, then check whether an even number of entires has been entered
if (modulo(size(param),2)) then
errmsg = msprintf(gettext("%s: Size of Options (list) should be even"), "fmincon");
error(errmsg);
end
-
-
- //Defining a function to calculate Gradient or Hessian if the respective user entry is OFF
- function [y,check] = gradhess(x,t)
- if t==1 then //To return Gradient
- if(execstr('y=numderivative(fun,x)','errcatch')==10000)
- y=0;
- check=1;
- else
- y=numderivative(fun,x);
- if (isreal(y)==%F) then
- y=0;
- check=1;
- else
- check=0;
- end
- end
- elseif t==2 then //To return Hessian
- if(execstr('[grad,y]=numderivative(fun,x)','errcatch')==10000)
- y=0;
- check=1;
- else
- [grad,y]=numderivative(fun,x);
- if (isreal(y)==%F) then
- y=0;
- check=1;
- else
- check=0;
- end
- end
- elseif t==3 then //To return Gradient
- if(execstr('y=numderivative(addnlc,x)','errcatch')==10000)
- y=0;
- check=1;
- else
- y=numderivative(addnlc,x);
- if (isreal(y)==%F) then
- y=0;
- check=1;
- else
- check=0;
- end
- end
- elseif t==4 then //To return Hessian
- if(execstr('[grad,y]=numderivative(addnlc,x)','errcatch')==10000)
- y=0;
- check=1;
- else
- [grad,y]=numderivative(addnlc,x);
- if (isreal(y)==%F) then
- y=0;
- check=1;
- else
- check=0;
- end
- end
- end
- endfunction
//To set default values for options, if user doesn't enter options
options = list("MaxIter", [3000], "CpuTime", [600]);
@@ -669,31 +590,65 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
//To check the user entry for options and storing it
for i = 1:(size(param))/2
select convstr(param(2*i-1),'l')
- case "maxiter" then
- options(2*i) = param(2*i); //Setting the maximum number of iterations as per user entry
+ case "maxiter" then
+ if (type(param(2*i))~=1) then
+ errmsg = msprintf(gettext("%s: Value for Maximum Iteration should be a Constant"), "fmincon");
+ error(errmsg);
+ else
+ options(2*i) = param(2*i); //Setting the maximum number of iterations as per user entry
+ end
case "cputime" then
- options(2*i) = param(2*i); //Setting the maximum CPU time as per user entry
+ if (type(param(2*i))~=1) then
+ errmsg = msprintf(gettext("%s: Value for Maximum Cpu-time should be a Constant"), "fmincon");
+ error(errmsg);
+ else
+ options(2*i) = param(2*i); //Setting the maximum CPU time as per user entry
+ end
case "gradobj" then
- flag1=1;
- fGrad=param(2*i);
+ if (type(param(2*i))==10) then
+ if (convstr(param(2*i))=="off") then
+ flag1 =0;
+ else
+ errmsg = msprintf(gettext("%s: Unrecognized String [%s] entered for the option- %s."), "fmincon",param(2*i), param(2*i-1));
+ error(errmsg);
+ end
+ else
+ flag1 = 1;
+ fGrad = param(2*i);
+ end
case "hessian" then
- flag2=1;
- lHess=param(2*i);
+ if (type(param(2*i))==10) then
+ if (convstr(param(2*i))=="off") then
+ flag2 =0;
+ else
+ errmsg = msprintf(gettext("%s: Unrecognized String [%s] entered for the option- %s."), "fmincon",param(2*i), param(2*i-1));
+ error(errmsg);
+ end
+ else
+ flag2 = 1;
+ lHess = param(2*i);
+ end
case "gradcon" then
- flag3=1;
- cGrad=param(2*i);
+ if (type(param(2*i))==10) then
+ if (convstr(param(2*i))=="off") then
+ flag3 =0;
+ else
+ errmsg = msprintf(gettext("%s: Unrecognized String [%s] entered for the option- %s."), "fmincon",param(2*i), param(2*i-1));
+ error(errmsg);
+ end
+ else
+ flag3 = 1;
+ cGrad = param(2*i);
+ end
else
errmsg = msprintf(gettext("%s: Unrecognized parameter name %s."), "fmincon", param(2*i-1));
error(errmsg);
end
end
- //To check for correct input of Gradient and Hessian functions from the user
+ //To check for correct input of Objective Gradient function from the user
if (flag1==1) then
- if (type(fGrad) ~= 11 & type(fGrad) ~= 13) then
- errmsg = msprintf(gettext("%s: Expected function for Gradient of Objective"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", fGrad, "fGrad", 10, "function");
if(execstr('sample_fGrad=fGrad(x0)','errcatch')==21)
errmsg = msprintf(gettext("%s: Gradient function of Objective and x0 did not match "), "fmincon");
@@ -704,32 +659,16 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
if (size(sample_fGrad,1)==s(2) & size(sample_fGrad,2)==1) then
elseif (size(sample_fGrad,1)==1 & size(sample_fGrad,2)==s(2)) then
- elseif (size(sample_fGrad,1)~=1 & size(sample_fGrad,2)~=1) then
- errmsg = msprintf(gettext("%s: Wrong Input for Objective Gradient function(10th Parameter)---->Vector function is Expected"), "fmincon");
+ else
+ errmsg = msprintf(gettext("%s: Wrong Input for Objective Gradient function(3rd Parameter)---->Row Vector function of size [1 X %d] is Expected"), "fmincon",s(2));
error(errmsg);
- end
-
- function [y,check] = fGrad1(x)
- if(execstr('y=fGrad(x)','errcatch')==32 | execstr('y=fGrad(x)','errcatch')==27)
- y = 0;
- check=1;
- else
- y=fGrad(x);
- if (isreal(y)==%F) then
- y = 0;
- check=1;
- else
- check=0;
- end
- end
- endfunction
-
+ end
end
+
+ //To check for correct input of Lagrangian Hessian function from the user
if (flag2==1) then
- if (type(lHess) ~= 11 & type(lHess) ~= 13) then
- errmsg = msprintf(gettext("%s: Expected function for Hessian of Objective"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", lHess, "lHess", 10, "function");
+
if(execstr('sample_lHess=lHess(x0,1,1:no_nlc)','errcatch')==21)
errmsg = msprintf(gettext("%s: Hessian function of Objective and x0 did not match "), "fmincon");
error(errmsg);
@@ -738,29 +677,12 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
if(size(sample_lHess,1)~=s(2) | size(sample_lHess,2)~=s(2)) then
errmsg = msprintf(gettext("%s: Wrong Input for Objective Hessian function(10th Parameter)---->Symmetric Matrix function is Expected "), "fmincon");
error(errmsg);
- end
-
- function [y,check] = lHess1(x,obj,lambda)
- if(execstr('y=lHess(x,obj,lambda)','errcatch')==32 | execstr('y=lHess(x,obj,lambda)','errcatch')==27)
- y = 0;
- check=1;
- else
- y=lHess(x,obj,lambda);
- if (isreal(y)==%F) then
- y = 0;
- check=1;
- else
- check=0;
- end
- end
- endfunction
-
+ end
end
+
+ //To check for correct input of Constraint Gradient function from the user
if (flag3==1) then
- if (type(cGrad) ~= 11 & type(cGrad) ~= 13) then
- errmsg = msprintf(gettext("%s: Expected function for Gradient of Constraint function"), "fmincon");
- error(errmsg);
- end
+ Checktype("fmincon", cGrad, "cGrad", 10, "function");
if(execstr('[sample_cGrad,sample_ceqg]=cGrad(x0)','errcatch')==21)
errmsg = msprintf(gettext("%s: Gradient function of Constraint and x0 did not match "), "fmincon");
@@ -790,36 +712,149 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
errmsg = msprintf(gettext("%s: Wrong Input for Constraint Gradient function(10th Parameter) (Refer Help)"), "fmincon");
error(errmsg);
end
-
- function [y,check] = addcGrad1(x)
- if(execstr('y=addcGrad(x)','errcatch')==32 | execstr('y=addcGrad(x)','errcatch')==27)
- y = 0;
- check=1;
- else
- y=addcGrad(x);
- if (isreal(y)==%F) then
+ end
+
+ //Defining an inbuilt Objective gradient function
+ function [y,check] = fGrad1(x)
+ if flag1==1 then
+ if(execstr('y=fGrad(x)','errcatch')==32 | execstr('y=fGrad(x)','errcatch')==27)
+ y = 0;
+ check=1;
+ else
+ y=fGrad(x);
+ if (isreal(y)==%F) then
+ y = 0;
+ check=1;
+ else
+ check=0;
+ end
+ end
+ else
+ if(execstr('y=numderivative(fun,x)','errcatch')==10000)
+ y=0;
+ check=1;
+ else
+ y=numderivative(fun,x);
+ if (isreal(y)==%F) then
+ y=0;
+ check=1;
+ else
+ check=0;
+ end
+ end
+ end
+ endfunction
+
+ //Defining an inbuilt Lagrangian Hessian function
+ function [y,check] = lHess1(x,obj,lambda)
+ if flag2==1 then
+ if(execstr('y=lHess(x,obj,lambda)','errcatch')==32 | execstr('y=lHess(x,obj,lambda)','errcatch')==27)
y = 0;
check=1;
- else
- check=0;
+ else
+ y=lHess(x,obj,lambda);
+ if (isreal(y)==%F) then
+ y = 0;
+ check=1;
+ else
+ check=0;
+ end
end
+ else
+ if(execstr('[grad,y]=numderivative(fun,x)','errcatch')==10000)
+ check1=1;
+ else
+ [grad,y1]=numderivative(fun,x);
+ if (isreal(y1)==%F) then
+ check1=1;
+ else
+ check1=0;
+ end
+ end
+
+
+ if check1==0 then
+ if no_nlc~=0 then
+ if(execstr('[grad,y]=numderivative(addnlc,x)','errcatch')==10000)
+ check2=1;
+ else
+ [grad,y2]=numderivative(addnlc,x);
+ if (isreal(y2)==%F) then
+ check2=1;
+ else
+ check2=0;
+ end
+ end
+ if check2==0 then
+ y2=matrix(y2, no_nlc*s(2)*s(2),1)
+ for i = 1:s(2)*s(2)
+ y(i)=0;
+ for j = 1:no_nlc
+ y(i)= y(i) + lambda(j)*y2((i-1)*no_nlc+j);
+ end
+ end
+
+ for i = 1:s(2)*s(2)
+ y(i) = y(i)+ obj*y1(i);
+ end
+ check=0;
+ else
+ check=1;
+ end
+ else
+ check=0;
+ for i = 1:s(2)*s(2)
+ y(i) = obj*y1(i);
+ end
+ end
+ else
+ check=1;
+ y=[0];
+ end
end
- endfunction
- end
-
- //To Convert the Gradient and Hessian into Error Debugable form
+ endfunction
+ //Defining an inbuilt Constraint gradient function
+ function [y,check] = addcGrad1(x)
+ if flag3==1 then
+ if(execstr('y=addcGrad(x)','errcatch')==32 | execstr('y=addcGrad(x)','errcatch')==27)
+ y = 0;
+ check=1;
+ else
+ y=addcGrad(x);
+ if (isreal(y)==%F) then
+ y = 0;
+ check=1;
+ else
+ check=0;
+ end
+ end
+ else
+ if(execstr('y=numderivative(addnlc,x)','errcatch')==10000)
+ y=0;
+ check=1;
+ else
+ y=numderivative(addnlc,x);
+ if (isreal(y)==%F) then
+ y=0;
+ check=1;
+ else
+ check=0;
+ end
+ end
+ end
+ endfunction
+
+ //Creating a Dummy Variable for IPopt use
+ empty=[0];
- //Dummy variable which is used by Ipopt
- empty=0;
-
//Calling the Ipopt function for solving the above problem
- [xopt,fopt,status,iter,cpu,obj_eval,dual,lambda1,zl,zu,gradient,hessian1] = solveminconp (f,gradhess,A,b,Aeq,beq,lb,ub,no_nlc,no_nlic,addnlc1,flag1,fGrad1,flag2,lHess1,flag3,addcGrad1,x0,options,empty)
+ [xopt,fopt,status,iter,cpu,obj_eval,dual,lambda1,zl,zu,gradient,hessian1] = solveminconp(f,A,b,Aeq,beq,lb,ub,no_nlc,no_nlic,addnlc1,fGrad1,lHess1,addcGrad1,x0,options,empty)
//Calculating the values for the output
xopt = xopt';
exitflag = status;
- output = struct("Iterations", [],"Cpu_Time",[],"Objective_Evaluation",[],"Dual_Infeasibility",[]);
+ output = struct("Iterations", [],"Cpu_Time",[],"Objective_Evaluation",[],"Dual_Infeasibility",[],"Message","");
output.Iterations = iter;
output.Cpu_Time = cpu;
output.Objective_Evaluation = obj_eval;
@@ -857,7 +892,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
lambda.ineqlin (j) = lambda1(i)
j= j+1;
end
- lambda.ineqlin = lambda.ineqlin';
+ lambda.ineqlin = lambda.ineqlin'
end
//Converting hessian of order (1 x (numberOfVariables)^2) received from Ipopt to order (numberOfVariables x numberOfVariables)
@@ -872,7 +907,7 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
if( status~=0 & status~=1 & status~=2 & status~=3 & status~=4 & status~=7 ) then
xopt=[];
fopt=[];
- output = struct("Iterations", [],"Cpu_Time",[]);
+ output = struct("Iterations", [],"Cpu_Time",[],"Message","");
output.Iterations = iter;
output.Cpu_Time = cpu;
lambda = struct("lower",[],"upper",[],"ineqlin",[],"eqlin",[],"ineqnonlin",[],"eqnonlin",[]);
@@ -886,37 +921,54 @@ function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fmincon (varargin
case 0 then
printf("\nOptimal Solution Found.\n");
+ output.Message="Optimal Solution Found";
case 1 then
printf("\nMaximum Number of Iterations Exceeded. Output may not be optimal.\n");
+ output.Message="Maximum Number of Iterations Exceeded. Output may not be optimal";
case 2 then
printf("\nMaximum CPU Time exceeded. Output may not be optimal.\n");
+ output.Message="Maximum CPU Time exceeded. Output may not be optimal";
case 3 then
printf("\nStop at Tiny Step\n");
+ output.Message="Stop at Tiny Step";
case 4 then
printf("\nSolved To Acceptable Level\n");
+ output.Message="Solved To Acceptable Level";
case 5 then
printf("\nConverged to a point of local infeasibility.\n");
+ output.Message="Converged to a point of local infeasibility";
case 6 then
printf("\nStopping optimization at current point as requested by user.\n");
+ output.Message="Stopping optimization at current point as requested by user";
case 7 then
printf("\nFeasible point for square problem found.\n");
+ output.Message="Feasible point for square problem found";
case 8 then
printf("\nIterates diverging; problem might be unbounded.\n");
+ output.Message="Iterates diverging; problem might be unbounded";
case 9 then
printf("\nRestoration Failed!\n");
+ output.Message="Restoration Failed!";
case 10 then
printf("\nError in step computation (regularization becomes too large?)!\n");
+ output.Message="Error in step computation (regularization becomes too large?)!";
case 11 then
printf("\nProblem has too few degrees of freedom.\n");
+ output.Message="Problem has too few degrees of freedom";
case 12 then
printf("\nInvalid option thrown back by Ipopt\n");
+ output.Message="Invalid option thrown back by Ipopt";
case 13 then
printf("\nNot enough memory.\n");
+ output.Message="Not enough memory";
case 15 then
printf("\nINTERNAL ERROR: Unknown SolverReturn value - Notify Ipopt Authors.\n");
+ output.Message="INTERNAL ERROR: Unknown SolverReturn value - Notify Ipopt Authors";
else
printf("\nInvalid status returned. Notify the Toolbox authors\n");
+ output.Message="Invalid status returned. Notify the Toolbox authors";
break;
end
+
endfunction