00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00014
00015
00016
00017
00018
00020
00021 #include "ModkSeparator.defs"
00022 #include "ModkSeparator.hh"
00023 #include <time.h>
00024 #include <stdio.h>
00025 #include <math.h>
00026 #include <stdlib.h>
00027
00028
00029
00030
00031
00032
00033
00034 ModkSeparator::ModkSeparator(int modfac, int nConstraints, int nVariables,
00035 int* nNonZeros, int** NonZeroIndices,
00036 double** NonZeroCoeff, double* rhs,
00037 double* lpSolution):
00038 ordGen(false),
00039 constraintsGen(false),
00040 congruenceSystemGen(0),
00041 colsGen(false),
00042 generatedCutsGen(false)
00043 {
00044 K = modfac;
00045
00046
00047 if( K == 2 ){
00048 Ord = new int[2];
00049 ordGen = true;
00050 Ord[0]=0;
00051 Ord[1]=1;
00052 }
00053 else if( K == 3 ){
00054 Ord = new int[3];
00055 ordGen = true;
00056 Ord[0]=0;
00057 Ord[1]=1;
00058 Ord[2]=2;
00059 }
00060 else if( K == 5 ){
00061 Ord = new int[5];
00062 ordGen = true;
00063 Ord[0]=0;
00064 Ord[1]=1;
00065 Ord[2]=3;
00066 Ord[3]=2;
00067 Ord[4]=4;
00068 }
00069 else if( K == 7 ){
00070 Ord = new int[7];
00071 ordGen = true;
00072 Ord[0] = 0;
00073 Ord[1] = 1;
00074 Ord[2] = 4;
00075 Ord[3] = 5;
00076 Ord[4] = 2;
00077 Ord[5] = 3;
00078 Ord[6] = 6;
00079 }
00080 else{
00081 printf("%d is not a supported modulo factor\n", K);
00082 exit(1);
00083 }
00084
00085
00086 if(IMPORTANTVARIABLES){
00087 int curVar = 0;
00088
00089
00090
00091 int* count = new int[nVariables];
00092 for(int i=0; i<nVariables; i++)
00093 count[i] = 0;
00094
00095 for(int i=0; i<nConstraints; i++)
00096 for(int j=0; j<nNonZeros[i]; j++)
00097 count[NonZeroIndices[i][j]]++;
00098
00099 for(int i=0; i<nVariables; i++){
00100 if(count[i] > count[curVar]){
00101 curVar = i;
00102 }
00103 }
00104
00105
00106 constraintsGen = true;
00107 NonZeroI = new int*[nConstraints];
00108 NonZeroC = new double*[nConstraints];
00109 nNonZ = new int[nConstraints];
00110 RHS = new double[nConstraints];
00111
00112 int maxCon = nVariables + MAXGAP;
00113 int* taken = new int[nConstraints];
00114 for(int i=0; i<nConstraints; i++)
00115 taken[i] = 0;
00116 bool doBreak = false;
00117 bool doContinue = false;
00118 nCon = 0;
00119 int nConOld = nCon;
00120
00121 while(1){
00122 nConOld = nCon;
00123
00124 for(int i=0; i<nConstraints; i++){
00125 if(nCon == maxCon){
00126
00127 doBreak = true;
00128 break;
00129 }
00130 if(taken[i] == 0){
00131 for(int j=0; j<nNonZeros[i]; j++){
00132 if(NonZeroIndices[i][j] == curVar){
00133
00134
00135 taken[i] = 1;
00136 nNonZ[nCon] = nNonZeros[i];
00137 NonZeroI[nCon] = NonZeroIndices[i];
00138 NonZeroC[nCon] = NonZeroCoeff[i];
00139 RHS[nCon] = rhs[i];
00140 nCon++;
00141 }
00142 }
00143 }
00144 }
00145 if(doBreak)
00146 break;
00147
00148 if(nConOld == nCon){
00149 doContinue = false;
00150 for(int i=0; i<nVariables; i++){
00151 if(count[i] != -1){
00152
00153 curVar = i;
00154 count[curVar] = -1;
00155 doContinue = true;
00156 break;
00157 }
00158 }
00159 if(doContinue)
00160 continue;
00161 else{
00162
00163 break;
00164 }
00165 }
00166
00167 count[curVar] = -1;
00168 for(int i=0; i<nVariables; i++){
00169 if(count[i] != -1)
00170 count[i] = 0;
00171 }
00172
00173 for(int i=0; i<nCon; i++){
00174 for(int j=0; j<nNonZ[i]; j++){
00175 if(count[NonZeroI[i][j]] >= 0)
00176 count[NonZeroI[i][j]]++;
00177 }
00178 }
00179
00180
00181 for(int i=0; i<nVariables; i++){
00182 if(count[i] > count[curVar]){
00183 curVar = i;
00184 }
00185 }
00186
00187
00188 }
00189
00190
00191 delete [] count;
00192 delete [] taken;
00193
00194 }
00195 else{
00196
00197
00198 NonZeroI = NonZeroIndices;
00199 NonZeroC = NonZeroCoeff;
00200 RHS = rhs;
00201 nNonZ = nNonZeros;
00202 nCon = nConstraints;
00203 }
00204
00205 nVar = nVariables;
00206 lpSol = lpSolution;
00207
00208 nRows = nVar + 1;
00209 nCols = nCon + 1;
00210
00211
00212
00213 CongruenceSystem = new int*[nRows];
00214 congruenceSystemGen = nRows;
00215
00216 for(int i=0; i<=nVar; i++)
00217 CongruenceSystem[i] = new int[nCols];
00218
00219
00220
00221 int entry=0;
00222
00223 for(int i=0; i<nCon; i++){
00224 for(int j=0; j<nVar; j++){
00225 entry = 0;
00226 for(int k=0; k<(int) nNonZ[i]; k++){
00227 if(NonZeroI[i][k] == j){
00228 entry = (K + (((int) NonZeroC[i][k]) % K)) % K;
00229 }
00230 }
00231 CongruenceSystem[j][i] = entry;
00232 }
00233 }
00234
00235 for(int i=0; i<nRows-1; i++)
00236 CongruenceSystem[i][nCols-1] = 0;
00237 CongruenceSystem[nRows-1][nCols-1] = K-1;
00238
00239 for(int i=0; i<nCon; i++){
00240 CongruenceSystem[nVar][i] = (K + (((int) RHS[i]) % K)) % K;
00241 }
00242 }
00243
00244
00245
00246
00247
00248 ModkSeparator::~ModkSeparator()
00249 {
00250
00251 if (ordGen) {
00252 delete [] Ord;
00253 }
00254 if (constraintsGen) {
00255
00256
00257 delete [] NonZeroI;
00258 delete [] NonZeroC;
00259 delete [] nNonZ;
00260 delete [] RHS;
00261 }
00262 if (congruenceSystemGen) {
00263 for (int i=0; i<congruenceSystemGen; i++) {
00264 delete [] CongruenceSystem[i];
00265
00266 }
00267 delete [] CongruenceSystem;
00268 }
00269 if (colsGen) {
00270 delete [] Cols;
00271 }
00272 if (generatedCutsGen) {
00273 for (int i=0; i<nGen; i++) {
00274 delete [] GeneratedCuts[i];
00275 }
00276 delete [] GeneratedCuts;
00277 }
00278 }
00279
00280
00281
00282
00283
00284 int ModkSeparator::separate(int nMaxCuts)
00285 {
00286 time_t begin = time(NULL);
00287
00288 bool doBreak = false;
00289 int* swap;
00290 int swap2=0;
00291 int idSize= ((nRows) < (nCols)) ? (nRows) : (nCols);
00292 int factor = 0;
00293
00294 Cols = new int[nCols];
00295 colsGen = true;
00296 for(int i=0; i<nCols; i++)
00297 Cols[i] = i;
00298
00302 if(RANDOMIZE){
00303 srand(time(NULL));
00304 for (int i=0; i<nCols; i++) {
00305 int r = rand() % nCols;
00306 int temp = Cols[i]; Cols[i] = Cols[r]; Cols[r] = temp;
00307 }
00308 }
00309
00310 swap = CongruenceSystem[0];
00311 CongruenceSystem[0] = CongruenceSystem[nRows-1];
00312 CongruenceSystem[nRows-1] = swap;
00313
00314
00315
00316
00317 int scalFac = 0;
00318
00319 for(int i=0; (i<nRows) && (i<nCols-1); i++){
00320
00321 if( CongruenceSystem[i][Cols[i]] > 1 ){
00322 scalFac = Ord[CongruenceSystem[i][Cols[i]]];
00323
00324
00325 for(int j=i; j<nCols; j++){
00326 CongruenceSystem[i][Cols[j]] = (scalFac*CongruenceSystem[i][Cols[j]]) % K;
00327 }
00328 }
00329
00330
00331
00332 if( CongruenceSystem[i][Cols[i]] == 0 ){
00333
00334 for(int k=i+1; k<nRows; k++){
00335
00336 if( CongruenceSystem[k][Cols[i]] != 0 ){
00337
00338 swap = CongruenceSystem[i];
00339 CongruenceSystem[i] = CongruenceSystem[k];
00340 CongruenceSystem[k] = swap;
00341 break;
00342
00343 }
00344 }
00345 }
00346
00347
00348 if( CongruenceSystem[i][Cols[i]] == 0 ){
00349
00350
00351 for(int k=i+1; k<nCols-1; k++){
00352
00353 if( CongruenceSystem[i][Cols[k]] != 0 ){
00354 swap2 = Cols[k];
00355
00356 Cols[k] = Cols[i];
00357 Cols[i] = swap2;
00358 break;
00359
00360 }
00361 }
00362 }
00363
00364
00365 if( CongruenceSystem[i][Cols[i]] == 0 ){
00366
00367 doBreak = false;
00368 for(int k=i+1; k<nRows; k++){
00369
00370
00371 for(int l=i+1; l<nCols-1; l++){
00372 if( CongruenceSystem[k][Cols[l]] != 0){
00373
00374 swap = CongruenceSystem[i];
00375
00376 CongruenceSystem[i] = CongruenceSystem[k];
00377 CongruenceSystem[k] = swap;
00378 swap2 = Cols[i];
00379
00380 Cols[i] = Cols[l];
00381 Cols[l] = swap2;
00382 doBreak = true;
00383
00384 break;
00385 }
00386
00387 if(doBreak) break;
00388 }
00389 }
00390 }
00391
00392
00393
00394 if( CongruenceSystem[i][Cols[i]] > 1 ){
00395
00396 scalFac = Ord[CongruenceSystem[i][Cols[i]]];
00397 for(int j=i; j<nCols; j++){
00398
00399 CongruenceSystem[i][Cols[j]] = (scalFac*CongruenceSystem[i][Cols[j]]) % K;
00400 }
00401 }
00402
00403
00404
00405 if( CongruenceSystem[i][Cols[i]] == 1 ){
00406
00407
00408 for(int k=i+1; k<nRows; k++){
00409 if( CongruenceSystem[k][Cols[i]] != 0 ){
00410
00411 factor = CongruenceSystem[k][Cols[i]];
00412 for(int j=i; j<nCols; j++){
00413
00414 CongruenceSystem[k][Cols[j]] = (K + ((CongruenceSystem[k][Cols[j]] -
00415 factor*CongruenceSystem[i][Cols[j]])%K))%K;
00416 }
00417 }
00418 }
00419 }
00420 else{
00421 idSize = i;
00422 break;
00423 }
00424 }
00425
00426
00427
00428
00429 for(int i=idSize-1; i>=0; i--){
00430 for(int k=i-1; k>=0; k--){
00431 if( CongruenceSystem[k][Cols[i]] != 0 ){
00432 factor = CongruenceSystem[k][Cols[i]];
00433 for(int j=i; j<nCols; j++){
00434
00435 CongruenceSystem[k][Cols[j]] = (K + ((CongruenceSystem[k][Cols[j]] -
00436 factor*CongruenceSystem[i][Cols[j]])%K))%K;
00437 }
00438 }
00439 }
00440 }
00441
00442
00443
00444 for(int i=idSize; i<nRows; i++)
00445 if(CongruenceSystem[i][nCols-1] != 0){
00446
00447
00448 return 0;
00449 }
00450
00451
00452
00453 GeneratedCuts = new int*[nMaxCuts];
00454 generatedCutsGen = true;
00455
00456 nGen = 0;
00457
00458
00459
00460
00461
00462 int checksum = 0;
00463 bool feasible = false;
00464
00465
00466
00467
00468 GeneratedCuts[nGen] = new int[nCols-1];
00469 generatedCutsGen = true;
00470 for(int l=0; l<nCols-1; l++){
00471 if(l<idSize)
00472 GeneratedCuts[nGen][Cols[l]] =
00473 CongruenceSystem[l][Cols[nCols-1]];
00474 else
00475 GeneratedCuts[nGen][Cols[l]] = 0;
00476 }
00477 nGen++;
00478
00479
00480
00481 bool cutGood = true;
00482
00483
00484 int curSum;
00485
00486 time_t START = time(NULL);
00487 time_t END;
00488 for(int i=1; i<pow(2,nCols-idSize); i++){
00489 END = time(NULL);
00490 if(nGen == nMaxCuts || difftime(END,START) > TIMEOUT)
00491 break;
00492
00493 GeneratedCuts[nGen] = new int[nCols-1];
00494 for(int l=0; l<nCols-1; l++)
00495 GeneratedCuts[nGen][l] = 0;
00496
00497 for(int j=0; j<idSize; j++){
00498
00499 curSum = 0;
00500 for(int k=idSize; k<nCols-1; k++){
00501 if( (i & (1 << (k-idSize))) ){
00502 curSum = curSum + CongruenceSystem[j][Cols[k]];
00503 GeneratedCuts[nGen][Cols[k]] = 1;
00504 }
00505 }
00506
00507 GeneratedCuts[nGen][Cols[j]] = (K + ((CongruenceSystem[j][Cols[nCols-1]] - curSum)%K )) %K;
00508 }
00509 cutGood = true;
00510
00516 if(MAXDENSITY){
00517 int nNZ = 0;
00518
00519 double* Cut = new double[nVar];
00520
00521 for(int j=0; j<nVar; j++)
00522 Cut[j] = 0;
00523
00524 double Multiplier = 1.0/((double) K);
00525
00526 for(int l=0; l<nCon; l++){
00527 if(GeneratedCuts[nGen][l] != 0){
00528 for(int j=0; j<nNonZ[l]; j++){
00529 Cut[NonZeroI[l][j]] = Cut[NonZeroI[l][j]] +
00530 Multiplier*GeneratedCuts[nGen][l]*NonZeroC[l][j];
00531 }
00532 }
00533 }
00534
00535 for(int l=0; l<nVar; l++){
00536 if(Cut[l] != 0)
00537 nNZ++;
00538 }
00539
00540 if(nNZ > MAXDENSITY){
00541 cutGood = false;
00542 }
00543
00544
00545 delete [] Cut;
00546 }
00547
00548 if(cutGood)
00549 nGen++;
00550 }
00551
00552
00553
00554 time_t end = time(NULL);
00555
00556
00557
00558 nCutsLeft = nGen;
00559 return nGen;
00560
00561 }
00562
00563
00564
00565
00566
00567 int* ModkSeparator::getCut() {
00568
00569 int nNonZerosCut=0;
00570 int k=-1;
00571 double* Cut = new double[nVar];
00572 double CutRHS = 0;
00573 int* sparseCut;
00574
00575 double entry=0;
00576
00577 if(nCutsLeft){
00578 nCutsLeft--;
00579 int cI = nGen-(nCutsLeft+1);
00580
00581
00582
00583 for(int j=0; j<nVar; j++)
00584 Cut[j] = 0;
00585
00586 double Multiplier = 1.0/((double) K);
00587
00588
00589 for(int i=0; i<nCon; i++){
00590
00591 if(GeneratedCuts[cI][i] != 0){
00592 for(int j=0; j<nNonZ[i]; j++){
00593 Cut[NonZeroI[i][j]] = Cut[NonZeroI[i][j]] +
00594 Multiplier*GeneratedCuts[cI][i]*NonZeroC[i][j];
00595 }
00596 CutRHS = CutRHS + Multiplier*GeneratedCuts[cI][i]*RHS[i];
00597 }
00598 }
00599
00600
00601 double lhs=0;
00602
00603
00604
00605 for(int i=0; i<nVar; i++)
00606 if(Cut[i] != 0)
00607 nNonZerosCut++;
00608 sparseCut = new int[(nNonZerosCut*2)+2];
00609 sparseCut[0] = nNonZerosCut;
00610 for(int i=0; i<nVar; i++){
00611 if(Cut[i] != 0){
00612 lhs = lhs + ((int) Cut[i])*lpSol[i];
00613 k++;
00614 sparseCut[k+1] = i;
00615 sparseCut[k+1+nNonZerosCut] = (int) Cut[i];
00616 }
00617 }
00618 sparseCut[(2*nNonZerosCut)+1] = (int) floor(CutRHS);
00619
00620
00621 delete [] Cut;
00622
00623 return sparseCut;
00624 }
00625 else {
00626
00627 delete [] Cut;
00628
00629 return NULL;
00630 }
00631 }
00632
00633
00634
00635
00636
00637
00638
00639
00640
00641
00642
00643 unsigned int ModkSeparator::hweight(unsigned int w)
00644 {
00645 unsigned int res = (w & 0x55555555) + ((w >> 1) & 0x55555555);
00646 res = (res & 0x33333333) + ((res >> 2) & 0x33333333);
00647 res = (res & 0x0F0F0F0F) + ((res >> 4) & 0x0F0F0F0F);
00648 res = (res & 0x00FF00FF) + ((res >> 8) & 0x00FF00FF);
00649 return (res & 0x0000FFFF) + ((res >> 16) & 0x0000FFFF);
00650 }
00651
00652