if (set_member(i, M->linset)) {      poly->EqualityIndex[i]=1;    }    for (j = 1; j <= M->colsize; j++) {      dd_set(poly->A[i-1][j-1], M->matrix[i-1][j-1]);      if (j==1 && dd_Nonzero(M->matrix[i-1][j-1])) poly->homogeneous = dd_FALSE;    }  /*of j*/  }  /*of i*/  dd_DoubleDescription2(poly, horder, err);_L99:  return poly; }void dd_MatrixIntegerFilter(dd_MatrixPtr M){   /* setting an almost integer to the integer. */  dd_rowrange i;  dd_colrange j;  mytype x;  dd_init(x);  for (i=0; i< M->rowsize; i++)    for (j=0; j< M->colsize; j++){       dd_SnapToInteger(x, M->matrix[i][j]);       dd_set(M->matrix[i][j],x);    }  dd_clear(x);}void dd_CopyRay(mytype *a, dd_colrange d_origsize, dd_RayPtr RR,   dd_RepresentationType rep, dd_colindex reducedcol){  long j,j1;  mytype b;  dd_init(b);  for (j = 1; j <= d_origsize; j++){    j1=reducedcol[j];    if (j1>0){      dd_set(a[j-1],RR->Ray[j1-1]);         /* the original column j is mapped to j1, and thus           copy the corresponding component */    } else {      dd_set(a[j-1],dd_purezero);          /* original column is redundant and removed for computation */    }  }  dd_set(b,a[0]);  if (rep==dd_Generator && dd_Nonzero(b)){    dd_set(a[0],dd_one);    for (j = 2; j <= d_origsize; j++)       dd_div(a[j-1],a[j-1],b);    /* normalization for generators */  }  dd_clear(b);}void dd_WriteRay(FILE *f, dd_colrange d_origsize, dd_RayPtr RR, dd_RepresentationType rep, dd_colindex reducedcol){  dd_colrange j;  static dd_colrange d_last=0;  static dd_Arow a;  if (d_last< d_origsize){    if (d_last>0) free(a);    dd_InitializeArow(d_origsize+1, &a);    d_last=d_origsize+1;  } dd_CopyRay(a, d_origsize, RR, rep, reducedcol);  for (j = 0; j < d_origsize; j++) dd_WriteNumber(f, a[j]);  fprintf(f, "\n");}void dd_WriteArow(FILE *f, dd_Arow a, dd_colrange d){  dd_colrange j;  for (j = 0; j < d; j++) dd_WriteNumber(f, a[j]);  fprintf(f, "\n");}void dd_WriteAmatrix(FILE *f, dd_Amatrix A, long rowmax, long colmax){  long i,j;  if (A==NULL){    fprintf(f, "WriteAmatrix: The requested matrix is empty\n");    goto _L99;  }  fprintf(f, "begin\n");#if defined GMPRATIONAL  fprintf(f, " %ld %ld rational\n",rowmax, colmax);#else  fprintf(f, " %ld %ld real\n",rowmax, colmax);#endif  for (i=1; i <= rowmax; i++) {    for (j=1; j <= colmax; j++) {      dd_WriteNumber(f, A[i-1][j-1]);    }    fprintf(f,"\n");  }  fprintf(f, "end\n");_L99:;}void dd_WriteBmatrix(FILE *f, dd_colrange d_size, dd_Bmatrix B){  dd_colrange j1, j2;  if (B==NULL){    fprintf(f, "WriteBmatrix: The requested matrix is empty\n");    goto _L99;  }  for (j1 = 0; j1 < d_size; j1++) {    for (j2 = 0; j2 < d_size; j2++) {      dd_WriteNumber(f, B[j1][j2]);    }  /*of j2*/    putc('\n', f);  }  /*of j1*/  putc('\n', f);_L99:;}void dd_WriteSetFamily(FILE *f, dd_SetFamilyPtr F){  dd_bigrange i;  if (F==NULL){    fprintf(f, "WriteSetFamily: The requested family is empty\n");    goto _L99;  }  fprintf(f,"begin\n");  fprintf(f,"  %ld    %ld\n", F->famsize, F->setsize);  for (i=0; i<F->famsize; i++) {    fprintf(f, " %ld %ld : ", i+1, set_card(F->set[i]));    set_fwrite(f, F->set[i]);  }  fprintf(f,"end\n");_L99:;}void dd_WriteSetFamilyCompressed(FILE *f, dd_SetFamilyPtr F){  dd_bigrange i,card;  if (F==NULL){    fprintf(f, "WriteSetFamily: The requested family is empty\n");    goto _L99;  }  fprintf(f,"begin\n");  fprintf(f,"  %ld    %ld\n", F->famsize, F->setsize);  for (i=0; i<F->famsize; i++) {    card=set_card(F->set[i]);    if (F->setsize - card >= card){      fprintf(f, " %ld %ld : ", i+1, card);      set_fwrite(f, F->set[i]);    } else {      fprintf(f, " %ld %ld : ", i+1, -card);      set_fwrite_compl(f, F->set[i]);    }  }  fprintf(f,"end\n");_L99:;}void dd_WriteMatrix(FILE *f, dd_MatrixPtr M){  dd_rowrange i, linsize;  if (M==NULL){    fprintf(f, "WriteMatrix: The requested matrix is empty\n");    goto _L99;  }  switch (M->representation) {    case dd_Inequality:      fprintf(f, "H-representation\n"); break;    case dd_Generator:      fprintf(f, "V-representation\n"); break;    case dd_Unspecified:      break;  }  linsize=set_card(M->linset);  if (linsize>0) {    fprintf(f, "linearity %ld ", linsize);    for (i=1; i<=M->rowsize; i++)       if (set_member(i, M->linset)) fprintf(f, " %ld", i);    fprintf(f, "\n");  }  dd_WriteAmatrix(f, M->matrix, M->rowsize, M->colsize);  if (M->objective!=dd_LPnone){    if (M->objective==dd_LPmax)      fprintf(f, "maximize\n");    else      fprintf(f, "minimize\n");          dd_WriteArow(f, M->rowvec, M->colsize);  }_L99:;}void dd_WriteLP(FILE *f, dd_LPPtr lp){  if (lp==NULL){    fprintf(f, "WriteLP: The requested lp is empty\n");    goto _L99;  }  fprintf(f, "H-representation\n");  dd_WriteAmatrix(f, lp->A, (lp->m)-1, lp->d);  if (lp->objective!=dd_LPnone){    if (lp->objective==dd_LPmax)      fprintf(f, "maximize\n");    else      fprintf(f, "minimize\n");          dd_WriteArow(f, lp->A[lp->objrow-1], lp->d);  }_L99:;}void dd_SnapToInteger(mytype y, mytype x){ /* this is broken.  It does nothing. */   dd_set(y,x);}void dd_WriteReal(FILE *f, mytype x){  long ix1,ix2,ix;  double ax;  ax=dd_get_d(x);    ix1= (long) (fabs(ax) * 10000. + 0.5);  ix2= (long) (fabs(ax) + 0.5);  ix2= ix2*10000;  if ( ix1 == ix2) {    if (dd_Positive(x)) {      ix = (long)(ax + 0.5);    } else {      ix = (long)(-ax + 0.5);      ix = -ix;    }    fprintf(f, " %2ld", ix);  } else    fprintf(f, " % .9E",ax);}void dd_WriteNumber(FILE *f, mytype x){#if defined GMPRATIONAL  mpz_t zn,zd;  mpz_init(zn); mpz_init(zd);  mpq_canonicalize(x);  mpq_get_num(zn,x);  mpq_get_den(zd,x);  fprintf(f," ");  if (mpz_sgn(zn)==0){    fprintf(f,"0");  } else if (mpz_cmp_ui(zd,1U)==0){    mpz_out_str(f,10,zn);  } else {    mpz_out_str(f,10,zn);fprintf(f,"/");mpz_out_str(f,10,zd);  }  mpz_clear(zn); mpz_clear(zd);#else  dd_WriteReal(f, x);#endif}void dd_WriteIncidence(FILE *f, dd_PolyhedraPtr poly){  dd_SetFamilyPtr I;  switch (poly->representation) {    case dd_Inequality:       fprintf(f, "ecd_file: Incidence of generators and inequalities\n");      break;    case dd_Generator:       fprintf(f, "icd_file: Incidence of inequalities and generators\n");      break;    default:      break;  }  I=dd_CopyIncidence(poly);  dd_WriteSetFamilyCompressed(f,I);  dd_FreeSetFamily(I);}void dd_WriteAdjacency(FILE *f, dd_PolyhedraPtr poly){  dd_SetFamilyPtr A;  switch (poly->representation) {    case dd_Inequality:       fprintf(f, "ead_file: Adjacency of generators\n");      break;    case dd_Generator:       fprintf(f, "iad_file: Adjacency of inequalities\n");      break;    default:      break;  }  A=dd_CopyAdjacency(poly);  dd_WriteSetFamilyCompressed(f,A);  dd_FreeSetFamily(A);}void dd_ComputeAinc(dd_PolyhedraPtr poly){/* This generates the input incidence array poly->Ainc, and   two sets: poly->Ared, poly->Adom. */  dd_bigrange k;  dd_rowrange i,m1;  dd_colrange j;  dd_boolean redundant;  dd_MatrixPtr M=NULL;  mytype sum,temp;  dd_init(sum); dd_init(temp);  if (poly->AincGenerated==dd_TRUE) goto _L99;  M=dd_CopyOutput(poly);  poly->n=M->rowsize;  m1=poly->m1;     /* this number is same as poly->m, except when      poly is given by nonhomogeneous inequalty:      !(poly->homogeneous) && poly->representation==Inequality,      it is poly->m+1.   See dd_ConeDataLoad.   */  poly->Ainc=(set_type*)calloc(m1, sizeof(set_type));  for(i=1; i<=m1; i++) set_initialize(&(poly->Ainc[i-1]),poly->n);  set_initialize(&(poly->Ared), m1);   set_initialize(&(poly->Adom), m1);   for (k=1; k<=poly->n; k++){    for (i=1; i<=poly->m; i++){      dd_set(sum,dd_purezero);      for (j=1; j<=poly->d; j++){        dd_mul(temp,poly->A[i-1][j-1],M->matrix[k-1][j-1]);        dd_add(sum,sum,temp);      }      if (dd_EqualToZero(sum)) {        set_addelem(poly->Ainc[i-1], k);      }    }    if (!(poly->homogeneous) && poly->representation==dd_Inequality){      if (dd_EqualToZero(M->matrix[k-1][0])) {        set_addelem(poly->Ainc[m1-1], k);  /* added infinity inequality (1,0,0,...,0) */      }    }  }  for (i=1; i<=m1; i++){    if (set_card(poly->Ainc[i-1])==M->rowsize){      set_addelem(poly->Adom, i);    }    }  for (i=m1; i>=1; i--){    if (set_card(poly->Ainc[i-1])==0){      redundant=dd_TRUE;      set_addelem(poly->Ared, i);    }else {      redundant=dd_FALSE;      for (k=1; k<=m1; k++) {        if (k!=i && !set_member(k, poly->Ared)  && !set_member(k, poly->Adom) &&             set_subset(poly->Ainc[i-1], poly->Ainc[k-1])){          if (!redundant){            redundant=dd_TRUE;          }          set_addelem(poly->Ared, i);        }      }    }  }  dd_FreeMatrix(M);  poly->AincGenerated=dd_TRUE;_L99:;  dd_clear(sum);  dd_clear(temp);}dd_boolean dd_InputAdjacentQ(dd_PolyhedraPtr poly,   dd_rowrange i1, dd_rowrange i2)/* Before calling this function, RedundantSet must be    a set of row indices whose removal results in a minimal   nonredundant system to represent the input polyhedron,   DominantSet must be the set of row indices which are   active at every extreme points/rays.*/{  dd_boolean adj=dd_TRUE;  dd_rowrange i;  static set_type common;  static long lastn=0;  if (poly->AincGenerated==dd_FALSE) dd_ComputeAinc(poly);  if (lastn!=poly->n){    if (lastn >0) set_free(common);    set_initialize(&common, poly->n);    lastn=poly->n;  }  if (set_member(i1, poly->Ared) || set_member(i2, poly->Ared)){    adj=dd_FALSE;    goto _L99;  }  if (set_member(i1, poly->Adom) || set_member(i2, poly->Adom)){      /* dominant inequality is considered adjacencent to all others. */    adj=dd_TRUE;    goto _L99;  }  set_int(common, poly->Ainc[i1-1], poly->Ainc[i2-1]);  i=0;  while (i<poly->m1 && adj==dd_TRUE){     i++;     if (i!=i1 && i!=i2 && !set_member(i, poly->Ared) &&        !set_member(i, poly->Adom) && set_subset(common,poly->Ainc[i-1])){      adj=dd_FALSE;    }  }_L99:;  return adj;} void dd_WriteInputIncidence(FILE *f, dd_PolyhedraPtr poly){  dd_SetFamilyPtr I;  if (poly->AincGenerated==dd_FALSE) dd_ComputeAinc(poly);  switch (poly->representation) {  case dd_Inequality:    fprintf(f,"icd_file: Incidence of inequalities and generators\n");    break;  case dd_Generator:   fprintf(f,"ecd_file: Incidence of generators and inequalities\n");    break;  default:    break;  }  I=dd_CopyInputIncidence(poly);  dd_WriteSetFamilyCompressed(f,I);  dd_FreeSetFamily(I);}void dd_WriteInputAdjacency(FILE *f, dd_PolyhedraPtr poly){  dd_SetFamilyPtr A;  if (poly->AincGenerated==dd_FALSE){    dd_ComputeAinc(poly);  }  switch (poly->representation) {  case dd_Inequality:    fprintf(f, "iad_file: Adjacency of inequalities\n");    break;  case dd_Generator:    fprintf(f, "ead_file: Adjacency of generators\n");    break;  default:    break;  }  A=dd_CopyInputAdjacency(poly);  dd_WriteSetFamilyCompressed(f,A);  dd_FreeSetFamily(A);}void dd_WriteProgramDescription(FILE *f){  fprintf(f, "* cddlib: a double description library:%s\n", dd_DDVERSION);  fprintf(f, "* compiled for %s arithmetic.\n", dd_ARITHMETIC);  fprintf(f,"* %s\n",dd_COPYRIGHT);}void dd_WriteTimes(FILE *f, time_t starttime, time_t endtime){  long ptime,ptime_sec,ptime_minu, ptime_hour;  /*    ptime=difftime(endtime,starttime);    This function is ANSI standard, but not available sometime */  ptime=(endtime - starttime);      /* This is to replace the line above, but it may not give     correct time in seconds */   ptime_hour=ptime/3600;  ptime_minu=(ptime-ptime_hour*3600)/60;  ptime_sec=ptime%60;  fprintf(f,"* Computation started at %s",asctime(localtime(&starttime)));  fprintf(f,"*             ended   at %s",asctime(localtime(&endtime)));  fprintf(f,"* Total processor time = %ld seconds\n",ptime);  fprintf(f,"*                      = %ld h %ld m %ld s\n",ptime_hour, ptime_minu, ptime_sec);}void dd_WriteDDTimes(FILE *f, dd_PolyhedraPtr poly){  dd_WriteTimes(f,poly->child->starttime,poly->child->endtime);}