diff --git a/problem/examples.go b/problem/examples.go index 44ac535..ed13b44 100644 --- a/problem/examples.go +++ b/problem/examples.go @@ -18,9 +18,9 @@ Description: x1 + 2 x2 + 2 x3 <= 4 3 x1 + 4 x3 <= 6 2 x1 + x2 + 4 x3 <= 8 - x1 >= 0 - x2 >= 0 - x3 >= 0 + x1 >= -1.0 + x2 >= -1.0 + x3 >= -1.0 */ func GetExampleProblem3() *OptimizationProblem { // Setup @@ -29,7 +29,7 @@ func GetExampleProblem3() *OptimizationProblem { // Create variables x := out.AddVariableVectorClassic( 3, - 0.0, + -1.0, symbolic.Infinity.Constant(), symbolic.Continuous, ) @@ -60,3 +60,99 @@ func GetExampleProblem3() *OptimizationProblem { // } return out } + +/* +GetExampleProblem4 +Description: + + Returns the LP from this youtube video: + https://www.youtube.com/watch?v=QAR8zthQypc&t=483s + It should look like this: + Maximize 4 x1 + 3 x2 + 5 x3 + Subject to + x1 + 2 x2 + 2 x3 <= 4 + 3 x1 + 4 x3 <= 6 + 2 x1 + x2 + 4 x3 <= 8 + x1 >= 0 + x2 >= 0 + x3 >= 0 +*/ +func GetExampleProblem4() *OptimizationProblem { + // Setup + out := NewProblem("TestProblem3") + + // Create variables + x := out.AddVariableVectorClassic( + 3, + 0.0, // Use this line to implement non-negativity constraints + symbolic.Infinity.Constant(), + symbolic.Continuous, + ) + + // Create Basic Objective + c := getKVector.From([]float64{4.0, 3.0, 5.0}) + out.SetObjective( + c.Transpose().Multiply(x), + SenseMaximize, + ) + + // Create Constraints (using one big matrix) + A := getKMatrix.From([][]float64{ + {1.0, 2.0, 2.0}, + {3.0, 0.0, 4.0}, + {2.0, 1.0, 4.0}, + }) + b := getKVector.From([]float64{4.0, 6.0, 8.0}) + out.Constraints = append(out.Constraints, A.Multiply(x).LessEq(b)) + + return out +} + +/* +GetExampleProblem5 +Description: + + Returns the LP from this youtube video: + https://www.youtube.com/watch?v=QAR8zthQypc&t=483s + It should look like this: + Maximize 4 x1 + 3 x2 + 5 x3 + Subject to + x1 + 2 x2 + 2 x3 <= 4 + 3 x1 + 4 x3 <= 6 + 2 x1 + x2 + 4 x3 <= 8 + x1 >= 0 + x2 >= 0 + x3 >= 0 +*/ +func GetExampleProblem5() *OptimizationProblem { + // Setup + out := NewProblem("TestProblem3") + + // Create variables + x := out.AddVariableVector(3) + + // Create Basic Objective + c := getKVector.From([]float64{4.0, 3.0, 5.0}) + out.SetObjective( + c.Transpose().Multiply(x), + SenseMaximize, + ) + + // Create Constraints (using one big matrix) + A := getKMatrix.From([][]float64{ + {1.0, 2.0, 2.0}, + {3.0, 0.0, 4.0}, + {2.0, 1.0, 4.0}, + }) + b := getKVector.From([]float64{4.0, 6.0, 8.0}) + out.Constraints = append(out.Constraints, A.Multiply(x).LessEq(b)) + + // Add non-negativity constraints + for _, varII := range x { + out.Constraints = append( + out.Constraints, + varII.GreaterEq(0.0), + ) + } + return out +} diff --git a/problem/optimization_problem.go b/problem/optimization_problem.go index c69d505..349e4c5 100644 --- a/problem/optimization_problem.go +++ b/problem/optimization_problem.go @@ -568,26 +568,6 @@ func (op *OptimizationProblem) LinearEqualityConstraintMatrices() (symbolic.KMat return COut2, dOut2, nil } -// func (op *OptimizationProblem) Simplify() OptimizationProblem { -// // Create a new optimization problem -// newProblem := NewProblem(op.Name + " (Simplified)") - -// // Add all variables to the new problem -// for _, variable := range op.Variables { -// newProblem.Variables = append(newProblem.Variables, variable) -// } - -// // Add all constraints to the new problem -// for _, constraint := range op.Constraints { -// newProblem.Constraints = append(newProblem.Constraints, constraint) -// } - -// // Set the objective of the new problem -// newProblem.Objective = op.Objective - -// return newProblem -// } - /* ToProblemWithAllPositiveVariables Description: @@ -602,6 +582,7 @@ Description: func (op *OptimizationProblem) ToProblemWithAllPositiveVariables() (*OptimizationProblem, error) { // Setup newProblem := NewProblem(op.Name + " (All Positive Variables)") + epsMagic := 1e-8 // TODO(Kwesi): Make this a parameter OR a constant in the package. // For each variable, let's create two new variables // and set the original variable to be the difference of the two @@ -610,20 +591,35 @@ func (op *OptimizationProblem) ToProblemWithAllPositiveVariables() (*Optimizatio // Setup xII := op.Variables[ii] + // Expression for the positive and negative parts + var expressionForReplacement symbolic.Expression = symbolic.K(0.0) + // Create the two new variables - newProblem.AddVariableClassic(0.0, symbolic.Infinity.Constant(), symbolic.Continuous) - nVariables := len(newProblem.Variables) - newProblem.Variables[nVariables-1].Name = xII.Name + " (+)" - variablePositivePart := newProblem.Variables[nVariables-1] + // - Positive Part + positivePartExists := xII.Upper >= 0 + positivePartExists = positivePartExists && !ConstraintIsRedundantGivenOthers(xII.LessEq(0.0-epsMagic), op.Constraints) + if positivePartExists { + newProblem.AddVariableClassic(0.0, symbolic.Infinity.Constant(), symbolic.Continuous) + nVariables := len(newProblem.Variables) + newProblem.Variables[nVariables-1].Name = xII.Name + " (+)" + variablePositivePart := newProblem.Variables[nVariables-1] + expressionForReplacement = expressionForReplacement.Plus(variablePositivePart) + } - newProblem.AddVariableClassic(0.0, symbolic.Infinity.Constant(), symbolic.Continuous) - nVariables = len(newProblem.Variables) - newProblem.Variables[nVariables-1].Name = xII.Name + " (-)" - variableNegativePart := newProblem.Variables[nVariables-1] + // - Negative Part + negativePartExists := xII.Lower < 0 + negativePartExists = negativePartExists && !ConstraintIsRedundantGivenOthers(xII.GreaterEq(0.0), op.Constraints) + if negativePartExists { + newProblem.AddVariableClassic(0.0, symbolic.Infinity.Constant(), symbolic.Continuous) + nVariables := len(newProblem.Variables) + newProblem.Variables[nVariables-1].Name = xII.Name + " (-)" + variableNegativePart := newProblem.Variables[nVariables-1] + + expressionForReplacement = expressionForReplacement.Minus(variableNegativePart) + } // Set the original variable to be the difference of the two new variables - mapFromOriginalVariablesToNewExpressions[xII] = - variablePositivePart.Minus(variableNegativePart) + mapFromOriginalVariablesToNewExpressions[xII] = expressionForReplacement } // Now, let's create the new constraints by replacing the variables in the @@ -697,8 +693,9 @@ func (problemIn *OptimizationProblem) ToLPStandardForm1() (*OptimizationProblem, } // Add all constraints to the new problem + problemWithPositivesAndCleanedConstraints := problemWithAllPositiveVariables.WithAllPositiveVariableConstraintsRemoved() slackVariables := []symbolic.Variable{} - for _, constraint := range problemWithAllPositiveVariables.Constraints { + for _, constraint := range problemWithPositivesAndCleanedConstraints.Constraints { // Create a new expression by substituting the variables according // to the map we created above oldLHS := constraint.Left() @@ -827,12 +824,70 @@ func (problemIn *OptimizationProblem) ToLPStandardForm1() (*OptimizationProblem, problemWithAllPositiveVariables.Objective.Sense, ) - // fmt.Printf("The slack variables are: %v\n", slackVariables) + // Simplify The Constraints If Possible + problemInStandardForm.SimplifyConstraints() // Return the new problem and the slack variables return problemInStandardForm, slackVariables, nil } +/* +WithAllPositiveVariableConstraintsRemoved +Description: + + Returns a new optimization problem that is the same as the original problem + but with all constraints of the following form removed: + x >= 0 + 0 <= x + Where x is a variable in the problem. + This is useful for removing redundant constraints that are already implied by the variable bounds. +*/ +func (op *OptimizationProblem) WithAllPositiveVariableConstraintsRemoved() *OptimizationProblem { + // Setup + newProblem := NewProblem(op.Name) + + // Copy the variables + for _, variable := range op.Variables { + newProblem.Variables = append(newProblem.Variables, variable) + } + + // Copy the constraints + for _, constraintII := range op.Constraints { + // Check if the constraint is a x >= 0 constraint + if symbolic.SenseGreaterThanEqual == constraintII.ConstrSense() { + lhsContains1Variable := len(constraintII.Left().Variables()) == 1 + rhs, rhsIsConstant := constraintII.Right().(symbolic.K) + if lhsContains1Variable && rhsIsConstant { + if float64(rhs) == 0.0 { + // If the constraint is of the form x >= 0, we can remove it + continue + } + } + } + + // Check if the constraint is a 0 <= x constraint + if symbolic.SenseLessThanEqual == constraintII.ConstrSense() { + rhsContains1Variable := len(constraintII.Left().Variables()) == 1 + lhs, lhsIsConstant := constraintII.Right().(symbolic.K) + if rhsContains1Variable && lhsIsConstant { + if float64(lhs) == 0.0 { + // If the constraint is of the form 0 <= x, we can remove it + continue + } + } + } + + // Otherwise, we can keep the constraint + newProblem.Constraints = append(newProblem.Constraints, constraintII) + } + + // Copy the objective + newProblem.Objective = op.Objective + + // Return the new problem + return newProblem +} + /* CheckIfLinear Description: diff --git a/testing/problem/optimization_problem_test.go b/testing/problem/optimization_problem_test.go index 427fc41..2be4cb7 100644 --- a/testing/problem/optimization_problem_test.go +++ b/testing/problem/optimization_problem_test.go @@ -1891,7 +1891,9 @@ Description: The problem will have: - a linear objective - 3 variables, - - and a single linear VECTORE equality constraint. + - and a single linear VECTOR equality constraint (n_ineq = 1). + Because each variable can be positive or negative, the resulting + linear equality constraint matrix should have 3 rows and 3*2+n_ineq columns. */ func TestOptimizationProblem_LinearEqualityConstraintMatrices7(t *testing.T) { // Constants @@ -1966,6 +1968,110 @@ func TestOptimizationProblem_LinearEqualityConstraintMatrices8(t *testing.T) { } } +/* +TestOptimizationProblem_LinearEqualityConstraintMatrices9 +Description: + + Tests the LinearEqualityConstraintMatrices function with a problem + that led to panics in the field. + The problem is Problem4 from our examples file. + The problem will have: + - a linear objective + - 3 variables (each labeled as positive via the LB input to AddVariableVectorClassic), + - and a single linear VECTOR equality constraint. +*/ +func TestOptimizationProblem_LinearEqualityConstraintMatrices9(t *testing.T) { + // Constants + p1 := problem.GetExampleProblem4() + + // Transform p1 into the standard form + p1Standard, _, err := p1.ToLPStandardForm1() + if err != nil { + t.Errorf("unexpected error: %v", err) + } + + // Attempt to Call LinearEqualityConstraintMatrices + A, b, err := p1Standard.LinearEqualityConstraintMatrices() + if err != nil { + t.Errorf("unexpected error: %v", err) + } + + // Check that the number of rows is as expected. + if A.Dims()[0] != 3 { + t.Errorf("expected the number of rows to be %v; received %v", + 3, A.Dims()[0]) + } + + // Check that the number of columns is as expected. + nVariables1 := len(p1.Variables) + nInequalityConstraints1 := p1.Constraints[0].Left().Dims()[0] + if A.Dims()[1] != nVariables1+nInequalityConstraints1 { + t.Errorf("expected the number of columns to be %v; received %v", + nVariables1+nInequalityConstraints1, A.Dims()[1]) + } + + // Check that the number of elements in b is as expected. + if len(b) != 3 { + t.Errorf("expected the number of elements in b to be %v; received %v", + 3, len(b)) + } +} + +/* +TestOptimizationProblem_LinearEqualityConstraintMatrices10 +Description: + + Tests the LinearEqualityConstraintMatrices function with a problem + that led to panics in the field. + The problem is Problem4 from our examples file. + The problem will have: + - a linear objective + - 3 variables (each labeled as positive implicitly via the final 3 constraints), + - and a single linear VECTOR equality constraint. + The resulting problem should have 6 constraints (3 coming from inequality constraint) + and 3 from the positivity constraints. + The resulting linear equality constraint matrix should have 3 rows and 3+3 columns + (3 variables, 3 inequality constraints, and 3 positivity constraints). +*/ +func TestOptimizationProblem_LinearEqualityConstraintMatrices10(t *testing.T) { + // Constants + p1 := problem.GetExampleProblem5() + + // Transform p1 into the standard form + p1Standard, slackVariables, err := p1.ToLPStandardForm1() + if err != nil { + t.Errorf("unexpected error: %v", err) + } + + // Attempt to Call LinearEqualityConstraintMatrices + A, b, err := p1Standard.LinearEqualityConstraintMatrices() + if err != nil { + t.Errorf("unexpected error: %v", err) + } + + // Check that the number of rows is as expected. + + if A.Dims()[0] != 3 { + t.Errorf("expected the number of rows to be %v; received %v", + 3, A.Dims()[0]) + } + + // Check that the number of columns is as expected. + nVariables1 := len(p1.Variables) + nSlackVariables1 := len(slackVariables) + expectedNumVariables := nVariables1 + nSlackVariables1 + if A.Dims()[1] != expectedNumVariables { + t.Errorf("expected the number of columns to be %v; received %v", + expectedNumVariables, A.Dims()[1]) + } + + // Check that the number of elements in b is as expected. + if len(b) != 3 { + t.Errorf("expected the number of elements in b to be %v; received %v", + 3, len(b)) + } +} + /* TestOptimizationProblem_ToProblemWithAllPositiveVariables1 Description: @@ -2017,6 +2123,62 @@ func TestOptimizationProblem_ToProblemWithAllPositiveVariables1(t *testing.T) { } } +/* +TestOptimizationProblem_ToProblemWithAllPositiveVariables2 +Description: + + Tests the ToProblemWithAllPositiveVariables function with a simple problem + that has: + - a constant objective + - 2 variables, + - and two scalar linear inequality constraints. + One of the variables is purely positive, while the other is purely negative. + The result should be a problem with 2 variables and 2 constraints. +*/ +func TestOptimizationProblem_ToProblemWithAllPositiveVariables2(t *testing.T) { + // Constants + p1 := problem.NewProblem("TestOptimizationProblem_ToProblemWithAllPositiveVariables2") + vv1 := p1.AddVariableVector(2) + // Add constraints + c1 := vv1.AtVec(0).GreaterEq(1.0) + c2 := vv1.AtVec(1).LessEq(-2.0) + + p1.Constraints = append(p1.Constraints, c1) + p1.Constraints = append(p1.Constraints, c2) + + // Create good objective + p1.Objective = *problem.NewObjective( + symbolic.K(3.14), + problem.SenseMaximize, + ) + + // Algorithm + p2, err := p1.ToProblemWithAllPositiveVariables() + if err != nil { + t.Errorf("unexpected error: %v", err) + } + + // Check that the number of variables is as expected. + if len(p2.Variables) != 2 { + t.Errorf("expected the number of variables to be %v; received %v", + 2, len(p2.Variables)) + } + + // Check that the number of constraints is as expected. + if len(p2.Constraints) != 2 { + t.Errorf("expected the number of constraints to be %v; received %v", + 2, len(p2.Constraints)) + } + + // Verify that the new constraints contain two variables in the left hand side + for _, c := range p2.Constraints { + if len(c.Left().Variables()) != 1 { + t.Errorf("expected the number of variables in the left hand side to be %v; received %v", + 1, len(c.Left().Variables())) + } + } +} + /* TestOptimizationProblem_ToLPStandardForm1_1 Description: @@ -2032,8 +2194,7 @@ func TestOptimizationProblem_ToLPStandardForm1_1(t *testing.T) { // Constants p1 := problem.NewProblem("TestOptimizationProblem_ToLPStandardForm1_1") v1 := p1.AddVariable() - p1.AddVariable() - c1 := v1.GreaterEq(1.0) + c1 := v1.GreaterEq(-1.0) p1.Constraints = append(p1.Constraints, c1) @@ -2055,7 +2216,7 @@ func TestOptimizationProblem_ToLPStandardForm1_1(t *testing.T) { expectedNumVariables += len(p1.Constraints) // slack variables if len(p2.Variables) != expectedNumVariables { t.Errorf("expected the number of variables to be %v; received %v", - 2, len(p2.Variables)) + expectedNumVariables, len(p2.Variables)) } // Check that the number of constraints is as expected. @@ -2393,7 +2554,7 @@ func TestOptimizationProblem_ToLPStandardForm1_7(t *testing.T) { // Check that the number of variables is as expected. expectedNumVariables := 0 - expectedNumVariables += 2 * len(p1.Variables) // original variables (positive and negative halfs) + expectedNumVariables += 3 // original variables (one with positive and negative halfs, the other with only the positive part) if len(p2.Variables) != expectedNumVariables { t.Errorf("expected the number of variables to be %v; received %v", expectedNumVariables, len(p2.Variables)) @@ -2505,7 +2666,7 @@ Description: C = [ 0 0 1 -1 0 0 0 -1 0 ] [ 0 0 0 0 1 -1 0 0 -1 ] and - b = [ 1 2 3 ] + b = [ -1 -2 -3 ] By creating a problem with 3 variables and 3 linear inequality constraints. The results should produce equality constraint matrix C with 3 rows and 6 columns and a vector b with 3 elements. @@ -2514,9 +2675,9 @@ func TestOptimizationProblem_ToLPStandardForm1_9(t *testing.T) { // Constants p1 := problem.NewProblem("TestOptimizationProblem_ToLPStandardForm1_9") vv1 := p1.AddVariableVector(3) - c1 := vv1.AtVec(0).GreaterEq(1.0) - c2 := vv1.AtVec(1).GreaterEq(2.0) - c3 := vv1.AtVec(2).GreaterEq(3.0) + c1 := vv1.AtVec(0).GreaterEq(-1.0) + c2 := vv1.AtVec(1).GreaterEq(-2.0) + c3 := vv1.AtVec(2).GreaterEq(-3.0) p1.Constraints = append(p1.Constraints, c1) p1.Constraints = append(p1.Constraints, c2)