Quadratic Programming in IronPython QuickStart Sample

Illustrates how to solve optimization problems a quadratic objective function and linear constraints using classes in the Numerics.NET.Optimization namespace in IronPython.

View this sample in: C# Visual Basic F#

```Python
import numerics

from System import Array

from Extreme.Mathematics import *

# The quadratic programming classes reside in their own namespace.
from Extreme.Mathematics.Optimization import *

# Illustrates solving quadratic programming problems
# using the classes in the Extreme.Mathematics.Optimization
# namespace of Numerics.NET.

# This QuickStart Sample illustrates the quadratic programming
# functionality by solving a portfolio optimization problem.

# Portfolio optimization is a common application of QP.
# For a collection of assets, the goal is to minimize
# the risk (variance of the return) while achieving
# a minimal return for a set maximum amount invested.

# The variables are the amounts invested in each asset.
# The quadratic term is the covariance matrix of the assets.
# THere is no linear term in this case.

# There are three ways to create a Quadratic Program.

# The first is in terms of matrices. The coefficients
# of the constraints and the quadratic terms are supplied 
# as matrices. The cost vector, right-hand side and 
# constraints on the variables are supplied as vectors.

# The linear term in the objective function:
c = Vector.CreateConstant(4, 0.0)
# The quaratic term in the objective function:
R = Matrix.CreateSymmetric(4, Array[float]([ \
    0.08,-0.05,-0.05,-0.05, \
    -0.05, 0.16,-0.02,-0.02, \
    -0.05,-0.02, 0.35, 0.06, \
    -0.05,-0.02, 0.06, 0.35 ]), MatrixTriangle.Upper)
# The coefficients of the constraints:
A = Matrix([ [1, 1, 1, 1], [-0.05, 0.2, -0.15, -0.30 ]])
# The right-hand sides of the constraints:
b = Vector([10000, -1000])

# We're now ready to call the constructor.
# The last parameter specifies the number of equality
# constraints.
qp1 = QuadraticProgram(c, R, A, b, 0)

# Now we can call the Solve method to run the Revised
# Simplex algorithm:
x = qp1.Solve()
print "Solution: {0:.1f}".format(x)
# The optimal value is returned by the Extremum property:
print "Optimal value:   {0:.1f}".format(qp1.OptimalValue)


# The second way to create a Quadratic Program is by constructing
# it by hand. We start with an 'empty' quadratic program.
qp2 = QuadraticProgram()

# Next, we add two variables: we specify the name, the cost, # and optionally the lower and upper bound.
qp2.AddVariable("X1", 0.0)
qp2.AddVariable("X2", 0.0)
qp2.AddVariable("X3", 0.0)
qp2.AddVariable("X4", 0.0)

# Next, we add constraints. Constraints also have a name.
# We also specify the coefficients of the variables, # the lower bound and the upper bound.
qp2.AddLinearConstraint("C1", Vector.Create(1.0, 1.0, 1.0, 1.0), ConstraintType.LessThanOrEqual, 10000)
qp2.AddLinearConstraint("C2", Vector.Create(0.05, -0.2, 0.15, 0.3), ConstraintType.GreaterThanOrEqual, 1000)
# If a constraint is a simple equality or inequality constraint, # you can supply a QuadraticProgramConstraintType value and the
# right-hand side of the constraint.

# Quadratic terms must be set individually.
# Each combination appears at most once.
qp2.SetQuadraticCoefficient("X1", "X1", 0.08)
qp2.SetQuadraticCoefficient("X1", "X2", -0.05 * 2)
qp2.SetQuadraticCoefficient("X1", "X3", -0.05 * 2)
qp2.SetQuadraticCoefficient("X1", "X4", -0.05 * 2)
qp2.SetQuadraticCoefficient("X2", "X2", 0.16)
qp2.SetQuadraticCoefficient("X2", "X3", -0.02 * 2)
qp2.SetQuadraticCoefficient("X2", "X4", -0.02 * 2)
qp2.SetQuadraticCoefficient("X3", "X3", 0.35)
qp2.SetQuadraticCoefficient("X3", "X4", 0.06 * 2)
qp2.SetQuadraticCoefficient("X4", "X4", 0.35)

# We can now solve the quadratic program:
x = qp2.Solve()
print "Solution: {0:.1f}".format(x)
print "Optimal value:   {0:.1f}".format(qp2.OptimalValue)


# Finally, we can create a quadratic program from an MPS file.
# The MPS format is a standard format.
qp3 = MpsReader.ReadQuadraticProgram(r"..\data\portfolio.qps")
# We can go straight to solving the quadratic program:
x = qp3.Solve()
print "Solution: {0:.1f}".format(x)
print "Optimal value:   {0:.1f}".format(qp3.OptimalValue)

```