diff options
Diffstat (limited to 'macros/fseminf.sci')
-rw-r--r-- | macros/fseminf.sci | 322 |
1 files changed, 322 insertions, 0 deletions
diff --git a/macros/fseminf.sci b/macros/fseminf.sci new file mode 100644 index 0000000..2ffcab5 --- /dev/null +++ b/macros/fseminf.sci @@ -0,0 +1,322 @@ +// Copyright (C) 2015 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Salman Anis, Harpreet Singh +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,exitflag,output,lambda,gradient,hessian] = fseminf (varargin) + // Solves a multi-variable constrainted optimization problem + // + // Calling Sequence + // xopt = fseminf(fun,x0,ntheta,seminfcon) + // xopt = fseminf(fun,x0,ntheta,seminfcon,A,b) + // xopt = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq) + // xopt = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub) + // xopt = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub,options) + // [xopt,fopt] = fseminf(.....) + // [xopt,fopt,exitflag]= fseminf(.....) + // [xopt,fopt,exitflag,output]= fseminf(.....) + // [xopt,fopt,exitflag,output,lambda]=fseminf(.....) + // [xopt,fopt,exitflag,output,lambda]=fseminf(.....) + // [xopt,fopt,exitflag,output,lambda]=fseminf(.....) + // + // Parameters + // fun : a function, representing the objective function of the problem + // x0 : a vector of doubles, containing the starting values of variables of size (1 X n) or (n X 1) where 'n' is the number of Variables + // ntheta : The number of semi-infinite constraints. + // seminfcon : a function that calculates the vector of nonlinear inequality constraints c, a vector of nonlinear equality constraints ceq, and ntheta semi-infinite constraints. See below for details. + // A : a matrix of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // b : a vector of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // Aeq : a matrix of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // beq : a vector of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // lb : a vector of double, contains lower bounds of the variables. + // ub : a vector of double, contains upper bounds of the variables. + // options : a list, containing the option for user to specify. See below for details. + + // xopt : a vector of double, the computed solution of the optimization problem. + // fopt : a double, the value of the function at x. + // exitflag : The exit status. See below for details. + // output : The structure consist of statistics about the optimization. See below for details. + // lambda : The structure consist of the Lagrange multipliers at the solution of problem. See below for details. + // + // Description + // Search the minimum of a constrained optimization problem specified by : + // Find the minimum of f(x) such that + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & f(x) \\ + // & \text{subject to} & A*x \leq b \\ + // & & Aeq*x \ = beq\\ + // & & lb \leq x \leq ub \\ + // & & c(x) \leq 0\\ + // & & ceq(x) \ = 0\\ + // & & K_i(x,w_i) \leq 0, 1 \leq i \leq n. + // \end{eqnarray} + // </latex> + // + // The routine calls Ipopt for solving the Constrained Optimization problem, Ipopt is a library written in C++. + // + // The options allows the user to set various parameters of the Optimization problem. + // It should be defined as type "list" and contains the following fields. + // <itemizedlist> + // <listitem>Syntax : options= list("MaxIter", [---], "CpuTime", [---], "GradObj", ---);</listitem> + // <listitem>MaxIter : a Scalar, containing the Maximum Number of Iteration that the solver should take.</listitem> + // <listitem>CpuTime : a Scalar, containing the Maximum amount of CPU Time that the solver should take.</listitem> + // <listitem>GradObj : a function, representing the gradient function of the Objective in Vector Form.</listitem> + // <listitem>Default Values : options = list("MaxIter", [3000], "CpuTime", [600]);</listitem> + // </itemizedlist> + // + // The exitflag allows to know the status of the optimization which is given back by Ipopt. + // <itemizedlist> + // <listitem>exitflag=0 : Optimal Solution Found </listitem> + // <listitem>exitflag=1 : Maximum Number of Iterations Exceeded. Output may not be optimal.</listitem> + // <listitem>exitflag=2 : Maximum amount of CPU Time exceeded. Output may not be optimal.</listitem> + // <listitem>exitflag=3 : Stop at Tiny Step.</listitem> + // <listitem>exitflag=4 : Solved To Acceptable Level.</listitem> + // <listitem>exitflag=5 : Converged to a point of local infeasibility.</listitem> + // </itemizedlist> + // + // For more details on exitflag see the ipopt documentation, go to http://www.coin-or.org/Ipopt/documentation/ + // + // The output data structure contains detailed informations about the optimization process. + // It has type "struct" and contains the following fields. + // <itemizedlist> + // <listitem>output.Iterations: The number of iterations performed during the search</listitem> + // <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> + // </itemizedlist> + // + // The lambda data structure contains the Lagrange multipliers at the end + // of optimization. In the current version the values are returned only when the the solution is optimal. + // It has type "struct" and contains the following fields. + // <itemizedlist> + // <listitem>lambda.lower: The Lagrange multipliers for the lower bound constraints.</listitem> + // <listitem>lambda.upper: The Lagrange multipliers for the upper bound constraints.</listitem> + // <listitem>lambda.eqlin: The Lagrange multipliers for the linear equality constraints.</listitem> + // <listitem>lambda.ineqlin: The Lagrange multipliers for the linear inequality constraints.</listitem> + // <listitem>lambda.eqnonlin: The Lagrange multipliers for the non-linear equality constraints.</listitem> + // <listitem>lambda.ineqnonlin: The Lagrange multipliers for the non-linear inequality constraints.</listitem> + // </itemizedlist> + // + // Examples + // function [y] = obj(x) + // y = (x-1)^2; + // endfunction + // function [c, ceq, K1, s] = seminfcon(x,s) + // // No finite nonlinear inequality and equality constraints + // c = []; + // ceq = []; + // // Sample set + // if isnan(s) + // // Initial sampling interval + // s = [0.01 0]; + // end + // t = 0:s(1):1; + // // Evaluate the semi-infinite constraint + // K1 = (x - 0.5) - (t - 0.5).^2; + // endfunction + // x = fseminf(obj,0.2,1,seminfcon) + // + // Authors + // Salman Anis, Harpreet Singh + + + //To check the number of input and output arguments + [lhs , rhs] = argn(); + + //To check the number of arguments given by the user + if ( rhs<4 | rhs==5 | rhs==7 | rhs>13 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while it should be 4,6,8,10,11"), "fseminf", rhs); + error(errmsg) + end + + //Storing the Input Parameters + _fun = varargin(1); + x0 = varargin(2); + ntheta = varargin(3); + seminfcon = varargin(4); + nbVar = size(x0,'*'); + + if(nbVar == 0) then + errmsg = msprintf(gettext("%s: Cannot determine the number of variables because input initial guess is empty"), "lsqcurvefit"); + error(errmsg); + end + + A = []; + b = []; + Aeq = []; + beq = []; + lb = []; + ub = []; + nlc = []; + param = list(); + + if (rhs>4) then + A = varargin(5); + b = varargin(6); + end + + if (rhs>6) then + Aeq = varargin(7); + beq = varargin(8); + end + + if (rhs>8) then + lb = varargin(9); + ub = varargin(10); + end + + if (rhs>10) then + param = varargin(10); + end + + if (size(lb,2)==0) then + lb = repmat(-%inf,nbVar,1); + end + + if (size(ub,2)==0) then + ub = repmat(%inf,nbVar,1); + end + + //Check type of variables + Checktype("fseminf", _fun, "fun", 1, "function") + Checktype("fseminf", x0, "x0", 2, "constant") + Checktype("fseminf", ntheta, "ntheta", 3, "constant") + Checktype("fseminf", seminfcon, "seminfcon", 4, ["function","constant"]) + Checktype("fseminf", A, "A", 5, "constant") + Checktype("fseminf", b, "b", 6, "constant") + Checktype("fseminf", Aeq, "Aeq", 7, "constant") + Checktype("fseminf", beq, "beq", 8, "constant") + Checktype("fseminf", lb, "lb", 9, "constant") + Checktype("fseminf", ub, "ub", 10, "constant") + Checktype("fseminf", param, "param", 10, "list") + + //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 + Checktype("fseminf", param(2*i), "maxiter", 10, "constant") + options(2) = param(2*i); //Setting the maximum number of iterations as per user entry + case "cputime" then + Checktype("fseminf", param(2*i), "cputime", 10, "constant") + options(4) = param(2*i); //Setting the maximum CPU time as per user entry + case "gradobj" then + Checktype("fseminf", param(2*i), "gradobj", 10, "string") + if(convstr(param(2*i),'l') == "on") then + function dy = graObj(x) + + endfunction + end + else + errmsg = msprintf(gettext("%s: Unrecognized parameter name %s."), "fmincon", param(2*i-1)); + error(errmsg); + end + end + + //To check and convert the 2nd Input argument (x0) to a row vector + if((size(x0,1)~=1) & (size(x0,2)~=1)) then + errmsg = msprintf(gettext("%s: Expected Vector for initial guess"), "fseminf"); + error(errmsg); + end + + if(size(x0,2)==1) then + x0=x0(:); + end + + //To check the match between fun (1st Parameter) and x0 (2nd Parameter) + if(execstr('init=_fun(x0)','errcatch')==21) then + errmsg = msprintf(gettext("%s: Objective function and x0 did not match"), "fseminf"); + error(errmsg); + end + + //Check the size of inequality constraint which should be equal to the number of variables + if ( size(A,2) ~= nbVar & size(A,2) ~= 0) then + errmsg = msprintf(gettext("%s: The number of columns in A must be the same as the number of elements of x0"), "fseminf"); + error(errmsg); + end + + nbConInEq = size(A,"r"); + + //Check the size of equality constraint which should be equal to the number of variables + if ( size(Aeq,2) ~= nbVar & size(Aeq,2) ~= 0 ) then + errmsg = msprintf(gettext("%s: The number of columns in Aeq must be the same as the number of elements of f"), "fseminf"); + error(errmsg); + end + + b = b(:); + beq = beq(:); + lb = lb(:); + ub = ub(:); + + //To check the contents of lb & ub + for i = 1:nbVar + if(ub(i)<lb(i)) then + errmsg = msprintf(gettext("%s: Problem has inconsistent variable bounds"), "fseminf"); + error(errmsg); + end + end + + if(typeof(seminfcon) == "function") + sample_S = %nan + if(execstr('[sample_c,sample_ceq,sample_K,sample_S] = seminfcon(x0,sample_S)','errcatch')==21) + errmsg = msprintf(gettext("%s: Semi-Infinite Constraint function & x0 did not match"), "fseminf"); + error(errmsg); + end + [sample_c, sample_ceq, sample_K, sample_S] = seminfcon(x0,sample_S); + + if (size(sample_c,1)~=1 & size(sample_c,1)~=0) then + errmsg = msprintf(gettext("%s: c in seminfcon should be a row vector or empty matrix"), "fseminf"); + error(errmsg) + end + + if (size(sample_ceq,1)~=1 & size(sample_ceq,1)~=0) then + errmsg = msprintf(gettext("%s: ceq in seminfcon should be a row vector or empty matrix"), "fseminf"); + error(errmsg) + end + + if (size(sample_K,"r")~=ntheta) then + errmsg = msprintf(gettext("%s: Number of rows in K should be equal to ntheta"), "fseminf"); + error(errmsg) + end + + if (size(sample_S,"c")~=2) then + errmsg = msprintf(gettext("%s: Number of columns in sampling interval should be equal to 2"), "fseminf"); + error(errmsg) + end + end + + ierr = execstr('init=_fun(x0)', "errcatch") + if ierr <> 0 then + lamsg = lasterror(); + lclmsg = "%s: Error while evaluating the function: ""%s""\n"; + error(msprintf(gettext(lclmsg), "fseminf", lamsg)); + end + + S = %nan; + + ierr = execstr('init=seminfcon(x0,S)', "errcatch") + if ierr <> 0 then + lamsg = lasterror(); + lclmsg = "%s: Error while evaluating the function: ""%s""\n"; + error(msprintf(gettext(lclmsg), "fseminf", lamsg)); + end + + function [c, ceq] = _seminfcon(x) + [c, ceq, K, S] = seminfcon(x,S) + K_max = max(K,"c"); + c= [c;K_max]; + ceq = ceq; + endfunction + + //Calling the fmincon function for solving the above problem + [xopt,fopt,exitflag,output,lambda,gradient] = fmincon(_fun,x0,A,b,Aeq,beq,lb,ub,_seminfcon,param) + +endfunction |