gqapsdp.gms : SDP Convexifications of the Generalized Quadratic Assignment Problem

Description

This model solves the Generalized Quadratic Assignment Problem
(GQAP) using different convexification methods. The GQAP describes
a broad class of quadratic integer programming problems, wherein
I pair-wise related equipments are assigned to N locations constrained
by the locations' ability to accommodate them.

Cplex has a build-in convexification method for MIQPs when the
quadratic terms involve binary variables only. The other methods
require the solution of a semidefinite program.  The communication
with the SDP solver is done through ASCII files.


References

  • Plateau M.C., Reformulations quadratiques convexes pour la programmation quadratique en variables 0-1. PhD thesis, Conservatoire National des Arts et Metiers, CEDRIC, 2006.
  • Guignard, M, Extension to Plateau Convexification Method for Nonconvex Quadratic 0-1 Programs. Tech. rep., The Wharton School, 2008.

Large Model of Type : RMIQCP


Category : GAMS Model library


Main file : gqapsdp.gms   includes :  runcsdp.inc  pb16x7.inc

$title SDP Convexifications of the Generalized Quadratic Assignment Problem (GQAPSDP,SEQ=339)

$onText
This model solves the Generalized Quadratic Assignment Problem
(GQAP) using different convexification methods. The GQAP describes
a broad class of quadratic integer programming problems, wherein
I pair-wise related equipments are assigned to N locations constrained
by the locations' ability to accommodate them.

Cplex has a build-in convexification method for MIQPs when the
quadratic terms involve binary variables only. The other methods
require the solution of a semidefinite program.  The communication
with the SDP solver is done through ASCII files.


Plateau M.C., Reformulations quadratiques convexes pour la
programmation quadratique en variables 0-1. PhD thesis,
Conservatoire National des Arts et Metiers, CEDRIC, 2006.

Guignard, M., Extension to Plateau Convexification Method for
Nonconvex Quadratic 0-1 Programs. Tech. rep., The Wharton School, 2008.

Keywords: relaxed mixed integer quadratic constraint programming, quadratic
          assignment problem, convexification methods, semidefinite programming
$offText

Set
   i 'equipments'
   j 'locations';

Alias (i,ip), (j,jp);

Parameter
   install(i,j) 'installation cost'
   flow(i,ip)   'volume of exchange'
   unit         'flow units'
   dist(j,jp)   'distance'
   a(i)         'size'
   b(j)         'capacity';

$if not set instance $set instance pb16x7.inc
$include %instance%

* Convexification options: cplex, u, uv, au, uav
$if not set opt    $set opt u
$if '%opt%'==u     $goto optok
$if '%opt%'==uv    $goto optok
$if '%opt%'==ua    $goto optok
$if '%opt%'==uav   $goto optok
$if '%opt%'==cplex $goto optok
$abort --opt=%opt% but should be u|uv|ua|uav|cplex
$label optok

$set dou 0
$set dov 0
$set doa 0
$if '%opt%'==cplex $goto nosdp

$eval imax card(i)
$eval jmax card(j)
$eval xmax %imax%*%jmax%
$eval dmax %imax%*%imax%*%jmax%

$set dou 1
$ifThen %opt%==u
$  set more ''
$elseIf %opt%==uv
$  set more ',v'
$  set dov 1
$elseIf %opt%==ua
$  set more ',d1*d%dmax%'
$  set doa 1
$elseIf %opt%==uav
$  set more ',d1*d%dmax%,v'
$  set dov 1
$  set doa 1
$else
$  abort missing test for opt=%opt%
$endIf

Set
   n      'SDP variables'       / #j,z,x1*x%xmax% /
   nx(n)  'original x'          /      x1*x%xmax% /
   nj(n)  'slack for capacity'  / #j              /
   m      'sdp constraints'     / #n,#i%more%     /
   mIx(m) 'identity x'          / #nx /
   mA(m)  'assignment'          / #i  /
   mC(m)  'capacity'            / #j  /
$ifThen %doa%==1
   mD(m)                        / d1*d%dmax%   /
   dmap(md,*,nx) 'd iij map'    / #md:(#i.#nx) /
$endIf
   xmap(*,*,j) 'x ij map'       / #nx:(#i.#j) /
   upper(n,n)  'diag and upper triangle'
   lower(n,n)  'lower triangle'
