← All NMath Code Examples
using System;
using CenterSpace.NMath.Core;
namespace CenterSpace.NMath.Examples.CSharp
{
/// <summary>
/// A .NET example in C# showing how to solve the Klee Minty cube linear
/// programming problem using both primal and dual simplex methods and
/// various pivoting strategies.
/// The Klee Minty cube is a unit cube with perturbed corners. Klee and
/// Minty showed that Dantzigs simplex algorithm can have worst-case
/// performance with this cube as the feasible region.
/// </summary>
public class PrimalDualSimplexExample
{
static void Main( string[] args )
{
int n = 200;
var lpProblem = KleeMintyCube( n );
Console.WriteLine( "Klee Minty Cube LP problem using Primal Simplex:" );
Console.WriteLine( "---------------------------------------------------" );
SolveKleeMintyProblemPrimal( lpProblem );
Console.WriteLine( "\n" );
Console.WriteLine( "Klee Minty Cube LP problem using Dual Simplex:" );
Console.WriteLine( "---------------------------------------------------" );
SolveKleeMintyProblemDual( lpProblem );
Console.WriteLine();
Console.WriteLine( "Press Enter Key" );
Console.Read();
}
/// <summary>
/// Solve the Klee Minty cube linear programming problem with the primal simplex
/// solver using different pivoting strategies.
/// </summary>
/// <param name="kleeMintyProblem">Klee Minty cube linear programming problem.</param>
static void SolveKleeMintyProblemPrimal( LinearProgrammingProblem kleeMintyProblem )
{
// Create a solver parameters object with the desired values.
// Here we set costing, or pivot strategy for selecting which
// variable enters the basis on a pivot, the maximum number of pivots
// to perform, tell the solver that we want to minimize, not maximize,
// the objective function. If the maximum number of pivots is exceeded
// the Result property of the solver will be SolveResult.SolverInterrupted.
var solverParams = new PrimalSimplexSolverParams
{
Costing = PrimalSimplexCosting.Partial, // Use partial pivoting.
MaxPivotCount = 5000, // Do not perform more than 5000 pivots
Minimize = true // minimize, not maximize the objective function.
};
// Create a primal simplex solver and solve with the above options.
var solver = new PrimalSimplexSolver();
solver.Solve( kleeMintyProblem, solverParams );
int partialPivotCount = solver.PivotCount;
Console.WriteLine( "Primal Simplex Solver, Partial pivot strategy result: " + solver.Result );
Console.WriteLine( "Partial pivot strategy pivot count = " + partialPivotCount );
Console.WriteLine();
// Now solve with best reduced cost pivoting strategy. This is the
// classic Dantzig rule.
solverParams.Costing = PrimalSimplexCosting.BestReducedCost;
solver.Solve( kleeMintyProblem, solverParams );
int bestReducedCostPivotCount = solver.PivotCount;
Console.WriteLine( "Primal Simplex Solver, Best Reduced Cost result: " + solver.Result );
Console.WriteLine( "Best reduced cost pivot strategy pivot count = " + bestReducedCostPivotCount );
Console.WriteLine();
// Is one pivoting strategy preferable?
if ( partialPivotCount < bestReducedCostPivotCount )
{
Console.WriteLine( "For Primal Simplex prefer steepest edge for Klee Minty problem" );
}
else if ( bestReducedCostPivotCount < partialPivotCount )
{
Console.WriteLine( "For Primal Simplex prefer best reduced cost for Klee Minty problem" );
}
else
{
Console.WriteLine( "For Primal Simplex both pivoting strategies yield the same number of pivots." );
}
}
/// <summary>
/// Solve the Klee Minty cube linear programming problem with the dual simplex
/// solver using different pivoting strategies.
/// </summary>
/// <param name="kleeMintyProblem">Klee Minty cube linear programming problem.</param>
static void SolveKleeMintyProblemDual( LinearProgrammingProblem kleeMintyProblem )
{
// Create a solver parameters object with the desired values. Note that
// partial pivoting is not an available option for the dual simplex
// solver.
var solverParams = new DualSimplexSolverParams
{
Costing = DualSimplexCosting.SteepestEdge, // Use steepest edge pivoting.
MaxPivotCount = 5000,
Minimize = true // minimize, not maximize the objective function.
};
// Create a dual simplex solver and solve with the above options.
var solver = new DualSimplexSolver();
solver.Solve( kleeMintyProblem, solverParams );
int steepestEdgePivotCount = solver.PivotCount;
Console.WriteLine( "Dual Simplex Solver, Steepest Edge result: " + solver.Result );
Console.WriteLine( "Steepest edge pivot strategy pivot count = " + steepestEdgePivotCount );
Console.WriteLine();
// Now solve with the best reduced cost pivot strategy. This is the
// classic Dantzig rule.
solverParams.Costing = DualSimplexCosting.BestReducedCost;
solver.Solve( kleeMintyProblem, solverParams );
int bestReducedCostPivotCount = solver.PivotCount;
Console.WriteLine( "Dual Simplex Solver, Best Reduced Cost result: " + solver.Result );
Console.WriteLine( "Best reduced cost pivot strategy pivot count = " + bestReducedCostPivotCount );
Console.WriteLine();
// Is one pivoting strategy preferable?
if ( steepestEdgePivotCount < bestReducedCostPivotCount )
{
Console.WriteLine( "For Dual Simplex prefer steepest edge for Klee Minty problem" );
}
else if ( bestReducedCostPivotCount < steepestEdgePivotCount )
{
Console.WriteLine( "For Dual Simplex prefer best reduced cost for Klee Minty problem" );
}
else
{
Console.WriteLine( "For Dual Simplex both pivoting strategies yield the same number of pivots." );
}
}
/// <summary>
/// Create the Klee Minty cube linear programming problem with the given dimension.
/// </summary>
/// <param name="n">The desired number of vertices.</param>
/// <returns><c>LinearProgramming</c> object describing the Klee Minty
/// cube linear programming problem.</returns>
static LinearProgrammingProblem KleeMintyCube( int n )
{
var objectiveCoeff = new DoubleVector( n );
objectiveCoeff[n - 1] = 1.0;
double eps = 1.0 / 300.0;
var problem = new LinearProgrammingProblem( objectiveCoeff );
problem.AddBounds( 0, 0.0, 1.0 );
for ( int i = 1; i < n; i++ )
{
// Lower bound constraint.
var lbConstraintCoeffs = new DoubleVector( n );
lbConstraintCoeffs[i] = 1.0;
lbConstraintCoeffs[i - 1] = -eps;
problem.AddLowerBoundConstraint( lbConstraintCoeffs, 0.0 );
// Upper bound constraint
var ubConstraintCoeffs = new DoubleVector( n );
ubConstraintCoeffs[i] = 1.0;
ubConstraintCoeffs[i - 1] = eps;
problem.AddUpperBoundConstraint( ubConstraintCoeffs, 1.0 );
}
return problem;
}
}
}
← All NMath Code Examples
