Euler
ReferenceSimplex and integer Simplex routines.
The Simplex algorithm is a core function of Euler. It can handle inequalities and equalities, restricted and unrestricted variables. There is a utility function to simplify the handling and the error checking.
Euler contains the optimization library LPSOLVE, which can be used alternatively. E.g., the function ilpsolve() is more efficient than the simple branch and bound method intsimplex() for integer programs.
function overwrite simplex (A:real, b:real column, c:real vector,
eq=-1, restrict=1, max:integer=0, min:integer=1, check:integer=1)
Minimize c.x with respect to A.x<=b and x>=0 or other restrictions.
Simplex method for maximization or minimization of a linear
function with linear constrains, and restricted or unrestricted
variables. The function calls the built-in Simplex algorithm.
Minimizes or maximizes c.x with the conditions A.x <= b, or A.x >=
b, or A.x == b, and optionally x>=0. The type of the condition is
contained in a column vector eq. It is in each row -1 for <=, 1 for
>=, 0 for equality. The column vector restrict determines, which
variables should be restricted to be non-negative (1), or
non-positive (-1), or unrestricted (0). If max=1, the function
maximizes the target function.
A : real matrix (nxm)
b : real column vector (nx1)
c : real row vector (1xm)
eq : real column vector with entries -1, 0, or 1 (nx1)
restrict : real row vector with entries 0, 1, or -1 (1xm)
max,min : flags to maximize and minimize (set either of them)
check: If true, an error will occur for unbounded or not feasible
problems. If false, you need to check the r flag for the result.
The return value is {x,r,i}, where x is the solution, and r has the
following meaning:
r=0 : success,
r=-1 : no feasible point,
r=1 : the problem is not bounded.
i is a column vector, which is 1 for each active side condition.
It must be remarked, that an equality condition generates two
conditions. So i might be longer than the original number of
conditions.
>A=[1,1;4,5]; b=[1000;4500]; // x+y<=1000, 4x+5y<=4500
>eq=[-1;-1]; // both equations <= (the default)
>c=[5,6]; // 4x+5y -> Max.
>restr=[1,1]; // x,y>=0 (the default)
>simplex(A,b,c,eq,restr,>max,>check) // check for errors
500
500
This function calls the built-in simplex() algorithm, which always
minimizes and always returns a result flag.
See:
intsimplex (Linear Programming)
function feasibleArea (A:real, b:real column, eq=-1, restrict=1)
Computes the corners of the set with A.x<=b A can have two columns only. Returns a 2xn matrix with rows x and y. These vectors can be used to plot the feasible area with the polt2d function and the parameter polygon=1. >A=[1,1;4,5]; b=[1000;4500]; // x+y<=1000, 4x+5y<=4500 >eq=[-1;-1]; // both equations <= (the default) >c=[5,6]; // 4x+5y -> Max. >restr=[1,1]; // x,y>=0 (the default) >X=feasibleArea(A,b,eq,restr); // determine polygon >plot2d(X[1],X[2],>filled,style="/"); >x=simplex(A,b,c,eq,restr,>max,>check); >plot2d(x[1],x[2],>points,>add);
function intsimplex (A:real, b:real column, c:real vector,
eq=-1, restrict=1, max:integer=0, min:integer=1, i="all",
check=false, ire:integer vector=none)
Minimize c.x=b under the conditions A.x<=b, x integer, x>=0.
This is the branch and bound algorithm for integer linear problems.
It is implemented in the Euler language. There is a faster and more
sophisticated algorithm ilpsolve().
Minimizes or maximizes c.x with the conditions A.x <= b, or A.x >=
b, or A.x == b, x integer. and optionally x>=0. The parameter i is a
row vector of indices of variables, which should be positive. If
i="all", all variables should be integer. This is the default.
A : real matrix (nxm)
b : real column vector (nx1)
c : real row vector (1xm)
eq : real column vector with entries -1, 0, or 1 (nx1)
restrict : vector with non-negative restricted variables.
max,min : flags to maximize and minimize
i : real row vector with indices of integer variables
ire : Alternative to i. Row vector with 0,1.
Note that i and restrict work in different ways. This is due to the
compatibility with ilpsolve().
check : If true, the function will throw an error, if the problem
is unbounded, or has no feasible point.
The function returns {x,r,a}, where x is the solution of the problem
if r=0. a is the optimal value.
>A=random(20,3)+1; b=ones(20,1)*1000; // random problem
>c=ones(1,3);
>intsimplex(A,b,c,>max,>check)
193
348
41
See:
ilpsolve (Linear Programming)
For the LPSOLVE package has been ported by Peter Notebaert for Euler. The full documentation is available on
defaultilpsolvedll:=0;
function startlpsolve ()
Start the LPSOLVE DLL.
function ilpsolve (A, b, c, eq = -1,
vlb = [], vub = [], i = "all",
scalemode = 0, keep = 0, max = 0, ire = none)
Minimize f.x subject a.x<=b, vlb<=x<=vub, and integer x.
This routine is contained in the lpsolve DLL, which is loaded, if
you load the lpsolve package with "load lpsolve". The routines have
been ported to Euler by Peter Notebaert.
The function solves the linear integer problem
max v = f'*x
a*x <> b
vlb <= x <= vub
x[i] are integer
Arguments: The first four arguments are required.
f: n vector of coefficients for a linear objective function.
a: m by n matrix representing linear constraints.
b: m vector of right sides for the inequality constraints.
e: m vector that determines the sense of the inequalities:
e(i) = -1 ==> Less Than
e(i) = 0 ==> Equals
e(i) = 1 ==> Greater Than
vlb: n vector of lower bounds. If empty or omitted,
then the lower bounds are set to zero.
vub: n vector of upper bounds. May be omitted or empty.
xint: vector of integer variables. May be omitted or empty.
scalemode: scale flag. Off when 0 or omitted.
keep: Flag for keeping the lp problem after it's been solved.
If omitted, the lp will be deleted when solved.
Output: The function returns {x,obj,duals}
x: Optimal value of the decision variables.
obj: Optimal value of the objective function.
duals: solution of the dual problem.
Loads the DLL once, when it is called.
>A=random(50,5)+1; b=ones(50,1)*1000; // random problem
>c=ones(1,5);
>load lpsolve;
>ilpsolve(A,b,c,i=[],>max,>check) // no integer problem
61.7789496937
156.76464272
76.0735082666
60.3330148438
232.179101403
>ilpsolve(A,b,c,>max,>check)
62
157
76
60
See:
intsimplex (Linear Programming),
lpsolve (LPSOLVE package)
This method minimizes f(x) under conditions g_i(x)<=0 using a global minimizer for f(x)-c*sum(log(-g_i(x)),i), and letting c to 0. The method for general convex f need the gradient and the Hessian for a good global minimizer.
function newtonbarrier (f$:string, df$:string, Hf$:string,
A:real, b:real column,
x:real vector, lambda:positive real=1,
c:positive real=0.1, cfactor:positive real=0.5,
history=false, eps=epsilon())
Newton-Barrier method starting from interior x. This function can minimize a convex function function f(x) subject to conditions Ax <= b. It needs the gradient df and the Hessian Hf of f as functions. The start point x must be an interior point with Ax < b. Note: f must be convex. f,df,Hf : Functions or call collections taking a row vector as a parameter. A,b : Conditions Ax<=b x : Start point with Ax<b lambda,c,cfactor : Parameters for the algorithms. lambda is the first step size for a minimization into a direction of descent. c is the starting value in the Newton-Barrier function, and cfactor is the factor at which c is reduced in each step. history : If true, the function returns a matrix with one step x in each column. Else it returns the last x only. eps : The absolute accuracy of the solution. The algorithm returns, when |x-xnew|<epsilon.