+ /*
+ * Capacity constraint block.
+ */
+ xbt_swag_foreach(cnst, cnst_list) {
+
+ CDEBUG2(surf_sdp_out, "0 %d 1 1 %f", block_num, -(cnst->bound));
+ addentry(constraints, &C, 0, block_num, 1, 1, -(cnst->bound),
+ C.blocks[block_num].blocksize);
+
+ elem_list = &(cnst->element_set);
+ xbt_swag_foreach(elem, elem_list) {
+ if (elem->variable->weight <= 0)
+ break;
+ matno = get_y(K, elem->variable->index);
+ CDEBUG3(surf_sdp_out, "%d %d 1 1 %f", matno, block_num,
+ -(elem->value));
+ addentry(constraints, &C, matno, block_num, 1, 1, -(elem->value),
+ C.blocks[block_num].blocksize);
+
+ }
+ block_num++;
+ }
+
+
+ /*
+ * Positivy constraint blocks.
+ */
+ for (i = 1; i <= pow(2, K); i++) {
+ matno = get_y(K, i);
+ CDEBUG2(surf_sdp_out, "%d %d 1 1 1", matno, block_num);
+ addentry(constraints, &C, matno, block_num, 1, 1, 1.0,
+ C.blocks[block_num].blocksize);
+ block_num++;
+ }
+ /*
+ * Latency constraint blocks.
+ */
+ xbt_swag_foreach(var, var_list) {
+ var->value = 0.0;
+ if (var->weight && var->bound > 0) {
+ matno = get_y(K, var->index);
+ CDEBUG3(surf_sdp_out, "%d %d 1 1 %f", matno, block_num, var->bound);
+ addentry(constraints, &C, matno, block_num, 1, 1, var->bound,
+ C.blocks[block_num].blocksize);
+ }
+ }
+
+ /*
+ * At this point, we'll stop to recognize whether any of the blocks
+ * are "hidden LP blocks" and correct the block type if needed.
+ */
+ for (i = 1; i <= nb_cnsts; i++) {
+ if ((C.blocks[i].blockcategory != DIAG) &&
+ (isdiag[i] == 1) && (C.blocks[i].blocksize > 1)) {
+ /*
+ * We have a hidden diagonal block!
+ */
+
+ blocksz = C.blocks[i].blocksize;
+ tempdiag = (double *) calloc((blocksz + 1), sizeof(double));
+ for (j = 1; j <= blocksz; j++)
+ tempdiag[j] = C.blocks[i].data.mat[ijtok(j, j, blocksz)];
+ free(C.blocks[i].data.mat);
+ C.blocks[i].data.vec = tempdiag;
+ C.blocks[i].blockcategory = DIAG;
+ };
+ };
+
+
+ /*
+ * Next, setup issparse and NULL out all nextbyblock pointers.
+ */
+ p = NULL;
+ for (i = 1; i <= k; i++) {
+ p = constraints[i].blocks;
+ while (p != NULL) {
+ /*
+ * First, set issparse.
+ */
+ if (((p->numentries) > 0.25 * (p->blocksize))
+ && ((p->numentries) > 15)) {
+ p->issparse = 0;
+ } else {
+ p->issparse = 1;
+ };
+
+ if (C.blocks[p->blocknum].blockcategory == DIAG)
+ p->issparse = 1;
+
+ /*
+ * Setup the cross links.
+ */
+
+ p->nextbyblock = NULL;
+ p = p->next;
+ };
+ };
+
+
+ /*
+ * Create cross link reference.
+ */
+ create_cross_link(constraints, nb_var);
+
+
+ /*
+ * Debuging print problem in SDPA format.
+ */
+ if (XBT_LOG_ISENABLED(surf_sdp, xbt_log_priority_debug)) {
+ DEBUG0("Printing SDPA...");
+ tmp = strdup("SURF-PROPORTIONNAL.sdpa");
+ write_prob(tmp, total_block_size, nb_var, C, a, constraints);
+ }
+
+ /*
+ * Initialize parameters.
+ */
+ DEBUG0("Initializing solution...");
+ initsoln(total_block_size, nb_var, C, a, constraints, &X, &y, &Z);
+
+
+
+ /*
+ * Call the solver.
+ */
+ DEBUG0("Calling the solver...");
+ stdout_sav = stdout;
+ stdout = fopen("/dev/null", "w");
+ ret =
+ easy_sdp(total_block_size, nb_var, C, a, constraints, 0.0, &X, &y,
+ &Z, &pobj, &dobj);
+ fclose(stdout);
+ stdout = stdout_sav;
+
+ switch (ret) {
+ case 0:
+ case 1:
+ DEBUG0("SUCCESS The problem is primal infeasible");
+ break;
+
+ case 2:
+ DEBUG0("SUCCESS The problem is dual infeasible");
+ break;
+
+ case 3:
+ DEBUG0
+ ("Partial SUCCESS A solution has been found, but full accuracy was not achieved. One or more of primal infeasibility, dual infeasibility, or relative duality gap are larger than their tolerances, but by a factor of less than 1000.");
+ break;
+
+ case 4:
+ DEBUG0("Failure. Maximum number of iterations reached.");
+ break;
+
+ case 5:
+ DEBUG0("Failure. Stuck at edge of primal feasibility.");
+ break;
+
+ }
+
+ if (XBT_LOG_ISENABLED(surf_sdp, xbt_log_priority_debug)) {
+ tmp = strdup("SURF-PROPORTIONNAL.sol");
+ write_sol(tmp, total_block_size, nb_var, X, y, Z);
+ }
+
+ /*
+ * Write out the solution if necessary.
+ */
+ xbt_swag_foreach(cnst, cnst_list) {
+
+ elem_list = &(cnst->element_set);
+ xbt_swag_foreach(elem, elem_list) {
+ if (elem->variable->weight <= 0)
+ break;
+
+ i = (int) get_y(K, elem->variable->index);
+ elem->variable->value = y[i];
+
+ }
+ }
+
+
+ /*
+ * Free up memory.
+ */
+ free_prob(total_block_size, nb_var, C, a, constraints, X, y, Z);
+
+ free(isdiag);
+ free(tempdiag);
+ free(tmp);
+
+ sys->modified = 0;
+
+ if (XBT_LOG_ISENABLED(surf_sdp, xbt_log_priority_debug)) {
+ lmm_print(sys);
+ }
+
+}
+
+
+/*
+ * Create the cross_link reference in order to have a byblock list.
+ */
+void create_cross_link(struct constraintmatrix *myconstraints, int k)
+{
+
+ int i, j;
+ int blk;
+ struct sparseblock *p;
+ struct sparseblock *q;
+
+ struct sparseblock *prev;
+
+ /*
+ * Now, cross link.
+ */
+ prev = NULL;
+ for (i = 1; i <= k; i++) {
+ p = myconstraints[i].blocks;
+ while (p != NULL) {
+ if (p->nextbyblock == NULL) {
+ blk = p->blocknum;
+
+ /*
+ * link in the remaining blocks.
+ */
+ for (j = i + 1; j <= k; j++) {
+ q = myconstraints[j].blocks;
+
+ while (q != NULL) {
+ if (q->blocknum == p->blocknum) {
+ if (p->nextbyblock == NULL) {
+ p->nextbyblock = q;
+ q->nextbyblock = NULL;
+ prev = q;
+ } else {
+ prev->nextbyblock = q;
+ q->nextbyblock = NULL;
+ prev = q;
+ }
+ break;
+ }
+ q = q->next;
+ }
+ }
+ }
+ p = p->next;
+ }
+ }
+}
+
+
+
+
+void addentry(struct constraintmatrix *constraints,
+ struct blockmatrix *C,
+ int matno,
+ int blkno, int indexi, int indexj, double ent, int blocksize)
+{
+ struct sparseblock *p;
+ struct sparseblock *p_sav;
+
+ p = constraints[matno].blocks;
+
+ if (matno != 0.0) {
+ if (p == NULL) {
+ /*
+ * We haven't yet allocated any blocks.
+ */
+ p = (struct sparseblock *) calloc(1, sizeof(struct sparseblock));
+
+ //two entries because this library ignores indices starting in zerox
+ p->constraintnum = matno;
+ p->blocknum = blkno;
+ p->numentries = 1;
+ p->next = NULL;
+
+ p->entries = calloc(p->numentries + 1, sizeof(double));
+ p->iindices = calloc(p->numentries + 1, sizeof(int));
+ p->jindices = calloc(p->numentries + 1, sizeof(int));
+
+ p->entries[p->numentries] = ent;
+ p->iindices[p->numentries] = indexi;
+ p->jindices[p->numentries] = indexj;
+
+ p->blocksize = blocksize;
+
+ constraints[matno].blocks = p;
+ } else {
+ /*
+ * We have some existing blocks. See whether this block is already
+ * in the chain.
+ */
+ while (p != NULL) {
+ if (p->blocknum == blkno) {
+ /*
+ * Found the right block.
+ */
+ p->constraintnum = matno;
+ p->blocknum = blkno;
+ p->numentries = p->numentries + 1;
+
+ p->entries =
+ realloc(p->entries, (p->numentries + 1) * sizeof(double));
+ p->iindices =
+ realloc(p->iindices, (p->numentries + 1) * sizeof(int));
+ p->jindices =
+ realloc(p->jindices, (p->numentries + 1) * sizeof(int));
+
+ p->entries[p->numentries] = ent;
+ p->iindices[p->numentries] = indexi;
+ p->jindices[p->numentries] = indexj;
+
+ return;
+ }
+ p_sav = p;
+ p = p->next;
+ }
+
+ /*
+ * If we get here, we have a non-empty structure but not the right block
+ * inside hence create a new structure.
+ */
+
+ p = (struct sparseblock *) calloc(1, sizeof(struct sparseblock));
+
+ //two entries because this library ignores indices starting in zerox
+ p->constraintnum = matno;
+ p->blocknum = blkno;
+ p->numentries = 1;
+ p->next = NULL;
+
+ p->entries = calloc(p->numentries + 1, sizeof(double));
+ p->iindices = calloc(p->numentries + 1, sizeof(int));
+ p->jindices = calloc(p->numentries + 1, sizeof(int));
+
+ p->entries[p->numentries] = ent;
+ p->iindices[p->numentries] = indexi;
+ p->jindices[p->numentries] = indexj;
+
+ p->blocksize = blocksize;
+
+ p_sav->next = p;
+ }
+ } else {
+ if (ent != 0.0) {
+ int blksz = C->blocks[blkno].blocksize;
+ if (C->blocks[blkno].blockcategory == DIAG) {
+ C->blocks[blkno].data.vec[indexi] = ent;
+ } else {
+ C->blocks[blkno].data.mat[ijtok(indexi, indexj, blksz)] = ent;
+ C->blocks[blkno].data.mat[ijtok(indexj, indexi, blksz)] = ent;
+ };
+ };
+
+ }