/***********************************************************************
* Xpress Optimizer Examples
* =========================
*
* Polygon_multimapdelta.cs
*`````````````````````
*
* (c) 2017 Fair Isaac Corporation
*
*
* Polygon example: maximise the area of an N sided polygon
*
* *** Demonstrating using a multi map (R^k->R^l) userfunction calculating its own derivatives ***
*
* Variables:
*
* rho : 0..N-1 ! Distance of vertex from the base point
* theta : 0..N-1 ! Angle from x-axis
*
* Objective:
* (sum (i in 1..N-2) (rho(i)*rho(i-1)*sin(theta(i)-theta(i-1)))) * 0.5
*
* Constraints:
* Vertices in increasing degree order:
* theta(i) >= theta(i-1) +.01 : i = 1..N-2
* Boundary conditions:
* theta(N-1) <= Pi
* 0.1 <= rho(i) <= 1 : i = 0..N-2
* Third side of all triangles <= 1
* rho(i)^2 + rho(j)^2 - rho(i)*rho(j)*2*cos(theta(j)-theta(i)) <= 1 : i in 0..N-3, j in i..N-2
***********************************************************************/
using System;
using System.IO;
using Optimizer;
using XSLPdnet;
namespace XSLPExamples
{
class XSLPExample
{
// userfunction of type "multimapdelta", implementing (x0,x1) -> (sin(x0-x1),cos(x0-x1)) alongside its partials
public int myTrigonometrics(double[] x, double[] deltas, double[] results)
{
int nInput = 2;
/*
* The order of return values for a multimapdelta function is as follows:
*
* Assuming f:R^k->R^l, there will be a total of (k+1)*l return values, loaded to the results array as
* f1(x), diff(f1(x), x1), diff(f1(x), x2) ... diff(f1(x), xk)
* f2(x), diff(f2(x), x1), diff(f2(x), x2) ... diff(f1(x), xk)
* ...
* fl(x), diff(fl(x), x1), diff(fl(x), x2) ... diff(fl(x), xk)
*/
results[0] = Math.Sin(x[0] - x[1]);
results[nInput+1] = Math.Cos(x[0] - x[1]);
if (deltas != null) // Delta may be used as a suggestion for a finite difference step size
// however it also indicates if a partial is requested, saving on effort in case only an evaluation is needed
{
if (deltas[0] != 0.0)
{
results[1+0] = Math.Cos(x[0] - x[1]);
results[nInput + 1 + 1 + 0] = -Math.Sin(x[0] - x[1]);
}
if (deltas[1] != 0.0)
{
results[1+1] = -Math.Cos(x[0] - x[1]);
results[nInput + 1 + 1 + 1] = Math.Sin(x[0] - x[1]);
}
}
return (0);
}
public void RunExp1()
{
/* Example 1 loads a simple problem, and setups the SLP coefficient using a string. */
XPRS.Init("");
try
{
XSLP sProb = new XSLP();
try
{
// Tell Optimizer to call OptimizerMsg whenever a message is output
sProb.XPRS.MessageCallbacks += new MessageCallback(this.OptimizerMsg);
System.Console.WriteLine("####################");
System.Console.WriteLine("# Polygon example #");
System.Console.WriteLine("####################");
System.Console.WriteLine(sProb.getBanner());
// Number of side of the Polygon
int nSide = 5;
/* Create rows */
int nRow = nSide - 2 + (nSide - 1) * (nSide - 2) / 2 + 1;
double[] RHS = new double[nRow];
for (int i = 0; i < nRow; i++) RHS[i] = 0;
char[] RowType = new char[nRow + 1];
string[] RowNames = new string[nRow];
nRow = 0;
// Objective constraint
RowType[nRow] = 'E'; /* OBJEQ */
RowNames[nRow] = "OBJEQ";
nRow++;
// Vertices in increasing degree order
for (int i = 1; i < nSide - 1; i++)
{
RowType[nRow] = 'G'; /* T2T1 .. T4T3 */
RowNames[nRow] = "T" + (i + 1) + "T" + i;
nRow++;
RHS[i] = 0.001;
}
// Third side of all triangles <= 1
for (int i = 1; i < nSide - 1; i++)
{
for (int j = i + 1; j < nSide; j++)
{
RowType[nRow] = 'L';
RHS[nRow] = 1.0;
RowNames[nRow] = "V" + i + "V" + j;
nRow++;
}
}
RowType[nRow] = '\0';
/* Columns and linear coefficients*/
int nCol = (nSide - 1) * 2 + 2;
int nElement = 0;
double[] OBJ = new double[nCol];
double[] Lower = new double[nCol];
double[] Upper = new double[nCol];
for (int i = 0; i < nCol; i++)
{
OBJ[i] = 0; /* objective function */
Lower[i] = 0; /* lower bound normally zero */
Upper[i] = XPRS.PLUSINFINITY; /* upper bound infinity */
}
// Arrays to load linear coefficients
int[] ColStart = new int[nRow + 1];
double[] Element = new double[10 * nCol];
int[] RowIndex = new int[10 * nCol];
string[] ColNames = new string[nCol];
/* OBJX */
nCol = 0;
ColStart[nCol] = nElement;
OBJ[nCol] = 1.0;
ColNames[nCol] = "OBJX";
Lower[nCol++] = XPRS.MINUSINFINITY; /* free column */
Element[nElement] = -1.0;
RowIndex[nElement++] = 0;
/* THETA1 - THETA 4 */
int iRow = 0;
for (int i = 1; i < nSide; i++)
{
ColNames[nCol] = "THETA" + i;
ColStart[nCol++] = nElement;
if (i < nSide - 1)
{
Element[nElement] = -1;
RowIndex[nElement++] = iRow + 1;
}
if (i > 1)
{
Element[nElement] = 1;
RowIndex[nElement++] = iRow;
}
iRow++;
}
Upper[nCol - 1] = 3.1415926;
/* Rho's */
for (int i = 1; i < nSide; i++)
{
Lower[nCol] = 0.01; /* lower bound */
Upper[nCol] = 1;
ColNames[nCol] = "RHO" + i;
ColStart[nCol++] = nElement;
}
ColStart[nCol] = nElement;
sProb.XPRS.LoadLP("Polygon", nCol, nRow, RowType, RHS, null, OBJ, ColStart, null, RowIndex, Element, Lower, Upper);
sProb.XPRS.AddNames(1, RowNames, 0, nRow - 1);
sProb.XPRS.AddNames(2, ColNames, 0, nCol - 1);
/* Add the user function*/
sProb.AddUserFunction("myTrigonometrics", 2 /* #inputs arguments */, 2 /* #output arguments*/, 0, myTrigonometrics);
/* Build up nonlinear coefficients */
/* Area */
string CoefBuffer = "0.5 * ( ";
for (int i = 1; i < nSide - 1; i++)
{
if (i > 1)
{
CoefBuffer += " + ";
}
// Expression using userfunction "myTrigonometrics" in place of the builtin sine function
CoefBuffer += string.Format("RHO{0} * RHO{1} * myTrigonometrics ( THETA{0} , THETA{1} : 1 )", i + 1, i);
}
CoefBuffer += " )";
sProb.ChgStringFormula(0, CoefBuffer);
/* Distances */
for (int i = 1; i < nSide - 1; i++)
{
for (int j = i + 1; j < nSide; j++)
{
// Expression using userfunction "myTrigonometrics" in place of the builtin cosine function
CoefBuffer = string.Format("RHO{0} ^ 2 + RHO{1} ^ 2 - 2 * RHO{0} * RHO{1} * myTrigonometrics ( THETA{0} , THETA{1} : 2 )", j, i);
sProb.ChgStringFormula(iRow, CoefBuffer);
iRow++;
}
}
sProb.SetIntControl(XSLPconst.JACOBIAN, 2); // Note: multmapdelta functions should use symbolic derivatives
sProb.Maxim("");
}
finally
{
/* Cleanup */
sProb.destroyProb();
}
}
finally
{
XPRS.Free();
}
}
public void OptimizerMsg(XPRSprob prob, object data, string message, int len, int msglvl)
{
if (message != null)
{
Console.WriteLine(message);
}
}
public static void Main(string[] args)
{
XSLPExample e = new XSLPExample();
e.RunExp1();
Console.ReadLine();
}
}
}
|
