/home/coin/SVN-release/OS-2.4.0/OS/src/OSModelInterfaces/OSnl2osil.cpp

Go to the documentation of this file.
00001 /* $Id: OSnl2osil.cpp 4292 2011-09-21 05:47:18Z kmartin $ */
00019 #include <iostream>
00020 #include "OSiLWriter.h"
00021 #include "OSnl2osil.h"
00022 #include "OSErrorClass.h"
00023 
00024 #include "CoinHelperFunctions.hpp"
00025 
00026 
00027 #include "nlp.h"
00028 #include "getstub.h"
00029 #include "r_opn.hd" /* for N_OPS */
00030 #include "opcode.hd"
00031 
00032 
00033 #ifdef HAVE_CMATH
00034 # include <cmath>
00035 #else
00036 # ifdef HAVE_CMATH_H
00037 #  include <cmath.h>
00038 # endif
00039 #endif
00040 
00041 #include <sstream>
00042 
00043 
00044 #define R_OPS  ((ASL_fg*)asl)->I.r_ops_
00045 #define OBJ_DE ((ASL_fg*)asl)->I.obj_de_
00046 #define VAR_E  ((ASL_fg*)asl)->I.var_e_
00047 #define CON_DE ((ASL_fg*)asl)->I.con_de_
00048 
00049 efunc *r_ops_int[N_OPS];
00050 
00051 #include <asl.h>
00052 
00053 using std::cerr;
00054 using std::cout;
00055 using std::endl;
00056 
00057 #ifdef HAVE_STDINT_H
00058 #include <stdint.h>
00059 #endif
00060 
00061 
00062 
00063 #define AMPLDEBUG
00064 
00065 OSnl2osil::OSnl2osil(std::string nlfilename)
00066 {
00067     //Initialize the AMPL library
00068     asl = ASL_alloc( ASL_read_fg);
00069     stub = &nlfilename[ 0];
00070     //cout << "READING FILE " << stub << endl;
00071     //Initialize the nl file reading
00072     nl = jac0dim(stub, (fint)strlen(stub));
00073     //Prepare *columnwise* parsing of nl file
00074     A_vals = (real *)Malloc(nzc*sizeof(real));
00075     //read nl file as a linear program
00076 #ifdef AMPLDEBUG
00077     cout << "number of nonzeros =   " << nzc << endl;
00078     cout << "number of variable =   " << n_var << endl;
00079     cout << "number of constraints =   " << n_con << endl;
00080     cout << "number of objs =   " << n_obj << endl;
00081     cout << "number of ranges =   " << nranges << endl;
00082     cout << "number of equations =   " << n_eqn << endl;
00083     //cout << "objective function type" << *objtype << endl;
00084 #endif
00085 #ifdef AMPLDEBUG
00086     cout << "Start f_read()" << endl;
00087 #endif
00088     //X0 = (real *)Malloc( n_var*sizeof(real));
00089     if(N_OPS > 0)
00090     {
00091         for(int i = 0; i < N_OPS; i++)
00092         {
00093             r_ops_int[i] = (efunc*)(unsigned long)i;
00094 
00095         }
00096         R_OPS = r_ops_int;
00097         want_derivs = 0;
00098         fg_read(nl, 0);
00099         R_OPS = 0;
00100     }
00101     numkount = 0;
00102 }
00103 
00104 OSnl2osil::~OSnl2osil()
00105 {
00106     osinstance->instanceData->linearConstraintCoefficients->start->bDeleteArrays = false;
00107     osinstance->instanceData->linearConstraintCoefficients->rowIdx->bDeleteArrays = false;
00108     osinstance->instanceData->linearConstraintCoefficients->value->bDeleteArrays = false;
00109     delete osinstance;
00110     osinstance = NULL;
00111  //   free( X0);
00112     free( A_vals);
00113     ASL_free(&asl);
00114 
00115 }
00116 
00117 OSnLNode* OSnl2osil::walkTree (expr *e)
00118 {
00119     OSnLNode *nlNodePoint;
00120     OSnLNodeVariable *nlNodeVariablePoint;
00121     OSnLNodeNumber *nlNodeNumberPoint;
00122     efunc *op;
00123     expr **ep;
00124     int opnum;
00125     int i;
00126     //std::cout << "Variable Index " << varIdx << std::endl;
00127     int j = ((expr_v *)e - VAR_E) - osinstance->getVariableNumber() ;
00128     //std::cout << "GET OPERATOR NUMBER" << std::endl;
00129     op = e->op;
00130     opnum = Intcast op;
00131     //std::cout << "OPERATOR NUMBER = " << opnum << std::endl;
00132     //Printf ("op %d  optype %d  ", opnum, optype[opnum]);
00133     try
00134     {
00135         switch( opnum)
00136         {
00137         case OPPLUS:
00138             //cout << "FOUND  PLUS NODE"  << endl;
00139             nlNodePoint = new OSnLNodePlus();
00140             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00141             nlNodePoint->m_mChildren[1] = walkTree (e->R.e);
00142             op_type.push_back( "PLUS");
00143             return nlNodePoint;
00144 
00145         case OPSUMLIST:
00146             i = 0;
00147             //cout << "INSIDE SUM OPERATOR" << endl;
00148             nlNodePoint = new OSnLNodeSum();
00149             nlNodePoint->inumberOfChildren = e->R.ep - e->L.ep;
00150             nlNodePoint->m_mChildren = new OSnLNode*[ e->R.ep - e->L.ep];
00151             for (ep = e->L.ep; ep < e->R.ep; *ep++)
00152                 nlNodePoint->m_mChildren[i++] = walkTree ( *ep);
00153             return nlNodePoint;
00154 
00155         case MAXLIST:
00156             i = 0;
00157             //cout << "INSIDE MAX OPERATOR" << endl;
00158             nlNodePoint = new OSnLNodeMax();
00159             nlNodePoint->inumberOfChildren = e->R.ep - e->L.ep;
00160             nlNodePoint->m_mChildren = new OSnLNode*[ e->R.ep - e->L.ep];
00161             for (ep = e->L.ep; ep < e->R.ep; *ep++)
00162                 nlNodePoint->m_mChildren[i++] = walkTree ( *ep);
00163             return nlNodePoint;
00164 
00165         case OPMINUS:
00166             nlNodePoint = new OSnLNodeMinus();
00167 
00168             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00169             nlNodePoint->m_mChildren[1] = walkTree (e->R.e);
00170             op_type.push_back( "MINUS");
00171 
00172             return nlNodePoint;
00173 
00174         case OPUMINUS:
00175             nlNodePoint = new OSnLNodeNegate();
00176             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00177             return nlNodePoint;
00178 
00179         case OPMULT:
00180             //cout << "FOUND MULT NODE"  << endl;
00181             nlNodePoint = new OSnLNodeTimes();
00182             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00183             nlNodePoint->m_mChildren[1] = walkTree (e->R.e);
00184             op_type.push_back( "TIMES");
00185             return nlNodePoint;
00186 
00187         case OPDIV:
00188             nlNodePoint = new OSnLNodeDivide();
00189             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00190             nlNodePoint->m_mChildren[1] = walkTree (e->R.e);
00191             return nlNodePoint;
00192 
00193         case OPPOW:
00194             //cout << "FOUND OPPOW NODE"  << endl;
00195             nlNodePoint = new OSnLNodePower();
00196             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00197             nlNodePoint->m_mChildren[1] = walkTree (e->R.e);
00198             op_type.push_back( "POWER");
00199             return nlNodePoint;
00200 
00201 
00202         case OP1POW:
00203             //cout << "FOUND OP1POW NODE"  << endl;
00204             //cout << "OP1POW EXPONENT =  "  << e->R.en->v<<  endl;
00205             nlNodePoint = new OSnLNodePower();
00206 
00207 
00208             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00209             nlNodeNumberPoint = new OSnLNodeNumber();
00210             nlNodeNumberPoint->value = e->R.en->v;
00211             nlNodePoint->m_mChildren[1] = nlNodeNumberPoint;
00212 
00213             op_type.push_back( "NUMBER");
00214             op_type.push_back( os_dtoa_format( numkount) );
00215             operand.push_back( e->R.en->v );
00216             numkount++;
00217 
00218             op_type.push_back( "POWER");
00219 
00220 
00221             return nlNodePoint;
00222 
00223         case OP2POW:
00224             //cout << "FOUND OP2POW NODE"  << endl;
00225             //nlNodePoint = new OSnLNodePower();
00226             //nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00227             //nlNodeNumberPoint = new OSnLNodeNumber();
00228             //nlNodeNumberPoint->value = 2;
00229             //nlNodePoint->m_mChildren[1] = nlNodeNumberPoint;
00230             nlNodePoint = new OSnLNodeSquare();
00231 
00232             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00233             op_type.push_back( "SQUARE");
00234             return nlNodePoint;
00235 
00236         case OPCPOW:
00237             //cout << "FOUND OPCPOW NODE"  << endl;
00238             //cout << "OPCPOW EXPONENT =  "  << e->R.en->v<<  endl;
00239             nlNodePoint = new OSnLNodePower();
00240             nlNodeNumberPoint = new OSnLNodeNumber();
00241             nlNodeNumberPoint->value = e->L.en->v;
00242             nlNodePoint->m_mChildren[0] = nlNodeNumberPoint;
00243             nlNodePoint->m_mChildren[1] = walkTree (e->R.e);
00244             return nlNodePoint;
00245 
00246         case OP_log:
00247             nlNodePoint = new OSnLNodeLn();
00248             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00249             op_type.push_back( "LOG");
00250             return nlNodePoint;
00251 
00252         case OP_sqrt:
00253             nlNodePoint = new OSnLNodeSqrt();
00254             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00255             return nlNodePoint;
00256 
00257         case OP_cos:
00258             nlNodePoint = new OSnLNodeCos();
00259             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00260             return nlNodePoint;
00261 
00262         case OP_sin:
00263             nlNodePoint = new OSnLNodeSin();
00264             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00265             return nlNodePoint;
00266 
00267         case OP_exp:
00268             nlNodePoint = new OSnLNodeExp();
00269             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00270             return nlNodePoint;
00271 
00272         case ABS:
00273             nlNodePoint = new OSnLNodeAbs();
00274             nlNodePoint->m_mChildren[0] = walkTree (e->L.e);
00275             return nlNodePoint;
00276 
00277         case OPNUM:
00278             //cout << "found a number node" << endl;
00279             nlNodeNumberPoint = new OSnLNodeNumber;
00280             //cout << "THE NUMBER" << (double) ((expr_n*)e)->v << endl;
00281             nlNodeNumberPoint->value = (double) ((expr_n*)e)->v;
00282             op_type.push_back( "NUMBER");
00283             op_type.push_back( os_dtoa_format(  numkount ) );
00284 
00285             operand.push_back( (double) ((expr_n*)e)->v   );
00286             numkount++;
00287 
00288 
00289             return nlNodeNumberPoint;
00290 
00291         case OPVARVAL:
00292             //cout << "found a variable node" << endl;
00293             // treat the common expression or defined variables
00294             if( j >= 0 )
00295             {
00296                 // process common expression
00297                 /*
00298                 std::cout << "como = "  << como  << std::endl;
00299                 std::cout << "comc = "  << comc  << std::endl;
00300                 std::cout << "comb = "  << comb  << std::endl;
00301                 std::cout << "como1 = "  << como1  << std::endl;
00302                 std::cout << "comc1 = "  << comc1  << std::endl;
00303                 std::cout << "ncom0 = "  << ncom0  << std::endl;
00304                 std::cout << "jjjjjjjjjjjjjjjjjj = "  << j  << std::endl;
00305                 */
00306 
00307                 // Orban: http://www.gerad.ca/~orban/drampl/def-vars.html
00308                 if(j < ncom0)
00309                 {
00310                     struct cexp *common = ((const ASL_fg *) asl) -> I.cexps_ + j ;
00311                     //walk the tree for the non-linear stuff
00312 
00313                     // now add the linear terms
00314                     int nlin = common -> nlin;
00315                     //std::cout << "Number of linear terms in common expression " << nlin << std::endl;
00316 
00317                     if( nlin > 0)
00318                     {
00319                         nlNodePoint = new OSnLNodeSum();
00320                         nlNodePoint->inumberOfChildren = nlin + 1;
00321                         nlNodePoint->m_mChildren = new OSnLNode*[ nlin + 1];
00322                         // we have linear variables
00323                         // get the index and coefficient
00324                         linpart *L = common -> L;
00325                         for(int kj = 0; kj < nlin; kj++)
00326                         {
00327 
00328                             // get the coefficient
00329                             //std::cout << "Linear coefficient  "  << L [kj].fac << std::endl;
00330 
00331                             // get the index
00332                             //std::cout  << "Variable index  "  << ((uintptr_t) (L [kj].v.rp) - (uintptr_t) VAR_E) / sizeof (expr_v) << std::endl;
00333                             // add an OSnLSumNode with the linear terms
00334                             nlNodePoint->m_mChildren[ kj] = new OSnLNodeVariable;
00335                             nlNodeVariablePoint = (OSnLNodeVariable*)nlNodePoint->m_mChildren[ kj];
00336                             nlNodeVariablePoint->coef = L [kj]. fac;
00337                             nlNodeVariablePoint->idx = ((uintptr_t) (L [kj].v.rp) - (uintptr_t) VAR_E) / sizeof (expr_v);
00338                         }
00339                         nlNodePoint->m_mChildren[ nlin] = walkTree(  common->e);
00340                         return nlNodePoint;
00341                     }
00342                     else return walkTree(  common->e);
00343                 }
00344                 else
00345                 {
00346                     struct cexp1 *common = ((const ASL_fg *) asl) -> I.cexps1_ + (j - ncom0);
00347                     //walk the tree for the non-linear stuff
00348 
00349                     // now add the linear terms
00350                     int nlin = common -> nlin;
00351                     //std::cout << "Number of linear terms in common expression " << nlin << std::endl;
00352 
00353                     if( nlin > 0)
00354                     {
00355                         nlNodePoint = new OSnLNodeSum();
00356                         nlNodePoint->inumberOfChildren = nlin + 1;
00357                         nlNodePoint->m_mChildren = new OSnLNode*[ nlin + 1];
00358                         // we have linear variables
00359                         // get the index and coefficient
00360                         linpart *L = common -> L;
00361                         for(int kj = 0; kj < nlin; kj++)
00362                         {
00363 
00364                             // get the coefficient
00365                             //std::cout << "Linear coefficient  "  << L [kj].fac << std::endl;
00366 
00367                             // get the index
00368                             //std::cout  << "Variable index  "  << ((uintptr_t) (L [kj].v.rp) - (uintptr_t) VAR_E) / sizeof (expr_v) << std::endl;
00369                             // add an OSnLSumNode with the linear terms
00370                             nlNodePoint->m_mChildren[ kj] = new OSnLNodeVariable;
00371                             nlNodeVariablePoint = (OSnLNodeVariable*)nlNodePoint->m_mChildren[ kj];
00372                             nlNodeVariablePoint->coef = L [kj]. fac;
00373                             nlNodeVariablePoint->idx = ((uintptr_t) (L [kj].v.rp) - (uintptr_t) VAR_E) / sizeof (expr_v);
00374                         }
00375                         nlNodePoint->m_mChildren[ nlin] = walkTree(  common->e);
00376                         return nlNodePoint;
00377                     }
00378                     else return walkTree(  common->e);
00379                 }
00380             }
00381             //if(e->a > osinstance->getVariableNumber() ) throw ErrorClass("OS cannot handle AMPL user defined variables, please reformulate");
00382             nlNodeVariablePoint = new OSnLNodeVariable;
00383             nlNodeVariablePoint->idx = e->a;
00384             nlNodeVariablePoint->coef = 1.0;
00385 
00386             op_type.push_back( "VARIABLE");
00387             op_type.push_back( os_dtoa_format(  e->a  ) );
00388 
00389 
00390             return nlNodeVariablePoint;
00391             break;
00392         default:
00393             std::ostringstream outStr;
00394             std::string error;
00395             outStr  << endl;
00396             outStr  << endl;
00397             error = "ERROR:  An unsupported operator found, AMPL operator number =  "  ;
00398             outStr << error;
00399             outStr << opnum;
00400             outStr << endl;
00401             error = outStr.str();
00402             throw ErrorClass( error );
00403         }//end switch
00404     }//end try
00405     catch(const ErrorClass& eclass)
00406     {
00407         throw;
00408     }
00409 }//walkTree
00410 
00411 
00412 bool OSnl2osil::createOSInstance()
00413 {
00414     osinstance = new OSInstance();
00415     int i, j;
00416     // put in instanceHeader information
00417     //
00418     osinstance->setInstanceDescription("Generated from AMPL nl file");
00419     //
00420     // put in instanceData information
00421     //
00422     // get the variable information
00423     //
00424     std::string colName;
00425     char vartype = 'C';
00426     osinstance->setVariableNumber( n_var);
00427 
00428 
00429 
00430 
00431     //first the nonlinear variables
00432     //welcome to the world of the ASL API
00433 
00434     int lower;
00435     int upper;
00436     lower = 0;
00437     upper = nlvb - nlvbi;
00438 #ifdef AMPLDEBUG
00439     std::cout << "LOWER = " << lower << std::endl;
00440     std::cout << "UPPER = " << upper << std::endl;
00441 #endif
00442     for(i = lower; i < upper; i++)  //continuous in an objective and in a constraint
00443     {
00444         vartype = 'C';
00445         osinstance->addVariable(i, var_name(i),
00446                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00447                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00448                                 vartype);
00449     }
00450 
00451     lower = nlvb - nlvbi;
00452     upper = (nlvb - nlvbi) + nlvbi;
00453     upper = nlvb; //
00454 #ifdef AMPLDEBUG
00455     std::cout << "LOWER = " << lower << std::endl;
00456     std::cout << "UPPER = " << upper << std::endl;
00457 #endif
00458     for(i = lower; i < upper; i++)  //integer in an objective and in a constraint
00459     {
00460         vartype = 'I';
00461         osinstance->addVariable(i, var_name(i),
00462                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00463                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00464                                 vartype);
00465     }
00466 
00467     lower = nlvb;
00468     upper = nlvb + (nlvc - (nlvb + nlvci)) ;
00469     upper = nlvc - nlvci;
00470 #ifdef AMPLDEBUG
00471     std::cout << "LOWER = " << lower << std::endl;
00472     std::cout << "UPPER = " << upper << std::endl;
00473 #endif
00474     for(i = lower; i < upper; i++)  //continuous just in constraints
00475     {
00476         vartype = 'C';
00477         osinstance->addVariable(i, var_name(i),
00478                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00479                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00480                                 vartype);
00481     }
00482 
00483 
00484     lower = nlvc - nlvci;
00485     upper = nlvc - nlvci + nlvci;
00486     upper = nlvc;
00487 #ifdef AMPLDEBUG
00488     std::cout << "LOWER = " << lower << std::endl;
00489     std::cout << "UPPER = " << upper << std::endl;
00490 #endif
00491     for(i = lower; i < upper; i++)  //integer just in constraints
00492     {
00493         vartype = 'I';
00494         osinstance->addVariable(i, var_name(i),
00495                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00496                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00497                                 vartype);
00498     }
00499 
00500     lower = nlvc;
00501     upper = nlvc + ( nlvo - (nlvc + nlvoi) );
00502     upper = nlvo - nlvoi;
00503 #ifdef AMPLDEBUG
00504     std::cout << "LOWER = " << lower << std::endl;
00505     std::cout << "UPPER = " << upper << std::endl;
00506 #endif
00507     for(i = lower; i < upper; i++)  //continuous just in objectives
00508     {
00509         vartype = 'C';
00510         osinstance->addVariable(i, var_name(i),
00511                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00512                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00513                                 vartype);
00514     }
00515 
00516     lower = nlvo - nlvoi;
00517     upper = nlvo - nlvoi + nlvoi;
00518     upper = nlvo ;
00519 #ifdef AMPLDEBUG
00520     std::cout << "LOWER = " << lower << std::endl;
00521     std::cout << "UPPER = " << upper << std::endl;
00522 #endif
00523     for(i = lower; i < upper; i++)  //integer just in objectives
00524     {
00525         vartype = 'I';
00526         osinstance->addVariable(i, var_name(i),
00527                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00528                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00529                                 vartype);
00530     }
00531 
00532 
00533     //now the other variables
00534 
00535     lower = CoinMax(nlvc, nlvo);
00536     upper = CoinMax(nlvc, nlvo) + nwv;
00537 #ifdef AMPLDEBUG
00538     std::cout << "LOWER = " << lower << std::endl;
00539     std::cout << "UPPER = " << upper << std::endl;
00540 #endif
00541     for(i = lower; i < upper; i++)  //linear arc variables
00542     {
00543         vartype = 'C';
00544         osinstance->addVariable(i, var_name(i),
00545                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00546                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00547                                 vartype);
00548     }
00549 
00550 
00551     lower = CoinMax(nlvc, nlvo) + nwv;
00552     upper =  CoinMax(nlvc, nlvo) + nwv + (n_var - (CoinMax(nlvc, nlvo) + niv + nbv + nwv) );
00553     upper = n_var -  niv - nbv;
00554 #ifdef AMPLDEBUG
00555     std::cout << "LOWER = " << lower << std::endl;
00556     std::cout << "UPPER = " << upper << std::endl;
00557 #endif
00558     for(i = lower; i < upper; i++)  //other linear
00559     {
00560         vartype = 'C';
00561         osinstance->addVariable(i, var_name(i),
00562                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00563                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00564                                 vartype);
00565     }
00566 
00567 
00568     lower = n_var -  niv - nbv;
00569     upper = n_var -  niv - nbv + nbv;
00570     upper = n_var -  niv ;
00571 #ifdef AMPLDEBUG
00572     std::cout << "LOWER = " << lower << std::endl;
00573     std::cout << "UPPER = " << upper << std::endl;
00574 #endif
00575     for(i = lower; i < upper; i++)  //linear binary
00576     {
00577         vartype = 'B';
00578         osinstance->addVariable(i, var_name(i),
00579                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00580                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00581                                 vartype);
00582     }
00583 
00584 
00585 
00586     lower = n_var -  niv;
00587     upper = n_var -  niv  + niv;
00588     upper =   n_var;
00589 #ifdef AMPLDEBUG
00590     std::cout << "LOWER = " << lower << std::endl;
00591     std::cout << "UPPER = " << upper << std::endl;
00592 #endif
00593     for(i = lower; i < upper; i++)  //linear integer
00594     {
00595         vartype = 'I';
00596         osinstance->addVariable(i, var_name(i),
00597                                 LUv[2*i]   > -OSDBL_MAX ? LUv[2*i]   : -OSDBL_MAX,
00598                                 LUv[2*i+1] <  OSDBL_MAX ? LUv[2*i+1] :  OSDBL_MAX,
00599                                 vartype);
00600     }
00601 
00602 
00603     // end of variables -- thank goodness!!!
00604 
00605 
00606 
00607 
00608     //
00609     //
00610     //(expr_v *)e;
00611     //
00612     //
00613     // now create the objective function
00614     // in the nl file, this is stored in dense form; convert to sparse.
00615     //
00616     double objWeight = 1.0;
00617     //  char    *objtype;       /* object type array: 0 == min, 1 == max */
00618     std::string objName="";
00619     SparseVector* objectiveCoefficients = NULL;
00620 
00621     osinstance->setObjectiveNumber( n_obj) ;
00622     for(i = 0; i < n_obj; i++)
00623     {
00624         int n_obj_coef = 0;
00625         for(og = Ograd[i]; og; og = og->next)
00626         {
00627             if (og->coef != 0) n_obj_coef++;
00628         }
00629         objectiveCoefficients = new SparseVector( n_obj_coef);
00630         int i_obj_coef = -1;
00631         for(og = Ograd[i]; og; og = og->next)
00632         {
00633             if (fabs(og->coef) > OS_EPS)
00634             {
00635                 i_obj_coef++;
00636                 objectiveCoefficients->values[i_obj_coef] = og->coef;
00637                 objectiveCoefficients->indexes[i_obj_coef] = og->varno;
00638             }
00639         }
00640         osinstance->addObjective(-n_obj + i, objName,
00641                                  (objtype[i] == 1)?"max":"min",
00642                                  objconst( i),  objWeight, objectiveCoefficients) ;
00643         delete objectiveCoefficients; // delete the temporary sparse vector
00644         objectiveCoefficients = NULL;
00645     }
00646     //
00647     // now fill in row information
00648     //
00649     osinstance->setConstraintNumber( n_con);
00650     // kipp -- important  -- figure out where the nl file stores a rhs constant
00651     double constant = 0.0;
00652     std::string rowName;
00653     for(i = 0; i < n_con; i++)
00654     {
00655         osinstance->addConstraint(i, con_name(i),
00656                                   LUrhs[2*i] > -OSDBL_MAX ? LUrhs[2*i] : -OSDBL_MAX,
00657                                   LUrhs[2*i+1] < OSDBL_MAX ? LUrhs[2*i+1] : OSDBL_MAX,
00658                                   constant);
00659     }
00660     //
00661     // Now the A-matrix
00662     //
00663     int colStart, colEnd, nCoefSqueezed;
00664     nCoefSqueezed = 0;
00665 
00666 #ifdef AMPLDEBUG
00667     cout << "A-matrix elements: ";
00668     for (i = 0; i < A_colstarts[ n_var]; i++)
00669         cout << A_vals[i] << " ";
00670     cout << endl;
00671     cout << "A-matrix rowinfo: ";
00672     for (i = 0; i < A_colstarts[ n_var]; i++)
00673         cout << A_rownos[i] << " ";
00674     cout << endl;
00675     cout << "A-matrix colstart: ";
00676     for (i = 0; i <= n_var; i++)
00677         cout << A_colstarts[i] << " ";
00678     cout << endl;
00679 #endif
00680 
00681     colEnd = 0;
00682     for (i = 0; i < n_var; i++)
00683     {
00684         colStart = colEnd;
00685         colEnd   = A_colstarts[i+1];
00686 #ifdef AMPLDEBUG
00687         cout << "col " << i << " from " << colStart << " to " << colEnd - 1 << endl;
00688 #endif
00689         for (j = colStart; j < colEnd; j++)
00690         {
00691             if (fabs(A_vals[ j]) > OS_EPS)
00692             {
00693                 A_vals[ j-nCoefSqueezed] = A_vals[ j];
00694                 A_rownos[ j-nCoefSqueezed] = A_rownos[j];
00695             }
00696             else
00697             {
00698 #ifdef AMPLDEBUG
00699                 cout << "squeeze out element " << j << endl;
00700 #endif
00701                 nCoefSqueezed++;
00702             }
00703         }
00704         A_colstarts[i+1] = A_colstarts[i+1] - nCoefSqueezed;
00705     }
00706 
00707 #ifdef AMPLDEBUG
00708     cout << "A-matrix elements: ";
00709     for (i = 0; i < A_colstarts[ n_var]; i++)
00710         cout << A_vals[i] << " ";
00711     cout << endl;
00712     cout << "A-matrix rowinfo: ";
00713     for (i = 0; i < A_colstarts[ n_var]; i++)
00714         cout << A_rownos[i] << " ";
00715     cout << endl;
00716     cout << "A-matrix colstart: ";
00717     for (i = 0; i <= n_var; i++)
00718         cout << A_colstarts[i] << " ";
00719     cout << endl;
00720     cout << "A-matrix nonzeroes: " << A_colstarts[ n_var] << "; nsqueezed: " << nCoefSqueezed << endl;
00721 #endif
00722 
00723     if(A_colstarts[ n_var] > 0)
00724     {
00725         osinstance->setLinearConstraintCoefficients(A_colstarts[ n_var],  true,
00726                 A_vals,   0,  A_colstarts[n_var] - 1,
00727                 A_rownos, 0,  A_colstarts[n_var] - 1,
00728                 A_colstarts,  0,  n_var);
00729     }
00730     /*  int valuesBegin = 0;
00731         int valuesEnd = A_colstarts[ n_var] - 1;
00732         int startsBegin = 0;
00733         int indexesBegin = 0;
00734         int indexesEnd = A_colstarts[n_var] - 1;
00735 
00736         // if A_vals has only zeros don't generate a linearConstraints section
00737         bool bNumAvalsIsPositive = false;
00738         i = valuesBegin;
00739         while( (i < valuesEnd) && (bNumAvalsIsPositive == false) ){
00740                 if(A_vals[ i] != 0) bNumAvalsIsPositive = true;
00741                 i++;
00742         }
00743         if(bNumAvalsIsPositive == true){
00744                 osinstance->setLinearConstraintCoefficients(nzc,  true,
00745                         A_vals, valuesBegin,  valuesEnd,
00746                         A_rownos,  indexesBegin,  indexesEnd,
00747                         A_colstarts,  startsBegin,  n_var);
00748         }
00749     */
00750 
00751     // Kipp: can AMPL identify QPs???
00752     //osinstance->setQuadraticTerms(numberOfQPTerms, VarOneIdx, VarTwoIdx, Coeff, begin, end);
00753     //loop over each row with a nonlinear term
00754     //
00755     if((nlc + nlo) > 0)
00756     {
00757         OSnLNode* m_treeRoot;
00758         //cout << nlc << endl;
00759         //cout << nlo << endl;
00760         osinstance->instanceData->nonlinearExpressions->numberOfNonlinearExpressions = nlc + nlo;
00761         osinstance->instanceData->nonlinearExpressions->nl = new Nl*[ nlc + nlo ];
00762         int iNLidx = 0;
00763         std::cout << "WALK THE TREE FOR NONLINEAR CONSTRAINT TERMS:  " << nlc << std::endl;
00764 
00765 
00766         if(nlc > 0)
00767         {
00768             while (iNLidx < nlc)
00769             {
00770                 m_treeRoot = walkTree ((CON_DE + iNLidx)->e);
00771                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx] = new Nl();
00772                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx]->idx = iNLidx;
00773                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx]->osExpressionTree = new OSExpressionTree();
00774                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx]->osExpressionTree->m_treeRoot = m_treeRoot;
00775                 iNLidx++;
00776                 //std::cout << m_treeRoot->getNonlinearExpressionInXML() << std::endl;
00777             }
00778         }
00779 
00780         std::cout << "WALK THE TREE FOR NONLINEAR OBJECTIVE TERMS  " << nlo  << std::endl;
00781         if(nlo > 0)
00782         {
00783             while ( iNLidx < nlc + nlo)
00784             {
00785                 m_treeRoot = walkTree ((OBJ_DE + iNLidx - nlc)->e);
00786                 //std::cout << "CREATING A NEW NONLINEAR TERM IN THE OBJECTIVE" << std::endl;
00787                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx] = new Nl();
00788                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx]->idx = -1 - (iNLidx - nlc);
00789                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx]->osExpressionTree = new OSExpressionTree();
00790                 osinstance->instanceData->nonlinearExpressions->nl[ iNLidx]->osExpressionTree->m_treeRoot = m_treeRoot;
00791                 iNLidx++;
00792                 //std::cout << m_treeRoot->getNonlinearExpressionInXML() << std::endl;
00793             }
00794         }
00795         //std::cout << "DONE WALKING THE TREE FOR NONLINEAR OBJECTIVE TERMS" << std::endl;
00796 
00797 
00798     }
00799 //      delete objectiveCoefficients;
00800 //      objectiveCoefficients = NULL;
00801     //
00802     // end loop of nonlinear rows
00803     //
00804 #ifdef AMPLDEBUG
00805     OSiLWriter osilwriter;
00806     std::cout << "WRITE THE INSTANCE" << std::endl;
00807     osilwriter.m_bWhiteSpace = true;
00808     std::cout << osilwriter.writeOSiL( osinstance) << std::endl;
00809     std::cout << "DONE WRITE THE INSTANCE" << std::endl;
00810 
00811     std::cout << osinstance->printModel() << std::endl;
00812 
00813 #endif
00814 
00815     return true;
00816 }

Generated on Thu Sep 22 03:06:03 2011 by  doxygen 1.4.7