FICO Xpress Optimization Examples Repository
 FICO Optimization Community FICO Xpress Optimization Home

Catenary - Solving a QCQP

Description
This model finds the shape of a hanging chain by minimizing its potential energy.

Source Files
By clicking on a file name, a preview is opened at the bottom of this page.

xbcatena.cs

/****************************************************************
BCL Example Problems
====================

file xbcatena.cs

QCQP test problem
Based on AMPL model catenary.mod
(Source: http://www.orfe.princeton.edu/~rvdb/ampl/nlmodels/ )

This model finds the shape of a hanging chain.

(c) 2008 Fair Isaac Corporation
authors: S.Heipcke, D.Brett, June 2008
****************************************************************/

using System;
using System.Text;
using System.IO;
using Optimizer;
using BCL;

namespace Examples
{

public class TestBCatena
{
const int N = 100;	               // Number of chainlinks
const int L = 1;	               // Difference in x-coordinates of endlinks

const double H = 2.0*L/N;	       // Length of each link

public static void Main()
{
XPRB.init();
int i;
XPRBvar[] x = new XPRBvar[N + 1];
XPRBvar[] y = new XPRBvar[N + 1];       // x-/y-coordinates of endpoints of chainlinks
XPRBexpr qe;
XPRBctr cobj;
XPRBprob prob = new XPRBprob("catenary");  // Initialize a new problem in BCL
XPRSprob xprob;

prob.setDictionarySize(BCLconstant.XPRB_DICT_NAMES,0);

/**** VARIABLES ****/
for(i=0;i<=N;i++)
x[i] = prob.newVar("x("+i+")", BCLconstant.XPRB_PL, -BCLconstant.XPRB_INFINITY,
BCLconstant.XPRB_INFINITY);
for(i=0;i<=N;i++)
y[i] = prob.newVar("y("+i+")", BCLconstant.XPRB_PL, -BCLconstant.XPRB_INFINITY,
BCLconstant.XPRB_INFINITY);
// Left anchor
x[0].fix(0);  y[0].fix(0);
// Right anchor
x[N].fix(L);  y[N].fix(0);

/****OBJECTIVE****/
/*  sum(j in 1..N) (y(j-1)+y(j))/2  */
qe = new XPRBexpr();
for(i=1;i<=N;i++) qe+= y[i-1]+y[i];
cobj = prob.newCtr("Obj", qe*0.5 );
prob.setObj(cobj);                     /* Set objective function */

/**** CONSTRAINTS ****/
/*  forall(j in 1..N) (x(j)-x(j-1))^2+(y(j)-y(j-1))^2 <= H^2  */
for(i=1;i<=N;i++)
prob.newCtr("Link_"+i, (x[i]-x[i-1]).sqr() + (y[i]-y[i-1]).sqr() <= H*H);

/****SOLVING + OUTPUT****/
prob.setSense(BCLconstant.XPRB_MINIM);   // Choose the sense of optimization

/* Problem printing and matrix output: */
/*
prob.print();
prob.exportProb(BCLconstant.XPRB_MPS, "caternary");
prob.exportProb(BCLconstant.XPRB_LP, "caternary");
*/

//Disable the convexivity check
xprob = prob.getXPRSprob();
xprob.IfCheckConvexity = 0;

//Solve the problem
prob.lpOptimize();

System.Console.WriteLine("Solution: " + prob.getObjVal());

for(i=0;i<=N;i++)
System.Console.WriteLine(i + ": " + x[i].getSol() + ", " + y[i].getSol());

return;
}
}
}