- Main Page
- Namespaces
- Classes
- Files
- Related Pages

- File List
- File Members

# /home/coin/SVN-release/OS-2.3.5/OS/src/OSModelInterfaces/OSgams2osil.cpp

Go to the documentation of this file.

    00001 // Copyright (C) GAMS Development and others 2007-2009
    00002 // All Rights Reserved.
    00003 // This code is published under the Common Public License.
    00004 //
    00005 // $Id: OSgams2osil.cpp 2898 2009-09-01 05:21:45Z kmartin $
    00006 //
    00007 // Author: Stefan Vigerske
    00008
    00009 #include "OSgams2osil.hpp"
    00010 #include "GAMSlinksConfig.h"
    00011 #include "GamsNLinstr.h"
    00012
    00013 #include "OSInstance.h"
    00014 #include "CoinHelperFunctions.hpp"
    00015
    00016 #include "gmomcc.h"
    00017 #include "gevmcc.h"
    00018
    00019 #include <sstream>
    00020
    00021 OSgams2osil::OSgams2osil(gmoHandle_t gmo_)
    00022 : gev(gmo_ ? (gevHandle_t)gmoEnvironment(gmo_) : NULL), gmo(gmo_), osinstance(NULL)
    00023 {
    00024         //assert(gev_ == NULL); // TODO better
    00025         assert(gmo_ == NULL); // TODO better
    00026 }
    00027
    00028 OSgams2osil::OSgams2osil(std::string gamsControlFile)
    00029 : gev(NULL), gmo(NULL), osinstance(NULL)
    00030 {
    00031         initGMO( gamsControlFile.c_str() ) ;
    00032 }
    00033
    00034 OSgams2osil::~OSgams2osil() {
    00035         delete osinstance;
    00036
    00037         //close output channels
    00038 //      gmoCloseGms(gmo);
    00039         gmoFree(&gmo);
    00040         gevFree(&gev);
    00041         gmoLibraryUnload();
    00042         gevLibraryUnload();
    00043 }
    00044
    00045 bool OSgams2osil::initGMO(const char* datfile) {
    00046         assert(gmo == NULL);
    00047         assert(gev == NULL);
    00048
    00049         char msg[1024];
    00050         int rc;
    00051
    00052         //if (!gmoCreate(&gmo, msg, sizeof(msg))) {
    00053   //            fprintf(stderr, "%s\n",msg);
    00054   //            return false;
    00055   //    }
    00056
    00057
    00058 //      if (!gmoCreateD(&gmo, GAMSIO_PATH, msg, sizeof(msg))) {
    00059                 if (!gmoCreate(&gmo, msg, sizeof(msg))) {
    00060                         fprintf(stderr, "%s\n",msg);
    00061                         return false;
    00062                 }
    00063 //      }
    00064
    00065         if (!gevCreate(&gev, msg, sizeof(msg))) {
    00066                 fprintf(stderr, "%s\n",msg);
    00067                 return EXIT_FAILURE;
    00068         }
    00069
    00070
    00071         gmoIdentSet(gmo, "OS link object");
    00072
    00073         // load control file
    00074   if ((rc = gevInitEnvironmentLegacy(gev, datfile))) {
    00075         fprintf(stderr, "Could not init gams environment: %s Rc = %d\n", datfile, rc);
    00076     gmoFree(&gmo);
    00077     gevFree(&gev);
    00078         return false;
    00079   }
    00080
    00081   if ((rc = gmoRegisterEnvironment(gmo, gev, msg))) {
    00082         gevLogStat(gev, "Could not register environment.");
    00083         gmoFree(&gmo);
    00084         gevFree(&gev);
    00085         return false;
    00086   }
    00087
    00088         // setup GAMS output channels
    00089 //      if ((rc = gmoOpenGms( gmo))) {
    00090 //              fprintf(stderr, "Could not open GAMS environment. Rc = %d\n", rc);
    00091 //              gmoFree(&gmo);
    00092 //              return false;
    00093 //      }
    00094
    00095         if ((rc = gmoLoadDataLegacy(gmo, msg))) {
    00096                 gevLogStat(gev, "Could not load model data.");
    00097 //              gmoCloseGms(gmo);
    00098                 gmoFree(&gmo);
    00099                 gevFree(&gev);
    00100                 return false;
    00101         }
    00102
    00103         gmoMinfSet(gmo, -OSDBL_MAX);
    00104         gmoPinfSet(gmo,  OSDBL_MAX);
    00105         gmoObjReformSet(gmo, 1);
    00106         gmoObjStyleSet(gmo, ObjType_Fun);
    00107         gmoIndexBaseSet(gmo, 0);
    00108
    00109         return true;
    00110 }
    00111
    00112 bool OSgams2osil::createOSInstance() {
    00113         assert(gmo != NULL);
    00114         assert(gev != NULL);
    00115
    00116         osinstance = new OSInstance();
    00117         int i, j;
    00118         char buffer[255];
    00119
    00120         // unfortunately, we do not know the model name
    00121         osinstance->setInstanceDescription("Generated from GAMS modeling object");
    00122         osinstance->setVariableNumber(gmoN(gmo));
    00123
    00124         char*        vartypes = new char[gmoN(gmo)];
    00125         std::string* varnames = new std::string[gmoN(gmo)];
    00126         for(i = 0; i < gmoN(gmo); ++i) {
    00127                 switch (gmoGetVarTypeOne(gmo, i)) {
    00128                         case var_X:
    00129                                 vartypes[i] = 'C';
    00130                                 break;
    00131                         case var_B:
    00132                                 vartypes[i] = 'B';
    00133                                 break;
    00134                         case var_I:
    00135                                 vartypes[i] = 'I';
    00136                                 break;
    00137                         default : {
    00138                                 // TODO: how to represent semicontinuous var. and SOS in OSiL ?
    00139                                 gevLogStat(gev, "Error: Unsupported variable type.");
    00140                                 return false;
    00141                         }
    00142                 }
    00143                 gmoGetVarNameOne(gmo, i, buffer);
    00144                 varnames[i] = buffer;
    00145         }
    00146
    00147         double* varlow = new double[gmoN(gmo)];
    00148         double* varup  = new double[gmoN(gmo)];
    00149         gmoGetVarLower(gmo, varlow);
    00150         gmoGetVarUpper(gmo, varup);
    00151
    00152         if (!osinstance->setVariables(gmoN(gmo), varnames, varlow, varup, vartypes))
    00153                 return false;
    00154
    00155         delete[] vartypes;
    00156         delete[] varnames;
    00157         delete[] varlow;
    00158         delete[] varup;
    00159
    00160         if (gmoModelType(gmo) == Proc_cns) { // no objective in constraint satisfaction models
    00161                 osinstance->setObjectiveNumber(0);
    00162         } else { // setup objective
    00163                 osinstance->setObjectiveNumber(1);
    00164
    00165                 SparseVector* objectiveCoefficients = new SparseVector(gmoObjNZ(gmo) - gmoObjNLNZ(gmo));
    00166
    00167                 int* colidx = new int[gmoObjNZ(gmo)];
    00168                 double* val = new double[gmoObjNZ(gmo)];
    00169                 int* nlflag = new int[gmoObjNZ(gmo)];
    00170                 int* dummy  = new int[gmoObjNZ(gmo)];
    00171
    00172                 if (gmoObjNZ(gmo)) nlflag[0] = 0; // workaround for gmo bug
    00173                 gmoGetObjSparse(gmo, colidx, val, nlflag, dummy, dummy);
    00174                 for (i = 0, j = 0; i < gmoObjNZ(gmo); ++i) {
    00175                         if (nlflag[i]) continue;
    00176                         objectiveCoefficients->indexes[j] = colidx[i];
    00177                         objectiveCoefficients->values[j]  = val[i];
    00178                         j++;
    00179                   assert(j <= gmoObjNZ(gmo) - gmoObjNLNZ(gmo));
    00180                 }
    00181           assert(j == gmoObjNZ(gmo) - gmoObjNLNZ(gmo));
    00182
    00183                 delete[] colidx;
    00184                 delete[] val;
    00185                 delete[] nlflag;
    00186                 delete[] dummy;
    00187
    00188                 std::string objname = "objective"; //TODO
    00189 //              if (dict.haveNames() && dict.getObjName(buffer, 256))
    00190 //                      objname = buffer;
    00191 //              std::cout << "gmo obj con: " << gmoObjConst(gmo) << std::endl;
    00192                 if (!osinstance->addObjective(-1, objname, gmoSense(gmo) == Obj_Min ? "min" : "max", gmoObjConst(gmo), 1., objectiveCoefficients)) {
    00193                         delete objectiveCoefficients;
    00194                         return false;
    00195                 }
    00196                 delete objectiveCoefficients;
    00197         }
    00198
    00199         osinstance->setConstraintNumber(gmoM(gmo));
    00200
    00201         double lb, ub;
    00202         for (i = 0;  i < gmoM(gmo);  ++i) {
    00203                 switch (gmoGetEquTypeOne(gmo, i)) {
    00204                         case equ_E:
    00205                                 lb = ub = gmoGetRhsOne(gmo, i);
    00206                                 break;
    00207                         case equ_L:
    00208                                 lb = -OSDBL_MAX;
    00209                                 ub =  gmoGetRhsOne(gmo, i);
    00210                                 break;
    00211                         case equ_G:
    00212                                 lb = gmoGetRhsOne(gmo, i);
    00213                                 ub = OSDBL_MAX;
    00214                                 break;
    00215                         case equ_N:
    00216                                 lb = -OSDBL_MAX;
    00217                                 ub =  OSDBL_MAX;
    00218                                 break;
    00219                         default:
    00220                                 gevLogStat(gev, "Error: Unknown row type. Exiting ...");
    00221                                 return false;
    00222                 }
    00223                 std::string conname;
    00224                 gmoGetEquNameOne(gmo, i, buffer);
    00225                 conname = buffer;
    00226                 if (!osinstance->addConstraint(i, conname, lb, ub, 0.))
    00227                         return false;
    00228         }
    00229
    00230         int nz = gmoNZ(gmo);
    00231         double* values  = new double[nz];
    00232         int* colstarts  = new int[gmoN(gmo)+1];
    00233         int* rowindexes = new int[nz];
    00234         int* nlflags    = new int[nz];
    00235
    00236         gmoGetMatrixCol(gmo, colstarts, rowindexes, values, nlflags);
    00237 //      for (i = 0; i < gmoNZ(gmo); ++i)
    00238 //              if (nlflags[i]) values[i] = 0.;
    00239         colstarts[gmoN(gmo)] = nz;
    00240
    00241         int shift = 0;
    00242         for (int col = 0; col < gmoN(gmo); ++col) {
    00243                 colstarts[col+1] -= shift;
    00244                 int k = colstarts[col];
    00245                 while (k < colstarts[col+1]) {
    00246                         values[k] = values[k+shift];
    00247                         rowindexes[k] = rowindexes[k+shift];
    00248                         if (nlflags[k+shift]) {
    00249                                 ++shift;
    00250                                 --colstarts[col+1];
    00251                         } else {
    00252                                 ++k;
    00253                         }
    00254                 }
    00255         }
    00256         nz -= shift;
    00257
    00258         if (!osinstance->setLinearConstraintCoefficients(nz, true,
    00259                 values, 0, nz-1,
    00260                 rowindexes, 0, nz-1,
    00261                 colstarts, 0, gmoN(gmo))) {
    00262                 delete[] nlflags;
    00263                 return false;
    00264         }
    00265
    00266         // values, colstarts, rowindexes are deleted by OSInstance
    00267         delete[] nlflags;
    00268
    00269         if (!gmoObjNLNZ(gmo) && !gmoNLNZ(gmo)) // everything linear -> finished
    00270                 return true;
    00271
    00272         osinstance->instanceData->nonlinearExpressions->numberOfNonlinearExpressions = gmoNLM(gmo) + (gmoObjNLNZ(gmo) ? 1 : 0);
    00273         osinstance->instanceData->nonlinearExpressions->nl = CoinCopyOfArrayOrZero((Nl**)NULL, osinstance->instanceData->nonlinearExpressions->numberOfNonlinearExpressions);
    00274         int iNLidx = 0;
    00275
    00276         int* opcodes = new int[gmoMaxSingleFNL(gmo)+1];
    00277         int* fields  = new int[gmoMaxSingleFNL(gmo)+1];
    00278         int constantlen = gmoNLConst(gmo);
    00279         double* constants = (double*)gmoPPool(gmo);
    00280         int codelen;
    00281
    00282         OSnLNode* nl;
    00283         if (gmoObjNLNZ(gmo)) {
    00284                 std::clog << "parsing nonlinear objective instructions" << std::endl;
    00285                 gmoDirtyGetObjFNLInstr(gmo, &codelen, opcodes, fields);
    00286
    00287                 nl = parseGamsInstructions(codelen, opcodes, fields, constantlen, constants);
    00288                 if (!nl) return false;
    00289
    00290                 double objjacval = gmoObjJacVal(gmo);
    00291                 std::clog << "obj jac val: " << objjacval << std::endl;
    00292                 if (objjacval == 1.) { // scale by -1/objjacval = negate
    00293                         OSnLNode* negnode = new OSnLNodeNegate;
    00294                         negnode->m_mChildren[0] = nl;
    00295                         nl = negnode;
    00296                 } else if (objjacval != -1.) { // scale by -1/objjacval
    00297                         OSnLNodeNumber* numbernode = new OSnLNodeNumber();
    00298                         numbernode->value = -1/objjacval;
    00299                         OSnLNodeTimes* timesnode = new OSnLNodeTimes();
    00300                         timesnode->m_mChildren[0] = nl;
    00301                         timesnode->m_mChildren[1] = numbernode;
    00302                         nl = timesnode;
    00303                 }
    00304                 assert(iNLidx < osinstance->instanceData->nonlinearExpressions->numberOfNonlinearExpressions);
    00305                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx] = new Nl();
    00306                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx]->idx = -1;
    00307                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx]->osExpressionTree = new OSExpressionTree();
    00308                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx]->osExpressionTree->m_treeRoot = nl;
    00309                 ++iNLidx;
    00310         }
    00311
    00312         for (i = 0; i < gmoM(gmo); ++i) {
    00313                 if (gmoDirtyGetRowFNLInstr(gmo, i, &codelen, opcodes, fields)) {
    00314                         std::clog << "got nonzero return at constraint " << i << std::endl;
    00315                 }
    00316                 if (!codelen) continue;
    00317                 std::clog << "parsing " << codelen << " nonlinear instructions of constraint " << osinstance->getConstraintNames()[i] << std::endl;
    00318                 nl = parseGamsInstructions(codelen, opcodes, fields, constantlen, constants);
    00319                 if (!nl) return false;
    00320                 assert(iNLidx < osinstance->instanceData->nonlinearExpressions->numberOfNonlinearExpressions);
    00321                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx] = new Nl();
    00322                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx]->idx = i; // correct that this is the con. number?
    00323                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx]->osExpressionTree = new OSExpressionTree();
    00324                 osinstance->instanceData->nonlinearExpressions->nl[iNLidx]->osExpressionTree->m_treeRoot = nl;
    00325                 ++iNLidx;
    00326         }
    00327         assert(iNLidx == osinstance->instanceData->nonlinearExpressions->numberOfNonlinearExpressions);
    00328
    00329         return true;
    00330 }
    00331
    00332 OSInstance* OSgams2osil::takeOverOSInstance() {
    00333         OSInstance* osinst = osinstance;
    00334         osinstance = NULL;
    00335         return osinst;
    00336 }
    00337
    00338 OSnLNode* OSgams2osil::parseGamsInstructions(int codelen, int* opcodes, int* fields, int constantlen, double* constants) {
    00339         std::vector<OSnLNode*> nlNodeVec;
    00340
    00341         const bool debugoutput = false;
    00342
    00343 //      for (int i=0; i<codelen; ++i)
    00344 //              std::clog << i << '\t' << GamsOpCodeName[opcodes[i+1]] << '\t' << fields[i+1] << std::endl;
    00345
    00346         nlNodeVec.reserve(codelen);
    00347
    00348         for (int i=0; i<codelen; ++i) {
    00349                 GamsOpCode opcode = (GamsOpCode)opcodes[i];
    00350                 int address = fields[i]-1;
    00351
    00352                 if (debugoutput) std::clog << '\t' << GamsOpCodeName[opcode] << ": ";
    00353 //              if (opcode == nlStore) {
    00354 //                      std::clog << "stop" << std::endl;
    00355 //                      break;
    00356 //              }
    00357                 switch(opcode) {
    00358                         case nlNoOp : { // no operation
    00359                                 if (debugoutput) std::clog << "ignored" << std::endl;
    00360                         } break;
    00361                         case nlPushV : { // push variable
    00362                                 address = gmoGetjSolver(gmo, address);
    00363                                 if (debugoutput) std::clog << "push variable " << osinstance->getVariableNames()[address] << std::endl;
    00364                                 OSnLNodeVariable *nlNode = new OSnLNodeVariable();
    00365                                 nlNode->idx=address;
    00366                                 nlNodeVec.push_back( nlNode );
    00367                         } break;
    00368                         case nlPushI : { // push constant
    00369                                 if (debugoutput) std::clog << "push constant " << constants[address] << std::endl;
    00370                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00371                                 nlNode->value = constants[address];
    00372                                 nlNodeVec.push_back( nlNode );
    00373                         } break;
    00374                         case nlStore: { // store row
    00375                                 if (debugoutput) std::clog << "ignored" << std::endl;
    00376                         } break;
    00377                         case nlAdd : { // add
    00378                                 if (debugoutput) std::clog << "add" << std::endl;
    00379                                 nlNodeVec.push_back( new OSnLNodePlus() );
    00380                         } break;
    00381                         case nlAddV: { // add variable
    00382                                 address = gmoGetjSolver(gmo, address);
    00383                                 if (debugoutput) std::clog << "add variable " << osinstance->getVariableNames()[address] << std::endl;
    00384                                 OSnLNodeVariable *nlNode = new OSnLNodeVariable();
    00385                                 nlNode->idx=address;
    00386                                 nlNodeVec.push_back( nlNode );
    00387                                 nlNodeVec.push_back( new OSnLNodePlus() );
    00388                         } break;
    00389                         case nlAddI: { // add immediate
    00390                                 if (debugoutput) std::clog << "add constant " << constants[address] << std::endl;
    00391                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00392                                 nlNode->value = constants[address];
    00393                                 nlNodeVec.push_back( nlNode );
    00394                                 nlNodeVec.push_back( new OSnLNodePlus() );
    00395                         } break;
    00396                         case nlSub: { // minus
    00397                                 if (debugoutput) std::clog << "minus" << std::endl;
    00398                                 nlNodeVec.push_back( new OSnLNodeMinus() );
    00399                         } break;
    00400                         case nlSubV: { // subtract variable
    00401                                 address = gmoGetjSolver(gmo, address);
    00402                                 if (debugoutput) std::clog << "substract variable " << osinstance->getVariableNames()[address] << std::endl;
    00403                                 OSnLNodeVariable *nlNode = new OSnLNodeVariable();
    00404                                 nlNode->idx=address;
    00405                                 nlNodeVec.push_back( nlNode );
    00406                                 nlNodeVec.push_back( new OSnLNodeMinus() );
    00407                         } break;
    00408                         case nlSubI: { // subtract immediate
    00409                                 if (debugoutput) std::clog << "substract constant " << constants[address] << std::endl;
    00410                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00411                                 nlNode->value = constants[address];
    00412                                 nlNodeVec.push_back( nlNode );
    00413                                 nlNodeVec.push_back( new OSnLNodeMinus() );
    00414                         } break;
    00415                         case nlMul: { // multiply
    00416                                 if (debugoutput) std::clog << "multiply" << std::endl;
    00417                                 nlNodeVec.push_back( new OSnLNodeTimes() );
    00418                         } break;
    00419                         case nlMulV: { // multiply variable
    00420                                 address = gmoGetjSolver(gmo, address);
    00421                                 if (debugoutput) std::clog << "multiply variable " << osinstance->getVariableNames()[address] << std::endl;
    00422                                 OSnLNodeVariable *nlNode = new OSnLNodeVariable();
    00423                                 nlNode->idx=address;
    00424                                 nlNodeVec.push_back( nlNode );
    00425                                 nlNodeVec.push_back( new OSnLNodeTimes() );
    00426                         } break;
    00427                         case nlMulI: { // multiply immediate
    00428                                 if (debugoutput) std::clog << "multiply constant " << constants[address] << std::endl;
    00429                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00430                                 nlNode->value = constants[address];
    00431                                 nlNodeVec.push_back( nlNode );
    00432                                 nlNodeVec.push_back( new OSnLNodeTimes() );
    00433                         } break;
    00434                         case nlDiv: { // divide
    00435                                 if (debugoutput) std::clog << "divide" << std::endl;
    00436                                 nlNodeVec.push_back( new OSnLNodeDivide() );
    00437                         } break;
    00438                         case nlDivV: { // divide variable
    00439                                 address = gmoGetjSolver(gmo, address);
    00440                                 if (debugoutput) std::clog << "divide variable " << osinstance->getVariableNames()[address] << std::endl;
    00441                                 OSnLNodeVariable *nlNode = new OSnLNodeVariable();
    00442                                 nlNode->idx=address;
    00443                                 nlNodeVec.push_back( nlNode );
    00444                                 nlNodeVec.push_back( new OSnLNodeDivide() );
    00445                         } break;
    00446                         case nlDivI: { // divide immediate
    00447                                 if (debugoutput) std::clog << "divide constant " << constants[address] << std::endl;
    00448                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00449                                 nlNode->value = constants[address];
    00450                                 nlNodeVec.push_back( nlNode );
    00451                                 nlNodeVec.push_back( new OSnLNodeDivide() );
    00452                         } break;
    00453                         case nlUMin: { // unary minus
    00454                                 if (debugoutput) std::clog << "negate" << std::endl;
    00455                                 nlNodeVec.push_back( new OSnLNodeNegate() );
    00456                         } break;
    00457                         case nlUMinV: { // unary minus variable
    00458                                 address = gmoGetjSolver(gmo, address);
    00459                                 if (debugoutput) std::clog << "push negated variable " << osinstance->getVariableNames()[address] << std::endl;
    00460                                 OSnLNodeVariable *nlNode = new OSnLNodeVariable();
    00461                                 nlNode->idx = address;
    00462                                 nlNode->coef = -1.;
    00463                                 nlNodeVec.push_back( nlNode );
    00464                         } break;
    00465                         case nlCallArg1 :
    00466                         case nlCallArg2 :
    00467                         case nlCallArgN : {
    00468                                 if (debugoutput) std::clog << "call function ";
    00469                                 GamsFuncCode func = GamsFuncCode(address+1); // here the shift by one was not a good idea
    00470                                 switch (func) {
    00471                                         case fnmin : {
    00472                                                 if (debugoutput) std::clog << "min" << std::endl;
    00473                                                 nlNodeVec.push_back( new OSnLNodeMin() );
    00474                                         } break;
    00475                                         case fnmax : {
    00476                                                 if (debugoutput) std::clog << "max" << std::endl;
    00477                                                 nlNodeVec.push_back( new OSnLNodeMax() );
    00478                                         } break;
    00479                                         case fnsqr : {
    00480                                                 if (debugoutput) std::clog << "square" << std::endl;
    00481                                                 nlNodeVec.push_back( new OSnLNodeSquare() );
    00482                                         } break;
    00483                                         case fnexp:
    00484                                         case fnslexp:
    00485                                         case fnsqexp: {
    00486                                                 if (debugoutput) std::clog << "exp" << std::endl;
    00487                                                 nlNodeVec.push_back( new OSnLNodeExp() );
    00488                                         } break;
    00489                                         case fnlog : {
    00490                                                 if (debugoutput) std::clog << "ln" << std::endl;
    00491                                                 nlNodeVec.push_back( new OSnLNodeLn() );
    00492                                         } break;
    00493                                         case fnlog10:
    00494                                         case fnsllog10:
    00495                                         case fnsqlog10: {
    00496                                                 if (debugoutput) std::clog << "log10 = ln * 1/ln(10)" << std::endl;
    00497                                                 nlNodeVec.push_back( new OSnLNodeLn() );
    00498                                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00499                                                 nlNode->value = 1./log(10.);
    00500                                                 nlNodeVec.push_back( nlNode );
    00501                                                 nlNodeVec.push_back( new OSnLNodeTimes() );
    00502                                         } break;
    00503                                         case fnlog2 : {
    00504                                                 if (debugoutput) std::clog << "log2 = ln * 1/ln(2)" << std::endl;
    00505                                                 nlNodeVec.push_back( new OSnLNodeLn() );
    00506                                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00507                                                 nlNode->value = 1./log(2.);
    00508                                                 nlNodeVec.push_back( nlNode );
    00509                                                 nlNodeVec.push_back( new OSnLNodeTimes() );
    00510                                         } break;
    00511                                         case fnsqrt: {
    00512                                                 if (debugoutput) std::clog << "sqrt" << std::endl;
    00513                                                 nlNodeVec.push_back( new OSnLNodeSqrt() );
    00514                                         } break;
    00515                                         case fnabs: {
    00516                                                 if (debugoutput) std::clog << "abs" << std::endl;
    00517                                                 nlNodeVec.push_back( new OSnLNodeAbs() );
    00518                                         } break;
    00519                                         case fncos: {
    00520                                                 if (debugoutput) std::clog << "cos" << std::endl;
    00521                                                 nlNodeVec.push_back( new OSnLNodeCos() );
    00522                                         } break;
    00523                                         case fnsin: {
    00524                                                 if (debugoutput) std::clog << "sin" << std::endl;
    00525                                                 nlNodeVec.push_back( new OSnLNodeSin() );
    00526                                         } break;
    00527                                         case fnpower:
    00528                                         case fnrpower: // x ^ y
    00529                                         case fncvpower: // constant ^ x
    00530                                         case fnvcpower: { // x ^ constant {
    00531                                                 if (debugoutput) std::clog << "power" << std::endl;
    00532                                                 nlNodeVec.push_back( new OSnLNodePower() );
    00533                                         } break;
    00534                                         case fnpi: {
    00535                                                 if (debugoutput) std::clog << "pi" << std::endl;
    00536                                                 nlNodeVec.push_back( new OSnLNodePI() );
    00537                                         } break;
    00538                                         case fndiv:
    00539                                         case fndiv0: {
    00540                                                 nlNodeVec.push_back( new OSnLNodeDivide() );
    00541                                         } break;
    00542                                         case fnslrec: // 1/x
    00543                                         case fnsqrec: { // 1/x
    00544                                                 if (debugoutput) std::clog << "divide" << std::endl;
    00545                                                 nlNodeVec.push_back( new OSnLNodeLn() );
    00546                                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00547                                                 nlNode->value = 1.;
    00548                                                 nlNodeVec.push_back( nlNode );
    00549                                                 nlNodeVec.push_back( new OSnLNodeDivide() );
    00550                                         } break;
    00551                                         case fnceil: case fnfloor: case fnround:
    00552                                         case fnmod: case fntrunc: case fnsign:
    00553                                         case fnarctan: case fnerrf: case fndunfm:
    00554                                         case fndnorm: case fnerror: case fnfrac: case fnerrorl:
    00555                             case fnfact /* factorial */:
    00556                             case fnunfmi /* uniform random number */:
    00557                             case fnncpf /* fischer: sqrt(x1^2+x2^2+2*x3) */:
    00558                             case fnncpcm /* chen-mangasarian: x1-x3*ln(1+exp((x1-x2)/x3))*/:
    00559                             case fnentropy /* x*ln(x) */: case fnsigmoid /* 1/(1+exp(-x)) */:
    00560                             case fnboolnot: case fnbooland:
    00561                             case fnboolor: case fnboolxor: case fnboolimp:
    00562                             case fnbooleqv: case fnrelopeq: case fnrelopgt:
    00563                             case fnrelopge: case fnreloplt: case fnrelople:
    00564                             case fnrelopne: case fnifthen:
    00565                             case fnedist /* euclidian distance */:
    00566                             case fncentropy /* x*ln((x+d)/(y+d))*/:
    00567                             case fngamma: case fnloggamma: case fnbeta:
    00568                             case fnlogbeta: case fngammareg: case fnbetareg:
    00569                             case fnsinh: case fncosh: case fntanh:
    00570                             case fnsignpower /* sign(x)*abs(x)^c */:
    00571                             case fnncpvusin /* veelken-ulbrich */:
    00572                             case fnncpvupow /* veelken-ulbrich */:
    00573                             case fnbinomial:
    00574                             case fntan: case fnarccos:
    00575                             case fnarcsin: case fnarctan2 /* arctan(x2/x1) */:
    00576                             case fnpoly: /* simple polynomial */
    00577                                         default : {
    00578                                                 if (debugoutput) std::cerr << "nr. " << func << " - unsuppored. Error." << std::endl;
    00579                                                 return NULL;
    00580                                         }
    00581                                 }
    00582                         } break;
    00583                         case nlMulIAdd: {
    00584                                 if (debugoutput) std::clog << "multiply constant " << constants[address] << " and add " << std::endl;
    00585                                 OSnLNodeNumber *nlNode = new OSnLNodeNumber();
    00586                                 nlNode->value = constants[address];
    00587                                 nlNodeVec.push_back( nlNode );
    00588                                 nlNodeVec.push_back( new OSnLNodeTimes() );
    00589                                 nlNodeVec.push_back( new OSnLNodePlus() );
    00590                         } break;
    00591                         case nlFuncArgN : {
    00592                                 if (debugoutput) std::clog << "ignored" << std::endl;
    00593                         } break;
    00594                         case nlArg: {
    00595                                 if (debugoutput) std::clog << "ignored" << std::endl;
    00596                         } break;
    00597                         case nlHeader: { // header
    00598                                 if (debugoutput) std::clog << "ignored" << std::endl;
    00599                         } break;
    00600                         case nlPushZero: {
    00601                                 if (debugoutput) std::clog << "push constant zero" << std::endl;
    00602                                 nlNodeVec.push_back( new OSnLNodeNumber() );
    00603                         } break;
    00604                         case nlStoreS: { // store scaled row
    00605                                 if (debugoutput) std::clog << "ignored" << std::endl;
    00606                         } break;
    00607                         // the following three should have been taken out by reorderInstr above; the remaining ones seem to be unused by now
    00608                         case nlPushS: // duplicate value from address levels down on top of stack
    00609                         case nlPopup: // duplicate value from this level to at address levels down and pop entries in between
    00610                         case nlSwap: // swap two positions on top of stack
    00611                         case nlAddL: // add local
    00612                         case nlSubL: // subtract local
    00613                         case nlMulL: // multiply local
    00614                         case nlDivL: // divide local
    00615                         case nlPushL: // push local
    00616                         case nlPopL: // pop local
    00617                         case nlPopDeriv: // pop derivative
    00618                         case nlUMinL: // push umin local
    00619                         case nlPopDerivS: // store scaled gradient
    00620                         case nlEquScale: // equation scale
    00621                         case nlEnd: // end of instruction list
    00622                         default: {
    00623                                 std::cerr << "not supported - Error." << std::endl;
    00624                                 return NULL;
    00625                         }
    00626                 }
    00627         }
    00628
    00629         if (!nlNodeVec.size()) return NULL;
    00630         // the vector is in postfix format - create expression tree and return it
    00631         return nlNodeVec[0]->createExpressionTreeFromPostfix(nlNodeVec);
    00632 }