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.

Large Model of Type : RMIQCP

Category : GAMS Model library

Main file : gqapsdp.gms   includes :  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.


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%

if(card(i)*card(j)>20, option limrow=0,limcol=0);

* 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%

$set dou 1
$ifThen %opt%==uv
$  set dov 1
$elseIf %opt%==ua
$  set doa 1
$elseIf %opt%==uav
$  set dov 1
$  set doa 1
$elseIf %opt%==u
$else
$  abort missing test for opt=%opt%
$endIf

Set n 'SDP variables' / 1*%xmax%, 0 /;
Alias(n, np, m);

Set map(n,i,j)  'mapping n to (i,j)';
map(n,i,j) = (ord(n) = (ord(i)-1) * %jmax% + ord(j));

Variables Y(n,n)     'PSDMATRIX'
          sdpz       'objective variable'
;

Equations sdpobj              'objective function'
          sdpassign(i)        'original assigment constraints'
          sdpassign1(i,ip,jp) 'assignment constraints multiplied by x(ip,jp)'
          sdpassign2          'sum of assignment constraints multiplied by themselve'
          sdpcapacity(j)      'location capacity in assigning equipment'
          sdpmatrixident(n)   'identities in SDP matrix'
;

sdpobj.. sdpz =e= sum((i,j,ip,jp,n,np)$(map(n,i,j) and map(np,ip,jp)), (Y(n,np)+Y(np,n))*(unit*flow(i,ip)*dist(j,jp)) / 2.0)
                + sum(map(n,i,j), install(i,j) * (Y(n,'0') + Y('0',n)) / 2.0)
                + eps*sum((n,np),Y(n,np));
                
sdpassign(i).. sum(map(n,i,j), Y(n,'0') + Y('0',n)) / 2.0 =e= 1;

sdpassign1(i,ip,jp).. sum((j,n,np)$(map(n,i,j) and map(np,ip,jp)), (Y(n,np) + Y(np,n)) / 2.0) =e= sum(map(np,ip,jp), Y(np,'0') + Y('0',np)) / 2.0;

sdpassign2.. sum(i, sum((j,jp,n,np)$(map(n,i,j) and map(np,i,jp)), Y(n,np) + Y(np,n)) / 2.0 
  - sum((j,n)$map(n,i,j), Y(n,'0') + Y('0',n)) + 1) =e= 0;

sdpcapacity(j)..  sum(map(n,i,j), a(i)*(Y(n,'0') + Y('0',n))/2.0) =l= b(j);

sdpmatrixident(n)$(not sameas(n,'0')).. Y(n,n) =e= (Y(n,'0') + Y('0',n)) / 2.0;

Y.fx('0','0') = 1.0;

model sdp / sdpobj, sdpassign, sdpcapacity,
$if %dou% == 1 sdpmatrixident
$if %doa% == 1 sdpassign1
$if %dov% == 1 sdpassign2
/;
sdp.dictfile = 1;

option lp = mosek;
solve sdp min sdpz using lp;

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

$if %dou%==1 u(i,j)        = -sum(n$map(n,i,j), sdpmatrixident.m(n));
$if %doa%==1 alpha(ip,i,j) = -sdpassign1.m(ip,i,j);
$if %dov%==1 v             = -sdpassign2.m;

$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), (u(i,j)*1.01)*(sqr(x(i,j)) - x(i,j)))
$if %doa%==1 +  sum((ip,i,j), 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 / defcost, assign, capacity /;

solve gqap using rmiqcp min z;