Changeset 473


Ignore:
Timestamp:
Oct 12, 2004 12:00:13 PM (15 years ago)
Author:
forrest
Message:

Scaling for binv

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/ClpSimplex.cpp

    r472 r473  
    65206520    abort();
    65216521  }
    6522   ClpFactorization * factorization = factorization_;
    65236522  CoinIndexedVector * rowArray0 = rowArray(0);
    65246523  CoinIndexedVector * rowArray1 = rowArray(1);
     
    65296528  columnArray0->clear();
    65306529  columnArray1->clear();
    6531   // put +1 in row
    6532   rowArray1->insert(row,1.0);
    6533   factorization->updateColumnTranspose(rowArray0,rowArray1);
     6530  // put +1 in row
     6531  // But swap if pivot variable was slack as clp stores slack as -1.0
     6532  int pivot = pivotVariable_[row];
     6533  double value;
     6534  // And if scaled then adjust
     6535  if (!rowScale_) {
     6536    if (pivot<numberColumns_)
     6537      value = 1.0;
     6538    else
     6539      value = -1.0;
     6540  } else {
     6541    if (pivot<numberColumns_)
     6542      value = columnScale_[pivot];
     6543    else
     6544      value = -1.0/rowScale_[pivot-numberColumns_];
     6545  }
     6546  rowArray1->insert(row,value);
     6547  factorization_->updateColumnTranspose(rowArray0,rowArray1);
    65346548  // put row of tableau in rowArray1 and columnArray0
    65356549  clpMatrix()->transposeTimes(this,1.0,
    65366550                            rowArray1,columnArray1,columnArray0);
    6537   memcpy(z,columnArray0->denseVector(),
    6538          numberColumns()*sizeof(double));
     6551  if (!rowScale_) {
     6552    memcpy(z,columnArray0->denseVector(),
     6553           numberColumns_*sizeof(double));
     6554  } else {
     6555    double * array = columnArray0->denseVector();
     6556    for (int i=0;i<numberColumns_;i++)
     6557      z[i] = array[i]/columnScale_[i];
     6558  }
    65396559  if (slack) {
    6540     int n = numberRows();
    6541     double * array = rowArray1->denseVector();
    6542     for (int i=0;i<n;i++) {
    6543       // clp stores slacks as -1.0  (Does not seem to matter - basics should be 1.0)
    6544       slack[i] =  array[i];
    6545     }
    6546     //memcpy(slack,rowArray1->denseVector(),
    6547     //   numberRows()*sizeof(double));
     6560    if (!rowScale_) {
     6561      memcpy(slack,rowArray1->denseVector(),
     6562             numberRows_*sizeof(double));
     6563    } else {
     6564      double * array = rowArray1->denseVector();
     6565      for (int i=0;i<numberRows_;i++)
     6566        slack[i] = array[i]*rowScale_[i];
     6567    }
    65486568  }
    65496569  // don't need to clear everything always, but doesn't cost
     
    65756595  rowArray1->clear();
    65766596  // put +1 in row
    6577   rowArray1->insert(row,1.0);
     6597  // But swap if pivot variable was slack as clp stores slack as -1.0
     6598  double value = (pivotVariable_[row]<numberColumns_) ? 1.0 : -1.0;
     6599  // What about scaling ?
     6600  rowArray1->insert(row,value);
    65786601  factorization->updateColumnTranspose(rowArray0,rowArray1);
    65796602  memcpy(z,rowArray1->denseVector(),numberRows()*sizeof(double));
     
    65896612    abort();
    65906613  }
    6591   ClpFactorization * factorization = factorization_;
    65926614  CoinIndexedVector * rowArray0 = rowArray(0);
    65936615  CoinIndexedVector * rowArray1 = rowArray(1);
     
    65966618  // get column of matrix
    65976619#ifndef NDEBUG
    6598   int n = numberColumns();
     6620  int n = numberColumns_+numberRows_;
    65996621  if (col<0||col>=n) {
    66006622    indexError(col,"getBInvACol");
    66016623  }
    66026624#endif
    6603   unpack(rowArray1,col);
    6604   factorization->updateColumn(rowArray0,rowArray1,false);
    6605   memcpy(vec,rowArray1->denseVector(),numberRows()*sizeof(double));
     6625  if (!rowScale_) {
     6626    if (col<numberColumns_) {
     6627      unpack(rowArray1,col);
     6628    } else {
     6629      rowArray1->insert(col-numberColumns_,1.0);
     6630    }
     6631  } else {
     6632    if (col<numberColumns_) {
     6633      unpack(rowArray1,col);
     6634      double multiplier = 1.0/columnScale_[col];
     6635      int number = rowArray1->getNumElements();
     6636      int * index = rowArray1->getIndices();
     6637      double * array = rowArray1->denseVector();
     6638      for (int i=0;i<number;i++) {
     6639        int iRow = index[i];
     6640        // make sure not packed
     6641        assert (array[iRow]);
     6642        array[iRow] *= multiplier;
     6643      }
     6644    } else {
     6645      rowArray1->insert(col-numberColumns_,rowScale_[col-numberColumns_]);
     6646    }
     6647  }
     6648  factorization_->updateColumn(rowArray0,rowArray1,false);
     6649  // But swap if pivot variable was slack as clp stores slack as -1.0
     6650  double * array = rowArray1->denseVector();
     6651  if (!rowScale_) {
     6652    for (int i=0;i<numberRows_;i++) {
     6653      double multiplier = (pivotVariable_[i]<numberColumns_) ? 1.0 : -1.0;
     6654      vec[i] = multiplier * array[i];
     6655    }
     6656  } else {
     6657    for (int i=0;i<numberRows_;i++) {
     6658      int pivot = pivotVariable_[i];
     6659      if (pivot<numberColumns_)
     6660        vec[i] = array[i] * columnScale_[pivot];
     6661      else
     6662        vec[i] = - array[i] / rowScale_[pivot-numberColumns_];
     6663    }
     6664  }
    66066665  rowArray1->clear();
    66076666}
     
    66276686#endif
    66286687  // put +1 in row
    6629   rowArray1->insert(col,1.0);
     6688  // But swap if pivot variable was slack as clp stores slack as -1.0
     6689  double value = (pivotVariable_[col]<numberColumns_) ? 1.0 : -1.0;
     6690  // What about scaling ?
     6691  rowArray1->insert(col,value);
    66306692  factorization->updateColumn(rowArray0,rowArray1,false);
    66316693  memcpy(vec,rowArray1->denseVector(),numberRows()*sizeof(double));
  • trunk/Test/unitTest.cpp

    r472 r473  
    660660    ClpSimplex model;
    661661    model.loadProblem(M, collb, colub, obj, rowlb, rowub);
    662     // For now - without scaling
    663     model.scaling(0);
    664662    model.dual(0,1); // keep factorization
    665663   
     
    668666    double * binvA = (double*) malloc((n_cols+n_rows) * sizeof(double));
    669667   
    670     printf("B-1 A");
     668    printf("B-1 A by row\n");
    671669    for(int i = 0; i < n_rows; i++){
    672670      model.getBInvARow(i, binvA,binvA+n_cols);
    673       printf("\nrow: %d -> ",i);
     671      printf("row: %d -> ",i);
    674672      for(int j=0; j < n_cols+n_rows; j++){
    675673        printf("%g, ", binvA[j]);
    676674      }
    677     }
    678     printf("\n");
    679     free(binvA);
     675      printf("\n");
     676    }
    680677    // See if can re-use factorization
    681678    model.primal(0,3); // keep factorization
     679    // And do by column
     680    printf("B-1 A by column\n");
     681    for(int i = 0; i < n_rows+n_cols; i++){
     682      model.getBInvACol(i, binvA);
     683      printf("column: %d -> ",i);
     684      for(int j=0; j < n_rows; j++){
     685        printf("%g, ", binvA[j]);
     686      }
     687      printf("\n");
     688    }
     689    free(binvA);
    682690    model.dual(0,2); // use factorization
    683691    model.dual(0,2); // hopefully will not use factorization
Note: See TracChangeset for help on using the changeset viewer.