Ignore:
Timestamp:
Jun 26, 2007 5:59:58 AM (14 years ago)
Author:
forrest
Message:

update branches/devel for threads

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/devel/Cbc/src/Cbc_ampl.cpp

    r602 r642  
    3636#include "CoinMpsIO.hpp"
    3737#include "CoinFloatEqual.hpp"
    38 
     38#ifdef COIN_HAS_CLP
     39#include "OsiClpSolverInterface.hpp"
     40#endif
    3941#include "Cbc_ampl.h"
    4042extern "C" {
     
    5557    printf("string %s\n",v);
    5658  // Say algorithm found
    57   strcpy(algFound,kw->desc);;
     59  strcpy(algFound,kw->desc);
    5860  return v;
    5961}
     
    124126        { "wantsol",    WS_val,         NULL, "write .sol file (without -AMPL)" }
    125127        };
    126 static Option_Info Oinfo = {"cbc", "Cbc 1.01", "cbc_options", keywds, nkeywds, 0, "",
    127                                 0,decodePhrase,0,0,0, 20060130 };
     128static Option_Info Oinfo = {"cbc", "Cbc 1.04", "cbc_options", keywds, nkeywds, 0, "",
     129                                0,decodePhrase,0,0,0, 20070606 };
    128130// strdup used to avoid g++ compiler warning
    129131 static SufDecl
     
    344346  else
    345347    fileName[0]='\0';
     348  int nonLinearType=-1;
     349  // testosi parameter - if >= 10 then go in through coinModel
     350  for (i=1;i<argc;i++) {
     351    if (!strncmp(argv[i],"testosi",7)) {
     352      char * equals = strchr(argv[i],'=');
     353      if (equals&&atoi(equals+1)>=10&&atoi(equals+1)<=20) {
     354        nonLinearType=atoi(equals+1);
     355        break;
     356      }
     357    }
     358  }
    346359  int saveArgc = argc;
    347360  if (info->numberRows!=-1234567)
     
    373386  csd = suf_iput("sstatus", ASL_Sufkind_var, info->columnStatus);
    374387  rsd = suf_iput("sstatus", ASL_Sufkind_con, info->rowStatus);
    375   // testosi parameter - if >= 10 then go in through coinModel
    376   int testOsi=-1;
    377   for (i=0;i<saveArgc;i++) {
    378     if (!strncmp(saveArgv[i],"testosi",7)) {
    379       testOsi = atoi(saveArgv[i]+8);
    380       break;
    381     }
    382   }
    383   if (!(nlvc+nlvo)&&testOsi<10) {
     388  if (!(nlvc+nlvo)&&nonLinearType<10) {
    384389    /* read linear model*/
    385390    f_read(nl,0);
     
    424429    obj = (double *) malloc(n_var*sizeof(double));
    425430    for (i=0;i<n_var;i++)
    426       obj[i]=0.0;;
     431      obj[i]=0.0;
    427432    if (n_obj) {
    428433      for (og = Ograd[0];og;og = og->next)
     
    459464    info->numberRows=n_con;
    460465    info->numberColumns=n_var;
    461     info->numberElements=nzc;;
     466    info->numberElements=nzc;
    462467    info->numberBinary=nbv;
    463468    info->numberIntegers=niv;
     
    503508    if (!strstr(fileName,".nl"))
    504509      strcat(fileName,".nl");
    505     CoinModel * model = new CoinModel(1,fileName,info);
    506     if (model->numberRows()>0)
     510    CoinModel * model = new CoinModel((nonLinearType>10) ? 2 : 1,fileName,info);
     511    if (model->numberRows()>0||model->numberColumns()>0)
    507512      *coinModel=(void *) model;
    508513    Oinfo.uinfo = tempBuffer;
     
    517522    info->numberRows=n_con;
    518523    info->numberColumns=n_var;
    519     info->numberElements=nzc;;
     524    info->numberElements=nzc;
    520525    info->numberBinary=nbv;
    521     info->numberIntegers=niv;
    522     if (nbv+niv+nlvbi+nlvci+nlvoi>0) {
     526    info->numberIntegers=niv+nlvci+nlvoi;
     527    // No idea if there are overlaps so compute
     528    int numberIntegers=0;
     529    for ( i=0; i<n_var;i++) {
     530      if (model->columnIsInteger(i))
     531          numberIntegers++;
     532    }
     533    info->numberIntegers=numberIntegers;
     534    if (numberIntegers>0) {
    523535      mip_stuff(); // get any extra info
    524536      if (info->cut)
     
    761773  ASL_pfgh * asl_;
    762774  double * non_const_x_;
     775  int * column_; // for jacobian
     776  int * rowStart_;
     777  double * gradient_;
     778  double * constraintValues_;
    763779  int nz_h_full_; // number of nonzeros in hessian
    764780  int nerror_;
    765781  bool objval_called_with_current_x_;
    766782  bool conval_called_with_current_x_;
     783  bool jacval_called_with_current_x_;
    767784} CbcAmplInfo;
    768785
    769 static bool get_nlp_info(void * amplInfo,int & n, int & m, int & nnz_jac_g,
    770                               int & nnz_h_lag)
    771 {
    772   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    773   ASL_pfgh* asl = info->asl_;
    774  
    775   n = n_var; // # of variables
    776   m = n_con; // # of constraints
    777   nnz_jac_g = nzc; // # of non-zeros in the jacobian
    778   nnz_h_lag = info->nz_h_full_; // # of non-zeros in the hessian
    779  
    780   return true;
    781 }
    782 
    783 static bool get_bounds_info(void * amplInfo,int  n, double * x_l,
    784                             double * x_u, int  m, double * g_l, double * g_u)
    785 {
    786   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    787   ASL_pfgh* asl = info->asl_;
    788   assert(n == n_var);
    789   assert(m == n_con);
    790   int i;
    791   for (i=0; i<n; i++) {
    792     x_l[i] = LUv[2*i];
    793     x_u[i] = LUv[2*i+1];
    794   }
    795  
    796   for (i=0; i<m; i++) {
    797     g_l[i] = LUrhs[2*i];
    798     g_u[i] = LUrhs[2*i+1];
    799   }
    800   return true;
    801 }
    802 
    803 
    804 bool get_constraints_linearity(void * amplInfo,int  n,
    805       int * const_types)
    806 {
    807   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    808   ASL_pfgh* asl = info->asl_;
    809   //check that n is good
    810   assert(n == n_con);
    811   // check that there are no network constraints
    812   assert(nlnc == 0 && lnc == 0);
    813   //the first nlc constraints are non linear the rest is linear
    814   int i;
    815   for (i=0; i<nlc; i++) {
    816     const_types[i]=1;
    817   }
    818   // the rest is linear
    819   for (i=nlc; i<n_con; i++)
    820     const_types[i]=0;
    821   return true;
    822 }
    823 #if 0
    824 bool get_starting_point(int  n, bool init_x, double * x, bool init_z,
    825                         double * z_L, double * z_U, int  m, bool init_lambda, double * lambda)
    826 {
    827   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    828   ASL_pfgh* asl = info->asl_;
    829   assert(n == n_var);
    830   assert(m == n_con);
    831   int i;
    832  
    833   if (init_x) {
    834     for (i=0; i<n; i++) {
    835       if (havex0[i]) {
    836         x[i] = X0[i];
    837       }
    838       else {
    839         x[i] = 0.0;
    840       }
    841     }
    842   }
    843  
    844   if (init_z) {
    845     for (i=0; i<n; i++) {
    846       z_L[i] = z_U[i] = 1.0;
    847     }
    848   }
    849  
    850   if (init_lambda) {
    851     for (i=0; i<m; i++) {
    852       if (havepi0[i]) {
    853         lambda[i] = pi0[i];
    854       }
    855       else {
    856         lambda[i] = 0.0;
    857       }
    858     }
    859   }
    860  
    861   return true;
    862 }
    863 #endif
    864 static bool internal_objval(CbcAmplInfo * info ,double & obj_val)
    865 {
    866   ASL_pfgh* asl = info->asl_;
    867   info->objval_called_with_current_x_ = false; // in case the call below fails
    868 
    869   if (n_obj==0) {
    870     obj_val = 0;
    871     info->objval_called_with_current_x_ = true;
    872     return true;
    873   }  else {
    874     double  retval = objval(0, info->non_const_x_, (fint*)info->nerror_);
    875     if (!info->nerror_) {
    876       obj_val = info->obj_sign_*retval;
    877       info->objval_called_with_current_x_ = true;
    878       return true;
    879     } else {
    880       abort();
    881     }
    882   }
    883  
    884   return false;
    885 }
    886 
    887 static bool internal_conval(CbcAmplInfo * info ,int  m, double * g)
    888 {
    889   ASL_pfgh* asl = info->asl_;
    890   info->conval_called_with_current_x_ = false; // in case the call below fails
    891  
    892   bool allocated = false;
    893   if (!g) {
    894     g = new double[m];
    895     allocated = true;
    896   }
    897  
    898   conval(info->non_const_x_, g, (fint*)info->nerror_);
    899 
    900   if (allocated) {
    901     delete [] g;
    902     g = NULL;
    903   }
    904  
    905   if (!info->nerror_) {
    906     info->conval_called_with_current_x_ = true;
    907     return true;
    908   } else {
    909     abort();
    910   }
    911   return false;
    912 }
    913 
    914 
    915 static bool apply_new_x(CbcAmplInfo * info  ,bool new_x, int  n, const double * x)
    916 {
    917   ASL_pfgh* asl = info->asl_;
    918  
    919   if (new_x) {
    920     // update the flags so these methods are called
    921     // before evaluating the hessian
    922     info->conval_called_with_current_x_ = false;
    923     info->objval_called_with_current_x_ = false;
    924 
    925     //copy the data to the non_const_x_
    926     if (!info->non_const_x_) {
    927       info->non_const_x_ = new double [n];
    928     }
    929 
    930     for (int  i=0; i<n; i++) {
    931       info->non_const_x_[i] = x[i];
    932     }
    933    
    934     // tell ampl that we have a new x
    935     xknowne(info->non_const_x_, (fint*)info->nerror_);
    936     return info->nerror_ ? false : true;
    937   }
    938  
    939   return true;
    940 }
    941 bool eval_f(void * amplInfo,int  n, const double * x, bool new_x, double & obj_value)
    942 {
    943   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    944   if (!apply_new_x(info,new_x, n, x)) {
    945     return false;
    946   }
    947  
    948   return internal_objval(info,obj_value);
    949 }
    950 
    951 bool eval_grad_f(void * amplInfo,int  n, const double * x, bool new_x, double * grad_f)
    952 {
    953   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    954   ASL_pfgh* asl = info->asl_;
    955   if (!apply_new_x(info,new_x, n, x)) {
    956     return false;
    957   }
    958   int i;
    959  
    960   if (n_obj==0) {
    961     for (i=0; i<n; i++) {
    962       grad_f[i] = 0.;
    963     }
    964   }
    965   else {
    966     objgrd(0, info->non_const_x_, grad_f, (fint*)info->nerror_);
    967     if (info->nerror_) {
    968       return false;
    969     }
    970    
    971     if (info->obj_sign_==-1) {
    972       for (i=0; i<n; i++) {
    973         grad_f[i] = -grad_f[i];
    974       }
    975     }
    976   }
    977   return true;
    978 }
    979 
    980 bool eval_g(void * amplInfo,int  n, const double * x, bool new_x, int  m, double * g)
    981 {
    982   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    983   ASL_pfgh* asl = info->asl_;
    984   assert(n == n_var);
    985   assert(m == n_con);
    986  
    987   if (!apply_new_x(info,new_x, n, x)) {
    988     return false;
    989   }
    990  
    991   return internal_conval(info,m, g);
    992 }
    993 
    994 bool eval_jac_g(void * amplInfo,int  n, const double * x, bool new_x,
    995                 int  m, int  nele_jac, int * iRow,
    996                 int  *jCol, double * values)
    997 {
    998   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    999   ASL_pfgh* asl = info->asl_;
    1000   assert(n == n_var);
    1001   assert(m == n_con);
    1002   int i;
    1003  
    1004   if (iRow && jCol && !values) {
    1005     // setup the structure
    1006     int  current_nz = 0;
    1007     for (i=0; i<n_con; i++) {
    1008       for (cgrad* cg=Cgrad[i]; cg; cg = cg->next) {
    1009         //iRow[cg->goff] = i + 1;
    1010         iRow[cg->goff] = i ;
    1011         jCol[cg->goff] = cg->varno + 1;
    1012         //                              iRow[current_nz] = i + 1;
    1013         //                              jCol[current_nz] = cg->varno+1;
    1014         current_nz++;
    1015       }
    1016     }
    1017     assert(current_nz == nele_jac);
    1018     return true;
    1019   }
    1020   else if (!iRow && !jCol && values) {
    1021     if (!apply_new_x(info,new_x, n, x)) {
    1022       return false;
    1023     }
    1024    
    1025     jacval(info->non_const_x_, values, (fint*)info->nerror_);
    1026     if (!info->nerror_) {
    1027       return true;
    1028     } else {
    1029       abort();
    1030     }
    1031   }
    1032   else {
    1033     assert(false && "Invalid combination of iRow, jCol, and values pointers");
    1034   }
    1035  
    1036   return false;
    1037 }
    1038 
    1039 bool eval_h(void * amplInfo,int  n, const double * x, bool new_x,
    1040             double  obj_factor, int  m, const double * lambda,
    1041             bool new_lambda, int  nele_hess, int * iRow,
    1042             int * jCol, double * values)
    1043 {
    1044   CbcAmplInfo * info = (CbcAmplInfo *) amplInfo;
    1045   ASL_pfgh* asl = info->asl_;
    1046   assert(n == n_var);
    1047   assert(m == n_con);
    1048   int i;
    1049  
    1050   if (iRow && jCol && !values) {
    1051     // setup the structure
    1052     int k=0;
    1053     for (int i=0; i<n; i++) {
    1054       for (int j=sputinfo->hcolstarts[i]; j<sputinfo->hcolstarts[i+1]; j++) {
    1055         //iRow[k] = i + 1;
    1056         iRow[k] = i;
    1057         jCol[k] = sputinfo->hrownos[j]+1;
    1058         k++;
    1059       }
    1060     }
    1061     assert(k==nele_hess);
    1062     return true;
    1063   }
    1064   else if (!iRow & !jCol && values) {
    1065     if (!apply_new_x(info,new_x, n, x)) {
    1066       return false;
    1067     }
    1068     if (!info->objval_called_with_current_x_) {
    1069       double  dummy;
    1070       internal_objval(info,dummy);
    1071       internal_conval(info,m,NULL);
    1072     }
    1073     if (!info->conval_called_with_current_x_) {
    1074       internal_conval(info,m,NULL);
    1075     }
    1076     // copy lambda to non_const_lambda - note, we do not store a copy like
    1077     // we do with x since lambda is only used here and not in other calls
    1078     double * non_const_lambda = new double [m];
    1079     for (i=0; i<m; i++) {
    1080       non_const_lambda[i] = lambda[i];
    1081     }
    1082    
    1083     real ow=info->obj_sign_*obj_factor;
    1084     sphes(values, -1, &ow, non_const_lambda);
    1085    
    1086     delete [] non_const_lambda;
    1087     return true;
    1088   }
    1089   else {
    1090     assert(false && "Invalid combination of iRow, jCol, and values pointers");
    1091   }
    1092  
    1093   return false;
    1094 }
    1095786void
    1096787CoinModel::gdb( int nonLinear, const char * fileName,const void * info)
     
    1193884    objective = (double *) malloc(n_var*sizeof(double));
    1194885    for (i=0;i<n_var;i++)
    1195       objective[i]=0.0;;
     886      objective[i]=0.0;
    1196887    if (n_obj) {
    1197888      for (og = Ograd[0];og;og = og->next)
     
    1228919    numberRows=n_con;
    1229920    numberColumns=n_var;
    1230     numberElements=nzc;;
     921    numberElements=nzc;
    1231922    numberBinary=nbv;
    1232923    numberIntegers=niv;
     
    1282973    objective = (double *) malloc(n_var*sizeof(double));
    1283974    for (i=0;i<n_var;i++)
    1284       objective[i]=0.0;;
     975      objective[i]=0.0;
    1285976    if (n_obj) {
    1286977      for (og = Ograd[0];og;og = og->next)
     
    13301021    numberRows=n_con;
    13311022    numberColumns=n_var;
    1332     numberElements=nzc;;
     1023    numberElements=nzc;
    13331024    numberBinary=nbv;
    13341025    numberIntegers=niv;
     
    13531044    info->conval_called_with_current_x_=false;
    13541045    info->non_const_x_=NULL;
     1046    info->jacval_called_with_current_x_=false;
     1047    info->rowStart_ = NULL;
     1048    info->column_ = NULL;
     1049    info->gradient_ = NULL;
     1050    info->constraintValues_ = NULL;
    13551051  } else if (nonLinear==2) {
    13561052    // General nonlinear!
    1357     ASL_pfgh* asl = (ASL_pfgh*)ASL_alloc(ASL_read_pfgh);
    1358     // was asl = ASL_alloc(ASL_read_fg);
     1053    //ASL_pfgh* asl = (ASL_pfgh*)ASL_alloc(ASL_read_pfgh);
     1054    asl = ASL_alloc(ASL_read_pfgh);
    13591055    nl = jac0dim(stub, (long) strlen(stub));
    13601056    free(stub);
     
    14241120    // Fix memory leak one day
    14251121    CbcAmplInfo * info = new CbcAmplInfo;
     1122    moreInfo_ = (void *) info;
    14261123    //amplGamsData_ = info;
    1427     info->asl_=asl;
     1124    info->asl_=(ASL_pfgh *) asl;
    14281125    // This is not easy to get from ampl so save
    14291126    info->nz_h_full_= sphsetup(-1, coeff_obj, mult_supplied, uptri);
     
    14331130    info->conval_called_with_current_x_=false;
    14341131    info->non_const_x_=NULL;
     1132    info->jacval_called_with_current_x_=false;
     1133    // Look at nonlinear
     1134    if (nzc) {
     1135      n_conjac[1]=nlc; // just nonlinear
     1136      int * rowStart = new int [nlc+1];
     1137      info->rowStart_ = rowStart;
     1138      // See how many
     1139      int  current_nz = 0;
     1140      for (int i=0; i<nlc; i++) {
     1141        for (cgrad* cg=Cgrad[i]; cg; cg = cg->next) {
     1142          current_nz++;
     1143        }
     1144      }
     1145      // setup the structure
     1146      int * column = new int [current_nz];
     1147      info->column_ = column;
     1148      current_nz = 0;
     1149      rowStart[0]=0;
     1150      for (int i=0; i<nlc; i++) {
     1151        for (cgrad* cg=Cgrad[i]; cg; cg = cg->next) {
     1152          cg->goff=current_nz;
     1153          //iRow[cg->goff] = i ;
     1154          //jCol[cg->goff] = cg->varno + 1;
     1155          column[cg->goff] = cg->varno ;
     1156          current_nz++;
     1157        }
     1158        rowStart[i+1]=current_nz;
     1159      }
     1160      info->gradient_ = new double [nzc];
     1161      info->constraintValues_ = new double [nlc];
     1162    }
    14351163    /* objective*/
    14361164    objective = (double *) malloc(n_var*sizeof(double));
    14371165    for (i=0;i<n_var;i++)
    1438       objective[i]=0.0;;
     1166      objective[i]=0.0;
    14391167    if (n_obj) {
    14401168      for (og = Ograd[0];og;og = og->next)
     
    14871215    numberRows=n_con;
    14881216    numberColumns=n_var;
    1489     numberElements=nzc;;
     1217    numberElements=nzc;
    14901218    numberBinary=nbv;
    14911219    numberIntegers=niv;
     
    15221250  for ( i=numberColumns-numberBinary-numberIntegers;
    15231251        i<numberColumns;i++) {
    1524       setColumnIsInteger(i,true);;
     1252    setColumnIsInteger(i,true);
    15251253  }
    15261254  // and non linear
    15271255  for (i=numberAllNonLinearBoth-numberIntegerNonLinearBoth;
    15281256       i<numberAllNonLinearBoth;i++) {
    1529     setColumnIsInteger(i,true);;
     1257    setColumnIsInteger(i,true);
    15301258  }
    15311259  for (i=numberAllNonLinearConstraints-numberIntegerNonLinearConstraints;
    15321260       i<numberAllNonLinearConstraints;i++) {
    1533     setColumnIsInteger(i,true);;
     1261    setColumnIsInteger(i,true);
    15341262  }
    15351263  for (i=numberAllNonLinearObjective-numberIntegerNonLinearObjective;
    15361264       i<numberAllNonLinearObjective;i++) {
    1537     setColumnIsInteger(i,true);;
     1265    setColumnIsInteger(i,true);
    15381266  }
    15391267  free(columnLower);
Note: See TracChangeset for help on using the changeset viewer.