Changeset 1057


Ignore:
Timestamp:
Jul 31, 2007 5:33:14 AM (12 years ago)
Author:
forrest
Message:

nonlinear plus fixes to perturbation

Location:
trunk/Clp/src
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/Clp/src/ClpSimplex.cpp

    r1056 r1057  
    15821582    //printf("ret after %d pivots\n",factorization_->pivots());
    15831583    return 1;
    1584   } else {
    1585     if (forceFactorization_>0&&
    1586         factorization_->pivots()==forceFactorization_) {
    1587       // relax
    1588       forceFactorization_ = (3+5*forceFactorization_)/4;
    1589       if (forceFactorization_>factorization_->maximumPivots())
    1590         forceFactorization_ = -1; //off
     1584  } else if (forceFactorization_>0&&
     1585             factorization_->pivots()==forceFactorization_) {
     1586    // relax
     1587    forceFactorization_ = (3+5*forceFactorization_)/4;
     1588    if (forceFactorization_>factorization_->maximumPivots())
     1589      forceFactorization_ = -1; //off
     1590    return 1;
     1591  } else if (numberIterations_>1000+10*(numberRows_+(numberColumns_>>2))) {
     1592    double random = CoinDrand48();
     1593    int maxNumber = (forceFactorization_<0) ? maximumPivots : CoinMin(forceFactorization_,maximumPivots);
     1594    if (factorization_->pivots()>=random*maxNumber) {
     1595      return 1;
     1596    } else if (numberIterations_>1000000+10*(numberRows_+(numberColumns_>>2))&&
     1597               numberIterations_<1001000+10*(numberRows_+(numberColumns_>>2))) {
    15911598      return 1;
    15921599    } else {
     
    15941601      return 0;
    15951602    }
     1603  } else {
     1604    // carry on iterating
     1605    return 0;
    15961606  }
    15971607}
  • trunk/Clp/src/ClpSimplexNonlinear.cpp

    r1056 r1057  
    33513351  int numberArtificials=0;
    33523352  // We could use nonlinearcost to do segments - maybe later
    3353 #define SEGMENTS 3 
     3353#define SEGMENTS 3
     3354  // Penalties may be adjusted by duals
     3355  // Both these should be modified depending on problem
     3356  // Possibly start with big bounds
     3357  double penalties[]={1.0e-2,1.0,1.0e9};
     3358  //double bounds[] = {1.0e-2,1.0,COIN_DBL_MAX};
     3359  double bounds[] = {1.0e-1,1.0e2,COIN_DBL_MAX};
    33543360  // see how many extra we need
    33553361  CoinBigIndex numberExtra=0;
     
    33853391  }
    33863392  int numberColumns2 = numberColumns_;
    3387   // Penalties may be adjusted by duals
    3388   // Both these should be modified depending on problem
    3389   double penalties[]={1.0e-2,1.0,1.0e9};
    3390   double bounds[] = {1.0e-2,1.0,COIN_DBL_MAX};
    33913393  if (numberArtificials) {
    33923394    numberArtificials *= SEGMENTS;
     
    34303432    delete [] addUpper;
    34313433    delete [] addCost;
     3434    //    newModel.primal(1);
    34323435  }
    34333436  // find nonlinear columns
     
    35263529  newModel.primal(1);
    35273530  // still infeasible
    3528   if (newModel.numberPrimalInfeasibilities()) {
     3531  if (newModel.problemStatus()==1) {
    35293532    delete [] listNonLinearColumn;
    35303533    return 0;
     3534  } else if (newModel.problemStatus()==2) {
     3535    // unbounded - add bounds
     3536    double * columnLower = newModel.columnLower();
     3537    double * columnUpper = newModel.columnUpper();
     3538    for (int i=0;i<numberColumns_;i++) {
     3539      columnLower[i]= CoinMax(-1.0e8,columnLower[i]);
     3540      columnUpper[i]= CoinMin(1.0e8,columnUpper[i]);
     3541    }
     3542    newModel.primal(1);
    35313543  }
    35323544  int numberRows = newModel.numberRows();
     
    35503562    double upper = columnUpper[iColumn];
    35513563    double lower = columnLower[iColumn];
    3552     if (solution[iColumn]<trueLower[jNon])
    3553       solution[iColumn]=trueLower[jNon];
    3554     else if (solution[iColumn]>trueUpper[jNon])
    3555       solution[iColumn]=trueUpper[jNon];
     3564    if (solution[iColumn]<lower)
     3565      solution[iColumn]=lower;
     3566    else if (solution[iColumn]>upper)
     3567      solution[iColumn]=upper;
    35563568#if 0
    35573569    double large = CoinMax(1000.0,10.0*fabs(solution[iColumn]));
     
    36983710    }
    36993711    double infPenalty=0.0;
     3712    // This penalty is for target drop
     3713    double infPenalty2=0.0;
    37003714    for (iConstraint=0;iConstraint<numberConstraints;iConstraint++) {
    37013715      ClpConstraint * constraint = constraints[iConstraint];
     
    37073721               newModel.primalRowSolution()[iRow],
    37083722               dualValue,offset);
     3723      double movement = newModel.primalRowSolution()[iRow];
     3724      movement = fabs((movement-functionValue)*dualValue);
     3725      infPenalty2 += movement;
    37093726      double infeasibility=0.0;
    37103727      if (functionValue<rowLower_[iRow]-1.0e-5) {
     
    38463863           numberColumns2+numberRows);
    38473864   
    3848     targetDrop=0.0;
     3865    targetDrop=infPenalty+infPenalty2;
    38493866    if (iPass) {
    38503867      // get reduced costs
  • trunk/Clp/src/ClpSimplexPrimal.cpp

    r1034 r1057  
    22172217  // printf("Pert type %d perturbation %g, maxF %g\n",type,perturbation,maximumFraction);
    22182218  if (type==1) {
     2219    double tolerance = 100.0*primalTolerance_;
    22192220    //double multiplier = perturbation*maximumFraction;
    22202221    for (iSequence=0;iSequence<numberRows_+numberColumns_;iSequence++) {
    22212222      if (getStatus(iSequence)==basic) {
    2222         double solutionValue = solution_[iSequence];
    22232223        double lowerValue = lower_[iSequence];
    22242224        double upperValue = upper_[iSequence];
    2225         double difference = upperValue-lowerValue;
    2226         difference = CoinMin(difference,perturbation);
    2227         difference = CoinMin(difference,fabs(solutionValue)+1.0);
    2228         double value = maximumFraction*(difference+1.0);
    2229         value = CoinMin(value,0.1);
    2230         value *= CoinDrand48();
    2231         if (solutionValue-lowerValue<=primalTolerance_) {
    2232           lower_[iSequence] -= value;
    2233         } else if (upperValue-solutionValue<=primalTolerance_) {
    2234           upper_[iSequence] += value;
    2235         } else {
     2225        if (upperValue>lowerValue+tolerance) {
     2226          double solutionValue = solution_[iSequence];
     2227          double difference = upperValue-lowerValue;
     2228          difference = CoinMin(difference,perturbation);
     2229          difference = CoinMin(difference,fabs(solutionValue)+1.0);
     2230          double value = maximumFraction*(difference+1.0);
     2231          value = CoinMin(value,0.1);
     2232          value *= CoinDrand48();
     2233          if (solutionValue-lowerValue<=primalTolerance_) {
     2234            lower_[iSequence] -= value;
     2235          } else if (upperValue-solutionValue<=primalTolerance_) {
     2236            upper_[iSequence] += value;
     2237          } else {
    22362238#if 0
    2237           if (iSequence>=numberColumns_) {
    2238             // may not be at bound - but still perturb (unless free)
    2239             if (upperValue>1.0e30&&lowerValue<-1.0e30)
    2240               value=0.0;
    2241             else
    2242               value = - value; // as -1.0 in matrix
    2243           } else {
    2244             value = 0.0;
     2239            if (iSequence>=numberColumns_) {
     2240              // may not be at bound - but still perturb (unless free)
     2241              if (upperValue>1.0e30&&lowerValue<-1.0e30)
     2242                value=0.0;
     2243              else
     2244                value = - value; // as -1.0 in matrix
     2245            } else {
     2246              value = 0.0;
     2247            }
     2248#else
     2249            value=0.0;
     2250#endif
    22452251          }
    2246 #else
    2247           value=0.0;
    2248 #endif
    2249         }
    2250         if (value) {
    2251           if (printOut)
    2252             printf("col %d lower from %g to %g, upper from %g to %g\n",
    2253                    iSequence,lower_[iSequence],lowerValue,upper_[iSequence],upperValue);
    2254           if (solutionValue) {
    2255             largest = CoinMax(largest,value);
    2256             if (value>(fabs(solutionValue)+1.0)*largestPerCent)
    2257               largestPerCent=value/(fabs(solutionValue)+1.0);
    2258           } else {
    2259             largestZero = CoinMax(largestZero,value);
    2260           }
     2252          if (value) {
     2253            if (printOut)
     2254              printf("col %d lower from %g to %g, upper from %g to %g\n",
     2255                     iSequence,lower_[iSequence],lowerValue,upper_[iSequence],upperValue);
     2256            if (solutionValue) {
     2257              largest = CoinMax(largest,value);
     2258              if (value>(fabs(solutionValue)+1.0)*largestPerCent)
     2259                largestPerCent=value/(fabs(solutionValue)+1.0);
     2260            } else {
     2261              largestZero = CoinMax(largestZero,value);
     2262            }
     2263          }
    22612264        }
    22622265      }
Note: See TracChangeset for help on using the changeset viewer.