C# Interior Point QP Example

using System;

using CenterSpace.NMath.Core;


namespace CenterSpace.NMath.Examples.CSharp
{
  /// <summary>
  /// A .NET example in C# showing how to use the interior point quadratic
  /// programming solver <c>InteriorPointQPSolver</c> to solve a quadratic
  /// programming problem.
  /// </summary>
  public class InteriorPointQPExample
  {
    static void Main( string[] args )
    {
      // Use the InteriorPointQPSolver class to solve a quadratic programming 
      // problem of the form
      // 0.5*xHx + xc 
      // Our problem will have 4 variables and eight constraints. In addition
      // the variables will be required to be nonnegative.

      // Construct the 4x4 H matrix and c vector for the problem.
      var HMatrixData = new double[,] { {3.0625,  0.0000,  0.0000, 0.0000},
                                        { 0.0000, 0.0405,  0.0000, 0.0000},
                                        { 0.0000,  0.0000, 0.0271, -0.0031},
                                        { 0.0000,  0.0000,  -0.0031, 0.0054} };
      var H = new DoubleMatrix( HMatrixData );
      var c = new DoubleVector( -2671.0, -135.0, -103.0, -19.0 );

      // Construct the quadratic programming problem from H and c.
      var problem = new QuadraticProgrammingProblem( H, c );

      // Our problem will have 8 linear constraints and will be of the form
      // Cx <= b
      // Where C is an 8x4 matrix and b is a vector of length 8.

      // Construct the C matrix and b vector.
      var CmatrixData = new double[,] 
      {  
        {-0.0401, -0.0162, -0.0039, 0.0002},
        {-0.1326, -0.0004, -0.0034, 0.0006},
        {1.5413, 0.0, 0.0, 0.0 },
        {0.0, 0.0203, 0.0, 0.0},
        {0.0, 0.0, 0.0136, -0.0075},
        {0.0, 0.0, -0.0016, 0.0027},
        {0.0160, 0.0004, 0.0005, 0.0002} 
      };
      var C = new DoubleMatrix( CmatrixData );
      var b = new DoubleVector( -92.6, -29.0, 2671.0, 135.0, 103.0, 19.0, 10.0 );

      // Now, add the constraints to the problem.
      for ( int i = 0; i < C.Rows; i++ )
      {
        problem.AddUpperBoundConstraint( C.Row( i ), b[i] );
      }

      // Add variable bounds. Variable values are required to be >= 0.
      for ( int i = 0; i < 4; i++ )
      {
        problem.AddLowerBound( i, 0.0 );
      }

      // Construct the solver.
      var solver = new InteriorPointQPSolver();

      // Set parameters for the solve.
      var solverParams = new InteriorPointQPSolverParams
      {
        KktForm = InteriorPointQPSolverParams.KktFormOption.Blended,
        Tolerance = 1e-6,
        MaxDenseColumnRatio = 0.9,
        PresolveLevel = InteriorPointQPSolverParams.PresolveLevelOption.Full,
        SymbolicOrdering = InteriorPointQPSolverParams.SymbolicOrderingOption.ApproximateMinDegree
      };

      solver.Solve( problem, solverParams );
      Console.WriteLine( "\nSolver Parameters:" );
      Console.WriteLine( solverParams.ToString() );
      Console.WriteLine( "\nResult = " + solver.Result );
      Console.WriteLine( "Optimal x = " + NMathFunctions.Round( solver.OptimalX, 0 ) );
      Console.WriteLine( "Optimal Function value = " + solver.OptimalObjectiveFunctionValue.ToString( "F3" ) );
      Console.WriteLine( "iterations = " + solver.IterationCount );

      Console.WriteLine();
      Console.WriteLine( "Press Enter Key" );
      Console.Read();
    }
  }
}