;

Alias (n,np,npp), (nx,nxp);
upper(n,np) = ord(n) <= ord(np);
lower(n,np) = ord(n) >  ord(np);

Parameter
   F(m,n,np) 'SDP constraint data'
   F0(n,np)  'SDP obj'
   c(m)      'SDP RHS data'
   qq(n,np)  'Q matrix'
   mLE(m)    'SDP constraints with =l=';

qq(nx,nxp) = sum((i,j,ip,jp)$(xmap(nx,i,j)*xmap(nxp,ip,jp)), unit*flow(i,ip)*dist(j,jp));
qq(nx,nxp) = (qq(nx,nxp) + qq(nxp,nx))/2;

* Objective
F0('z',nx)        = -sum(xmap(nx,i,j), install(i,j))/2;
F0(upper(nx,nxp)) = -qq(nx,nxp);

* Identity for x
F(mIx,nx,nx)$sameas(mIx,nx)  =  2;
F(mIx,'z',nx)$sameas(mIx,nx) = -1;

* Assignment
F(mA,'z',nx) = sum(xmap(nx,i,j)$sameas(mA,i), 1);
c(mA)        = 2;

* Capacity =l= constraint
F(mC,'z',nx)              = sum(xmap(nx,i,j)$sameas(mC,j), a(i));
F(mC,nj,nj)$sameas(mC,nj) = 1;
c(mC)                     = 2*sum(j$sameas(mC,j), b(j));
mLE(mC)                   = yes;

$ifThen %doa%==1
* Alpha
   F(mD,'z',nx)        = -sum(dmap(mD,mA,nx), 1);
   loop((dmap(mD,mA,nx),xmap(nxp,mA,j)), F(mD,nx,nxp) = 2;);
   F(mD,upper(nx,nxp)) = (F(mD,nxp,nx) + F(mD,nx,nxp))/2;
   F(mD,lower) = 0;
$endIf

$ifThen %dov%==1
* v
   F('v','z',nx)        = -1;
   F('v',upper(nx,nxp)) =  sum((i,j,jp)$(xmap(nx,i,j) and xmap(nxp,i,jp)), 1);
   c('v')               = -card(i);
$endIf

* z = 1
F('z','z','z') = 1;
c('z')         = 1;

Parameter
   xvec(m)       'SDP solution'
   u(i,j)        'dual of X_ijij = xij'
   alpha(ip,i,j) 'dual of assignment constraints'
   v             'dual of |Ax-b|';

$if not set skip $set skip 0
$ifThen %skip%==0
   execute_unload 'csdpin.gdx' n, m, mLE, c, F, F0;
   execute 'gams runcsdp.inc lo=%gams.lo% --strict=1';
   abort$errorLevel 'Problems running RunCSDP. Check listing file for details.';
$endIf
execute_load 'csdpout.gdx', xvec;

$if %dou%==1 u(i,j)        = sum(xmap(mIx,i,j), xvec(mIx));
$if %doa%==1 alpha(ip,i,j) = sum((dmap(mD,ip,nx),xmap(nx,i,j)), xvec(mD));
$if %dov%==1 v             = xvec('v');

$label nosdp
Variable
   x(i,j) 'equipment i is assigned to location j'
   z      'total costs';

Binary Variable x;

Equation
   defcost     'objective function'
   assign(i)   'assign each equipment to a location'
   capacity(j) 'location capacity in assigning equipment';

defcost.. z =e= sum((i,j,ip,jp), x(i,j)*(unit*flow(i,ip)*dist(j,jp))*x(ip,jp))
             +  sum((i,j), install(i,j)*x(i,j))
$if %dou%==1 +  sum((i,j), 2*(u(i,j)*1.01)*(sqr(x(i,j)) - x(i,j)))
$if %doa%==1 +  sum((ip,i,j), 2*alpha(ip,i,j)*x(i,j)*(sum(jp,x(ip,jp)) - 1))
$if %dov%==1 +  v*sum(i, sqr(sum(j,x(i,j)) - 1))
;

assign(i)..   sum(j, x(i,j))      =e= 1;

capacity(j).. sum(i, a(i)*x(i,j)) =l= b(j);

Model gqap / all /;

if(card(i)*card(j) > 20, option limRow = 0, limCol = 0;);

solve gqap using rmiqcp min z;