/* Copyright (C) 2000 Dale Schuurmans, Finnegan Southey */ /* This work is licensed under the Gnu General Public License (see gpl.txt). */ /* Additive updates, using our standard hinge penalty */ #include #include #include #include #include #include "cop.h" #include "satfront.h" #define swapint(x,y) {tempint=x; x=y; y=tempint;} #define swapreal(x,y) {tempreal=x; x=y; y=tempreal;} typedef double REAL; int tempint; REAL tempreal; int n, m, N, M; int *X, *A, *nCs, **Cs, **Ccoefs; int *B, *nXs, **Xs, **Xcoefs; int MaxCoef, MinB, MaxB; int steps, dsteps, considered, r; int T, R, W; int shift, S, COR; REAL alpha, rho, oneminusrho; int *incrsat, *decrsat, *incrvio, *decrvio; int *xvreducers, *minxs, *wherexvreduce, *cviol, *wherecviol, *V; int **numup, **numdn; REAL *Y, *primpenalty, *dualpenalty; REAL **Yup, **Ydn, **diff1ppen, **diff2ppen, *diffL; void preSolve(SATParams *params) { int minV, maxV, minval, maxval, c, s; double walkprob; T = params->intParams[0]; if (T == 0) T = 500000; R = params->intParams[1]; if (R == 0) R = 1; COR = params->intParams[2]; if (COR == 0) COR = 100; alpha = params->doubleParams[0]; if (alpha == 0) alpha = 1; /* rho = params->doubleParams[1]; if (rho == 0) rho = 0.995; oneminusrho = 1 - rho; */ walkprob = params->doubleParams[2]; if (walkprob == 0) walkprob = 0.002; W = RAND_MAX * walkprob; /*alpha = 1.0 + (REAL)par1 * n / (REAL)m;*/ if (MaxCoef > 1) { printf("Error: this implementation assumes {-1,0,+1} valued "); printf("coefficients\n"); exit(1); } minV = INT_MAX; maxV = INT_MIN; for (c = 0; c < m; c++) { if ( (nXs[c] + B[c]) % 2 ) { printf("Error: can only handle even constraint violation values\n"); exit(1); } minval = (-nXs[c] - B[c]) / 2; maxval = (nXs[c] - B[c]) / 2; if (minval < minV) minV = minval; if (maxval > maxV) maxV = maxval; } S = maxV - minV + 1; shift = -minV; incrsat = (int *) malloc(n * sizeof(int)); decrsat = (int *) malloc(n * sizeof(int)); incrvio = (int *) malloc(n * sizeof(int)); decrvio = (int *) malloc(n * sizeof(int)); xvreducers = (int *) malloc(n * sizeof(int)); minxs = (int *) malloc(n * sizeof(int)); wherexvreduce = (int *) malloc(n * sizeof(int)); cviol = (int *) malloc(m * sizeof(int)); wherecviol = (int *) malloc(m * sizeof(int)); V = (int *) malloc(m * sizeof(int)); Y = (REAL *) malloc(m * sizeof(REAL)); diff1ppen = (REAL **) malloc(2 * sizeof(REAL *)); diff2ppen = (REAL **) malloc(2 * sizeof(REAL *)); diffL = (REAL *) malloc(n * sizeof(REAL)); primpenalty = (REAL *) malloc(S * sizeof(REAL)); diff1ppen[0] = (REAL *) malloc(S * sizeof(REAL)); diff1ppen[1] = (REAL *) malloc(S * sizeof(REAL)); diff2ppen[0] = (REAL *) malloc(S * sizeof(REAL)); diff2ppen[1] = (REAL *) malloc(S * sizeof(REAL)); dualpenalty = (REAL *) malloc(S * sizeof(REAL)); numup = (int **) malloc(S * sizeof(int *)); numdn = (int **) malloc(S * sizeof(int *)); Yup = (REAL **) malloc(S * sizeof(REAL *)); Ydn = (REAL **) malloc(S * sizeof(REAL *)); for (s = 0; s < S; s++) { numup[s] = (int *) malloc(n * sizeof(int)); numdn[s] = (int *) malloc(n * sizeof(int)); Yup[s] = (REAL *) malloc(n * sizeof(REAL)); Ydn[s] = (REAL *) malloc(n * sizeof(REAL)); /*** Dual penalty ***/ dualpenalty[s] = (s <= shift) ? 0.0 : 1.0; } /*** Primal penalty ***/ for (s = 0; s < S; s++) primpenalty[s] = (s <= shift) ? 0.0 : 1.0; for (s = 0; s < S; s++) { if (s > 0) diff1ppen[1][s] = primpenalty[s-1] - primpenalty[s]; if (s < S-1) diff1ppen[0][s] = primpenalty[s+1] - primpenalty[s]; if (s > 0 && s < S - 1) { diff2ppen[0][s] = 2*primpenalty[s]-primpenalty[s+1]-primpenalty[s-1]; diff2ppen[1][s] = -diff2ppen[0][s]; } } } void runSolver(SATParams *params, SATResults *results) { REAL Yval = 1.0 / m; int x, c, s, i, j, t, d, nmxs, xflip, nxvr, numvio = 0; int Scurr, Sflip, agreexflip, agreex; REAL mindiff, ytotal; steps = considered = dsteps = 0; for (r = 0; r < R; r++) { /*******************************************/ for (x = 0; x < n; x++) X[x] = (mrand48() < 0) ? -1 : 1; steps++; considered++; numvio = 0; for (c = 0; c < m; c++) { V[c] = -B[c]; for (i = 0; i < nXs[c]; i++) V[c] += X[Xs[c][i]] == Xcoefs[c][i] ? 1 : -1; if (V[c] > 0) { cviol[numvio] = c; wherecviol[c] = numvio; numvio++; } V[c] /= 2; Y[c] = Yval; Scurr = V[c] + shift; } d = 1; nxvr = 0; for (x = 0; x < n; x++) { incrsat[x] = decrsat[x] = incrvio[x] = decrvio[x] = 0; diffL[x] = 0.0; for (s = 0; s < S; s++) { numup[s][x] = numdn[s][x] = 0; Yup[s][x] = Ydn[s][x] = 0.0; } for (j = 0; j < nCs[x]; j++) { c = Cs[x][j]; Scurr = V[c] + shift; if (X[x] == Ccoefs[x][j]) { numdn[Scurr][x]++; Ydn[Scurr][x] += Y[c]; diffL[x] += Y[c] * diff1ppen[1][Scurr]; if (V[c] > 0) { if (!decrvio[x]++) { xvreducers[nxvr] = x; wherexvreduce[x] = nxvr; nxvr++; } } else decrsat[x]++; } else { numup[Scurr][x]++; Yup[Scurr][x] += Y[c]; diffL[x] += Y[c] * diff1ppen[0][Scurr]; if (V[c] > 0) incrvio[x]++; else incrsat[x]++; } } } for (t = 0; t < T; t++) { /***************************/ if (numvio == 0) break; if (d % COR == 0) { /*** CORRECTION ****/ ytotal = 0.0; for (c = 0; c < m; c++) ytotal += Y[c]; for (c = 0; c < m; c++) { Y[c] /= ytotal; } d = 1; for (x = 0; x < n; x++) { diffL[x] = 0.0; for (s = 0; s < S; s++) Yup[s][x] = Ydn[s][x] = 0.0; for (j = 0; j < nCs[x]; j++) { c = Cs[x][j]; Scurr = V[c] + shift; if (X[x]==Ccoefs[x][j]) { Ydn[Scurr][x] += Y[c]; diffL[x] += Y[c] * diff1ppen[1][Scurr]; } else { Yup[Scurr][x] += Y[c]; diffL[x] += Y[c] * diff1ppen[0][Scurr]; } } } } mindiff = DBL_MAX; /*** Candidates ****/ nmxs = 0; for (i = 0; i < nxvr; i++) { x = xvreducers[i]; if (diffL[x] < mindiff) { mindiff = diffL[x]; minxs[nmxs = 0] = x; nmxs++; } else if (diffL[x] == mindiff) minxs[nmxs++] = x; } if (mindiff >= 0.0 && random() < W) { nmxs = 1; x = minxs[0] = lrand48()%n; diffL[x] = 0.0; if (decrvio[x] == 0) for (s = 0; s < S; s++) diffL[x] += Yup[s][x] * diff1ppen[0][s] + Ydn[s][x] * diff1ppen[1][s]; mindiff = -1.0; } if (mindiff >= 0.0) { /*** DUAL update ***/ dsteps++; for (c = 0; c < m; c++) Y[c] += alpha * dualpenalty[V[c]+shift]; for (i = 0; i < nxvr; i++) diffL[xvreducers[i]] = 0.0; for (s = 0; s < S; s++) { for (x = 0; x < n; x++) { Yup[s][x] += numup[s][x] * alpha * dualpenalty[s]; Ydn[s][x] += numdn[s][x] * alpha * dualpenalty[s]; if (decrvio[x] > 0) diffL[x] += Yup[s][x] * diff1ppen[0][s] + Ydn[s][x] * diff1ppen[1][s]; } } d++; } else { /*** PRIMAL update */ if (nmxs == 1) xflip = minxs[0]; else xflip = minxs[lrand48()%nmxs]; X[xflip] = -X[xflip]; diffL[xflip] = -diffL[xflip]; for(s = 0; s < S-1; s++) { swapint(numup[s][xflip], numdn[s+1][xflip]); swapreal(Yup[s][xflip], Ydn[s+1][xflip]); } swapint(incrsat[xflip],decrsat[xflip]); swapint(incrvio[xflip],decrvio[xflip]); if (incrvio[xflip] > 0 && decrvio[xflip] == 0) { nxvr--; xvreducers[wherexvreduce[xflip]] = xvreducers[nxvr]; wherexvreduce[xvreducers[nxvr]] = wherexvreduce[xflip]; } else if (incrvio[xflip] == 0 && decrvio[xflip] > 0) { xvreducers[nxvr] = xflip; wherexvreduce[xflip] = nxvr; nxvr++; } steps++; considered += n; for (j = 0; j < nCs[xflip]; j++) { c = Cs[xflip][j]; agreexflip = X[xflip] == Ccoefs[xflip][j]; Scurr = V[c] + shift; Sflip = Scurr + (agreexflip ? 1 : -1); if (Sflip > shift && Scurr <= shift) { cviol[numvio] = c; wherecviol[c] = numvio; numvio++; for (i = 0; i < nXs[c]; i++) { x = Xs[c][i]; if (X[x] == Xcoefs[c][i]) { if (!decrvio[x]++) { xvreducers[nxvr] = x; wherexvreduce[x] = nxvr; nxvr++; diffL[x] = 0.0; for (s = 0; s < S; s++) diffL[x] += Yup[s][x] * diff1ppen[0][s] + Ydn[s][x] * diff1ppen[1][s]; } decrsat[x]--; } else { incrvio[x]++; incrsat[x]--; } } } else if (Sflip <= shift && Scurr > shift) { numvio--; cviol[wherecviol[c]] = cviol[numvio]; wherecviol[cviol[numvio]] = wherecviol[c]; for (i = 0; i < nXs[c]; i++) { x = Xs[c][i]; if (X[x] == Xcoefs[c][i]) { if (!--decrvio[x]) { nxvr--; xvreducers[wherexvreduce[x]] = xvreducers[nxvr]; wherexvreduce[xvreducers[nxvr]] = wherexvreduce[x]; } decrsat[x]++; } else { incrvio[x]--; incrsat[x]++; } } } for (i = 0; i < nXs[c]; i++) if ((x = Xs[c][i]) != xflip) { if ((agreex = X[x]) == Xcoefs[c][i]) { numdn[Scurr][x]--; numdn[Sflip][x]++; Ydn[Scurr][x] -= Y[c]; Ydn[Sflip][x] += Y[c]; } else { numup[Scurr][x]--; numup[Sflip][x]++; Yup[Scurr][x] -= Y[c]; Yup[Sflip][x] += Y[c]; } if (decrvio[x] > 0) diffL[x] += Y[c] * ( agreex == agreexflip ? diff2ppen[0][Scurr] : diff2ppen[1][Sflip]); }; V[c] = Sflip - shift; } } /*******************/ } /***************************/ if (numvio == 0) break; } /*******************************************/ if (numvio > 0) r--; results->numFlips = steps; results->numReweights = dsteps; results->numRestarts = r + 1; results->numSatisfied = m - numvio; results->isSatisfied = numvio == 0; #include "copysoln.c" } void postSolve(SATParams *params) { int s; free(Y); free(V); free(primpenalty); free(dualpenalty); free(diffL); free(minxs); free(incrsat); free(decrsat); free(incrvio); free(decrvio); free(xvreducers); free(wherexvreduce); free(cviol); free(wherecviol); free(diff1ppen[0]); free(diff1ppen[1]); free(diff1ppen); free(diff2ppen[0]); free(diff2ppen[1]); free(diff2ppen); for (s = 0; s < S; s++) { free(numup[s]); free(numdn[s]); free(Yup[s]); free(Ydn[s]); } free(numup); free(numdn); free(Yup); free(Ydn); } SATParams *getSATParams() { SATParams *params; params = malloc(sizeof(SATParams)); params->numIntParams = 3; params->intParams = malloc(params->numIntParams * sizeof(int)); params->intParamNames = malloc(params->numIntParams * sizeof(char*)); params->intParamDescs = malloc(params->numIntParams * sizeof(char*)); params->intParamNames[0] = "-mf"; params->intParamDescs[0] = "max flips before restarting"; params->intParams[0] = 0; params->intParamNames[1] = "-mr"; params->intParamDescs[1] = "max restarts before aborting"; params->intParams[1] = 0; params->intParamNames[2] = "-cp"; params->intParamDescs[2] = "number of reweights between corrections"; params->intParams[2] = 0; params->numDoubleParams = 3; params->doubleParams = malloc(params->numDoubleParams * sizeof(double)); params->doubleParamNames = malloc(params->numDoubleParams * sizeof(char*)); params->doubleParamDescs = malloc(params->numDoubleParams * sizeof(char*)); params->doubleParamNames[0] = "-alpha"; params->doubleParamDescs[0] = "target post-flood gradient"; params->doubleParams[0] = 0; params->doubleParamNames[1] = "-rho"; params->doubleParamDescs[1] = "rate of weight shrinkage to mean for sat clauses"; params->doubleParams[1] = 0; params->doubleParamNames[2] = "-noise"; params->doubleParamDescs[2] = "probability of random walk when stuck"; params->doubleParams[2] = 0; params->numStringParams = 0; return params; } char *getSolverIdentity() { return "\$RCSfile: asg.c,v \$ \$Revision: 1.8 \$ \$Date: 2001/11/23 06:35:04 \$"; }