Get rid of `beta_ineq` and `A_ineq`

[a-torgovitsky]

Going forward, we are going to always assume that the underlying problem is standard form:

min c'x x >= 0 Ax = b

any LP can be written in this way

So we don't need (or want) inequalities, except for the bound constraints on x.

[conroylau]

Going forward, we are going to always assume that the underlying problem is standard form:

May I clarify with you that in getting rid of the inequalities (i.e. they are matrices like A_ineq and beta_ineq in the current module), do you mean that I should only allow users to pass through equality constraints to the module (i.e. if they have inequality constraints in their problem, they need to make it in standard form before using functions like estbounds)? Or do you want me to have a function that would transform the matrices that correspond to inequality constraints into standard form, being further passing them to the functions that do the tests?

Thank you!

[a-torgovitsky]

It's a good question.

When I wrote this, I assumed that automatically transforming a general form LP problem into standard form is a rather difficult programming problem. I might be wrong about that. When I reformulate linear programs I just do it by intuition/feel, but maybe there is a clean step-by-step algorithm that can be used? Could you research this and find out? If so, perhaps someone has already coded this in R and we can use their package directly.

[conroylau]

If we are primarily concerned with the inequality constraints in linear programs, I think reformulating them in standard form is mainly related to the change of the inequality signs and adding the slack or surplus variables. If this is the case, the simplex function from the boot package can be useful for us.

The detailed description of the function simplex and the object returned from the function can be found here. But for convenience, I am illustrating how it might be useful for us below.

To begin with, the arguments of the simplex function is as follows:

• a is a vector for the objective function
• A1 is an m1 x n matrix for <= type constraints and b1 is the corresponding RHS vector
• A2 is an m2 x n matrix for >= type constraints and b2 is the corresponding RHS vector
• A3 is an m3 x n matrix for = type constraints and b3 is the corresponding RHS vector
• maxi is set to TRUE for max. problems and FALSE for min. problems
• n.iter is the total number of iterations for each phase of the simplex method

All values in b1, b2 and b3 have to be nonnegative.

Suppose the linear program that we want to solve has the following constraints (where I set the coefficients randomly): x + 2y <= 7 3x + 4y >= 8 5x + 6y = 9 x, y >= 0

Then the following code can be used to transform the above into standard form:

library(boot)
simp = simplex(a = c(1, 1), 
               A1 = c(1, 2),
               b1 = 7,
               A2 = c(3, 4),
               b2 = 8,
               A3 = c(5, 6),
               b3 = 9,
               maxi =  FALSE, 
               n.iter = 0)
simp$A

In the simplex function, simp$A is returning the final constraints matrix.

• If there is a feasible solution, then it will have n+m1+m2 columns , where the final m1+m2 columns correspond to slack and surplus variables
• If there is no feasible solution, then there will be an additional m1+m2+m3 columns for the artificial variables used to solve the first phase of the problem, but I think these columns are irrelevant for our case

Since I set n.iter = 0, there is no iteration in the two-step simplex method. Thus, there is no feasible solution from the simplex function. In the above example, the following matrix is returned:

     [,1] [,2] [,3] [,4] [,5] [,6]
[1,]    1    2    1    0    0    0
[2,]    3    4    0   -1    1    0
[3,]    5    6    0    0    0    1

where

• Columns 1 and 2 correspond to the coefficients in the constraint
• Columns 3 and 4 correspond to the slack and surplus variables due to the <= and >= constraints respectively
• Columns 5 and 6 corresponds to the artificial variables that are added to the >= constraint and equality constraint respectively in phase one of the two-step simplex method

Thus, we can transform A_ineq and beta_ineq into standard form using the above function by setting a to be some vectors of length equals to the number of columns in A_ineq (I think the objective function does not matter here as we are primarily concerned with the constraints matrix) and setting n.iter = 0. Then, we can just extract the first n + m1 + m2 columns of the matrix as the updated constraints matrix in standard form. In the above example, we can extract the first four columns of simp$A to be the matrix we want.

[a-torgovitsky]

Ok that's great!

Am I understanding correctly that by setting n.iter = 0, the function doesn't actually try to solve the linear program? But still returns the problem in standard form?

I think what we want to do is:

1. Provide the user with a function called standard.form which is basically just a simple wrapper for simplex in the boot package, used as you did above.
2. Do some testing to make sure this is working as expected.
   • Start with an LP in general form.
   • Solve with Gurobi.
   • Convert to standard form with standard.form and solve with Gurobi.
   • Check that these match.
3. Add unit tests which start with a general form LP, manually convert to standard form, attempt to solve one of the simpler procedures in lpinfer (such as mincriterion), and compare the results with when the LP is automatically converted to standard form with standard.form.

[conroylau]

Yes, the simplex function will not solve the linear program if we set n.iter = 0.

Here is a zip file that contains the three files that correspond to the three points that you mentioned above:

• lpsf.R contains the function standard.form
• test_eg.R contains some linear programming examples that I have taken from some online notes. I have compared the solutions (both the optimal value and the optimal arguments) obtained from solving the LP in the general form and from the LP in the standard form (using standard.form). They give the same solution.
• test_lpsf.R contains an example where I randomly design inequality constraints for the mincriterion problem. The solution computed from using beta_shp_ineq and A_shp_ineq and the solution computed only using equality constraints are the same.

[a-torgovitsky]

Ok you should set it up as a unit test using this package

[conroylau]

I have updated the file for the unit test using the testthat package. The updated file can be found here.

In addition, I am thinking what would be the best way to organize the constraints around the lpmodel object (as per issue #21). Currently I can think of the following two ways:

• For all testing procedures, the input only accepts standard form constraints. If the user would like to pass some inequality constraints, they need to first convert them to standard form using the standard.form function. Then, they can pass the updated matrices to the testing procedure.
• Alternatively, add an argument like sense.shp that refers to the sense of the constraints. This is an optional argument to lpmodel. If sense.shp is NULL, it means that the constraints are all equality constraints. Otherwise, users may specify the inequality or equality sign for each row of the matrix. Inside the testing procedure (such as dkqs), use the standard.form function in the beginning of the procedure to make the constraints in standard form.

Do you think the above works, or would there be a better approach? Thank you!

[a-torgovitsky]

I would prefer the first approach, where all of our procedures only accept standard form, and the user has to convert by themselves.

I am worried that automatic conversions may lead to strange results in certain situations (perhaps the rows get switched, or something with the slack variables, etc.), so I think it is safer to have that done manually.

[conroylau]

I am dropping the matrices A_shp_ineq and beta_shp_ineq in the argument of the functions so that we only accept matrices of equality constraints in our module.

For the inequality constraints in the LP that are not coming from A_shp_ineq and beta_shp_ineq, such as step 2 of the estbounds procedure below, may I know that should I also make it in standard form so that all the LP in the module are in standard form for consistency? Thank you!

[a-torgovitsky]

No. The point is that we want the user to input the LP in a specified form so that we know exactly what we are getting from them. This has no bearing on how the problems are specified later inside the code. You can do it in whatever form is most natural, and Gurobi will take care of it.

[conroylau]

I see, no problem. I have kept how the constraints are set inside the code and have restricted the input of the LP to be matrices for equality constraints.