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

- File List
- File Members

# /home/coin/SVN-release/OS-2.3.5/ApplicationTemplates/osDip/solvers/OSDipBlockBearcatSolver.cpp

Go to the documentation of this file.

    00001 /* $Id: OSBlockBearcatSolver.cpp 3038 2009-11-07 11:43:44Z kmartin $ */
    00017 #include "OSDipBlockBearcatSolver.h"
    00018 #include "OSErrorClass.h"
    00019 #include "OSDataStructures.h"
    00020
    00021 #include <sstream>
    00022 using std::ostringstream;
    00023
    00024
    00025
    00026
    00027 OSDipBlockBearcatSolver::OSDipBlockBearcatSolver():
    00028         m_osinstance(NULL) {
    00029         std::cout << "INSIDE OSDipBlockBearcatSolver CONSTRUCTOR" << std::endl;
    00030 }//end OSDipBlockBearcatSolver
    00031
    00032
    00033 OSDipBlockBearcatSolver::OSDipBlockBearcatSolver( OSInstance *osinstance,  OSOption *osoption) {
    00034           std::cout << "INSIDE OSDipBlockBearcatSolver CONSTRUCTOR" << std::endl;
    00035           std::cout << "whichBlock = " << m_whichBlock<< std::endl;
    00036         try{
    00037                 m_osinstance = osinstance;
    00038                 m_numberOfVar = m_osinstance->getVariableNumber();
    00039
    00040                 m_osoption = osoption;
    00041
    00042                 m_demand = NULL;
    00043                 m_u = NULL;
    00044                 m_v = NULL;
    00045                 m_g = NULL;
    00046                 m_px = NULL;
    00047                 m_tx =NULL;
    00048                 m_varIdx = NULL;
    00049
    00050                 //now get data
    00051                 std::vector<SolverOption*> solverOptions;
    00052                 std::vector<SolverOption*>::iterator vit;
    00053                 std::vector<int >demand;
    00054
    00055
    00056                 solverOptions = m_osoption->getSolverOptions("Dip");
    00057                 //iterate over the vector
    00058
    00059                 int tmpVal;
    00060
    00061                 for (vit = solverOptions.begin(); vit != solverOptions.end(); vit++) {
    00062
    00063
    00064                         //std::cout << (*vit)->name << std::endl;
    00065
    00066                         //(*vit)->name.compare("initialCol") == 0
    00067                         //if(rowNames[ i3].find("routeCapacity(1)") == string::npos )
    00068
    00069                         if( (*vit)->name.find("numHubs") !=  std::string::npos){
    00070
    00071                                 std::istringstream buffer( (*vit)->value);
    00072                                 buffer >> m_numHubs;
    00073                                 std::cout << "numHubs = " << m_numHubs <<  std::endl;
    00074
    00075                         }else{
    00076
    00077                                 if((*vit)->name.find("numNodes") !=  std::string::npos){
    00078
    00079                                         std::istringstream buffer( (*vit)->value);
    00080                                         buffer >> m_numNodes;
    00081                                         std::cout << "numNodes = " <<  m_numNodes <<  std::endl;
    00082
    00083                                 }else{
    00084                                         if((*vit)->name.find("totalDemand") !=  std::string::npos){
    00085
    00086                                                 std::istringstream buffer( (*vit)->value);
    00087                                                 buffer >> m_totalDemand;
    00088                                                 std::cout << "m_totalDemand = " << m_totalDemand <<  std::endl;
    00089
    00090                                         }else{
    00091                                                 if((*vit)->name.find("minDemand") !=  std::string::npos){
    00092
    00093                                                         std::istringstream buffer( (*vit)->value);
    00094                                                         buffer >> m_minDemand;
    00095                                                         std::cout << "m_minDemand = " << m_minDemand <<  std::endl;
    00096
    00097                                                 }else{
    00098                                                         if( (*vit)->name.find("demand") !=  std::string::npos ){
    00099
    00100                                                                 std::istringstream buffer( (*vit)->value);
    00101                                                                 buffer >> tmpVal;
    00102                                                                 demand.push_back( tmpVal);
    00103                                                                 //std::cout << "demand = " << tmpVal <<  std::endl;
    00104
    00105                                                         }
    00106                                                 }
    00107                                         }
    00108                                 }
    00109                         }//end if on solver options
    00110
    00111                 }//end for loop on options
    00112
    00113                 //now fill in demand
    00114                 m_demand = new int[ m_numNodes];
    00115
    00116                 m_varIdx = new int[ m_numNodes];
    00117
    00118                 m_u = new double*[ m_numNodes];
    00119                 m_v = new double*[ m_numNodes];
    00120                 m_g = new double*[ m_numNodes];
    00121
    00122                 m_px = new int*[ m_numNodes];
    00123                 m_tx = new int*[ m_numNodes];
    00124
    00125                 if(m_numNodes != demand.size( ) ) throw ErrorClass("inconsistent number of demand nodes");
    00126                 int i;
    00127                 for (i = 0; i < m_numNodes; i++) {
    00128
    00129                         m_demand[ i] = demand[i];
    00130                         m_u[ i] = new double[ m_totalDemand + 1];
    00131                         m_v[ i] = new double[ m_totalDemand + 1];
    00132                         m_g[ i] = new double[ m_numNodes];
    00133
    00134                         m_px[ i] = new int[ m_totalDemand + 1];
    00135                         m_tx[ i] = new int[ m_totalDemand + 1];
    00136
    00137
    00138                 }
    00139
    00140
    00141         } catch (const ErrorClass& eclass) {
    00142
    00143                 throw ErrorClass(eclass.errormsg);
    00144
    00145         }
    00146
    00147 }//end OSDipBlockBearcatSolver
    00148
    00149 OSDipBlockBearcatSolver::~OSDipBlockBearcatSolver(){
    00150
    00151         std::cout << "INSIDE ~OSDipBlockBearcatSolver()" << std::endl;
    00152
    00153
    00154         std::vector<IndexValuePair*>::iterator  vit;
    00155
    00156         for (vit = m_primalVals.begin(); vit != m_primalVals.end(); vit++) {
    00157
    00158                 delete *vit;
    00159         }
    00160
    00161         //delete data structures for arrays used in calculating minimum reduced cost
    00162         int i;
    00163
    00164         for(i = 0; i < m_numNodes; i++){
    00165
    00166             delete[] m_u[i];
    00167                 m_u[i] = NULL;
    00168
    00169             delete[] m_v[i];
    00170                 m_v[i] = NULL;
    00171
    00172             delete[] m_g[i];
    00173                 m_g[i] = NULL;
    00174
    00175             delete[] m_px[i];
    00176                 m_px[i] = NULL;
    00177
    00178             delete[] m_tx[i];
    00179                 m_tx[i] = NULL;
    00180         }
    00181
    00182         delete[] m_u;
    00183         m_u= NULL;
    00184
    00185         delete[] m_v;
    00186         m_v= NULL;
    00187
    00188         delete[] m_g;
    00189         m_g= NULL;
    00190
    00191         delete[] m_px;
    00192         m_px= NULL;
    00193
    00194         delete[] m_tx;
    00195         m_tx= NULL;
    00196
    00197
    00198
    00199         if(m_demand != NULL) delete[] m_demand;
    00200
    00201         if(m_varIdx != NULL) delete[] m_varIdx;
    00202         //if(m_osrlreader != NULL) delete m_osrlreader;
    00203 }//end ~OSDipBlockBearcatSolver
    00204
    00205 void OSDipBlockBearcatSolver::solve(double *cost, std::vector<IndexValuePair*> *solIndexValPair, double *optVal){
    00206
    00207         try{
    00208
    00209                 int i;
    00210
    00211
    00212
    00213                 struct IndexValuePair *primalValPair;
    00214
    00215                 std::vector<IndexValuePair*>::iterator  vit;
    00216
    00217                 for (vit = m_primalVals.begin(); vit != m_primalVals.end(); vit++) {
    00218
    00219                         delete *vit;
    00220                 }
    00221                 m_primalVals.clear();
    00222
    00223
    00224                 //now, lets find the minimum reduced cost
    00225
    00226                 *optVal = OSDBL_MAX;
    00227
    00228                 //std::cout << "Number of Var " << m_numberOfVar  << std::endl;
    00229
    00230                 int l ;
    00231
    00232                 int bestl;
    00233                 int kountVar;
    00234                 double testVal;
    00235
    00236                 for(l = m_minDemand; l <= m_totalDemand; l++){
    00237
    00238                         testVal = qrouteCost(m_whichBlock,  l,  cost,  &kountVar);
    00239                         if( testVal < *optVal){
    00240                                 *optVal = testVal;
    00241                                 bestl = l;
    00242                         }
    00243
    00244
    00245
    00246                 }
    00247
    00248                 *optVal = qrouteCost(m_whichBlock,  bestl,  cost,  &kountVar);
    00249                 std::cout <<   "best reduced cost = " << *optVal << std::endl;
    00250
    00251                 std::map<int, int> indexCount;
    00252
    00253                 std::map<int, int>::iterator mit;
    00254
    00255
    00256                 for(i = 0; i < kountVar; i++){
    00257
    00258                         if( indexCount.find(  m_varIdx[ i]) == indexCount.end()  ){
    00259
    00260                                 indexCount[ m_varIdx[ i]] = 1;
    00261
    00262                         }else{
    00263
    00264                                 indexCount[ m_varIdx[ i]] += 1;
    00265
    00266                         }
    00267
    00268
    00269                 }
    00270
    00271                 for (mit = indexCount.begin(); mit != indexCount.end(); mit++){
    00272
    00273                         //std::cout << "Variable Index " << mit->first << " Count =  " << mit->second << std::endl;
    00274
    00275                         //std::cout << "Variable Name = " << m_osinstance->getVariableNames()[   mit->first]  << std::endl;
    00276
    00277                         // get unique indexes
    00278                         primalValPair = new IndexValuePair();
    00279
    00280                         primalValPair->value =  mit->second;
    00281                         primalValPair->idx =  mit->first;
    00282
    00283
    00284                         m_primalVals.push_back( primalValPair);
    00285                 }
    00286
    00287                 /*** old way could lead to duplicate indexes
    00288                 for (i = 0; i < kountVar; i++){
    00289
    00290                         //std::cout << "Variable Number = " << m_varIdx[ i] << std::endl;
    00291                         std::cout << "Variable Name = " << m_osinstance->getVariableNames()[  m_varIdx[ i]]  << std::endl;
    00292
    00293                         // get unique indexes
    00294                         primalValPair = new IndexValuePair();
    00295
    00296                         primalValPair->value = 1.0;
    00297                         primalValPair->idx = m_varIdx[ i];
    00298
    00299
    00300                         m_primalVals.push_back( primalValPair);
    00301
    00302
    00303                 }
    00304         */
    00305
    00306                  *solIndexValPair = m_primalVals;
    00307
    00308
    00309         } catch (const ErrorClass& eclass) {
    00310
    00311                 throw ErrorClass(eclass.errormsg);
    00312
    00313         }
    00314
    00315 }//end solve
    00316
    00317
    00318 void OSDipBlockBearcatSolver::solve(double *cost, std::string *osrl){
    00319
    00320
    00321         try{
    00322                 //not implemented
    00323
    00324         } catch (const ErrorClass& eclass) {
    00325
    00326                 throw ErrorClass(eclass.errormsg);
    00327
    00328         }
    00329
    00330 }//end solve
    00331
    00332
    00333
    00334 double OSDipBlockBearcatSolver::qrouteCost(const int& k, const int& l, double* c, int* kountVar){
    00335
    00336         //critical -- nodes 0, ..., m_numNodes - 1 are the hub nodes
    00337         // we are doing the calculation for hub k, k <= m_numNodes - 1
    00338         //l  should be the total demand on the network -- we are NOT using 0 based counting
    00339         double rcost;
    00340         rcost = OSDBL_MAX;
    00341
    00342         //start of the cost vector for hub k plus move to start of l demands
    00343         // now move the pointer up
    00344         //int startPnt = k*totalDemand*(m_numNodes*m_numNodes - m_numNodes) + (l - 1)*(m_numNodes*m_numNodes - m_numNodes);
    00345
    00346         int startPnt = (l - 1)*(m_numNodes*m_numNodes - m_numNodes);
    00347         c+=  startPnt ;
    00348
    00349         *kountVar = 0;
    00350         int bestLastNode;
    00351         int i;
    00352         int j;
    00353         int l1;
    00354         int l2;
    00355         //for(i = 0; i < 20; i++){
    00356         //      std::cout << "i =  " << i  <<  " c[i] = " << *(c + i) << std::endl ;
    00357         //}
    00358
    00359
    00360
    00361         // l  is the total demand on this network
    00362         //address of node (j, i) is j*(m_numNodes-1) + i when i < j
    00363         //address of node (j, i) is j*(m_numNodes-1) + i - 1 when i > j
    00364
    00365         //
    00366         // initialize
    00367
    00368         for(i = m_numHubs; i < m_numNodes; i++){
    00369
    00370
    00371                 for(l1 = m_minDemand; l1 <= l; l1++){  //l-1  is total demand on network
    00372
    00373                         m_u[i][l1] = OSDBL_MAX;
    00374                         m_v[i][l1] = OSDBL_MAX;
    00375
    00376                         if(l1 == *(m_demand + i) ){
    00377
    00378                                 m_px[i][l1] = k;
    00379                                 // want the cost for arc (k, i)
    00380                                 // note: k < i so use that formula
    00381                                 m_u[i][l1] = *(c + k*(m_numNodes - 1) + i - 1);  // put in correct cost
    00382
    00383
    00384                         }
    00385                 }
    00386         }
    00387         //end initialize
    00388         //
    00389
    00390         //if l = minDemand we visit only one node, let's get the reduced cost in this case
    00391         if(l == m_minDemand){
    00392
    00393                 for(i = m_numHubs; i < m_numNodes; i++){
    00394
    00395
    00396                         if(  m_u[i][l] + *(c + i*(m_numNodes-1) + k )  < rcost){
    00397
    00398                                 rcost = m_u[i][l] + *(c + i*(m_numNodes-1) + k );
    00399
    00400                                 //std::cout << " m_u[i][l2] = "  << m_u[i][l2] << std::endl;
    00401
    00402                                 //std::cout << " *(c + i*(m_numNodes-1) + k ) = "  << *(c + i*(m_numNodes-1) + k ) << std::endl;
    00403                                 //push node back
    00404                                 bestLastNode = i;
    00405                         }
    00406
    00407                 }
    00408
    00409                 //go from node bestLastNode to node k
    00410                 //node bestLastNode is a higher number than k
    00411                 *(m_varIdx + (*kountVar)++) = startPnt + bestLastNode*(m_numNodes - 1)  +  k ;
    00412                 *(m_varIdx + (*kountVar)++) = startPnt + k*(m_numNodes - 1)  + bestLastNode - 1;
    00413
    00414
    00415                 return rcost;
    00416         }//end if on l == minDemand
    00417
    00418 // now calculate values for demand 2 or greater
    00419         //address of node (j, i) is j*(m_numNodes-1) + i when i < j
    00420         //address of node (j, i) is j*(m_numNodes-1) + i - 1 when i > j
    00421         // we start l2 at 2 since demand must be at least 1
    00422         // change to min demand + 1
    00423         int lowerVal = m_minDemand + 1;
    00424         for(l2 = lowerVal; l2 <= l; l2++){// loop over possible demand values assuming we have already gone to at least one node
    00425
    00426                 for(i = m_numHubs; i < m_numNodes; i++) { //we are finding least cost to node i
    00427
    00428
    00429                         if( *(m_demand + i) < l2 ){ // kipp < OR <= ????
    00430
    00431                                 for(j = m_numHubs; j < i; j++){ //we are coming from node j
    00432
    00433
    00434
    00435                                         //calculate g  -- l2 - d[ i]  is the demand trough and including node j
    00436                                         l1 = l2 - *(m_demand + i);
    00437
    00438                                         if( m_px[j][ l1 ] != i ){//check to make sure we did not come into j from i
    00439
    00440
    00441                                                 m_g[j][i] = m_u[ j][ l1 ] + *(c + j*(m_numNodes-1) + i - 1) ;
    00442
    00443
    00444
    00445
    00446                                         }else{
    00447
    00448                                                 m_g[j][i] = m_v[ j][ l1] + *(c + j*(m_numNodes-1) + i - 1) ;
    00449
    00450
    00451
    00452                                         }
    00453
    00454                                         // update u and the pointer for p
    00455
    00456                                         if(m_g[j][i] < m_u[i][l2] ){
    00457
    00458                                                 m_u[i][l2] = m_g[j][i];
    00459                                                 m_px[i][l2] =  j;
    00460
    00461                                         }
    00462
    00463
    00464
    00465                                 }//end of first for loop on j, j < i
    00466
    00467
    00468                                 for(j = i + 1; j < m_numNodes; j++){ //we are coming from node j
    00469
    00470
    00471                                         //calculate g  -- l2 - d[ i]  is the demand trough and including node j
    00472                                         l1 = l2 - *(m_demand + i);
    00473
    00474                                         if( m_px[j][ l1 ] != i ){//check to make sure we did not come into j from i
    00475
    00476
    00477                                                 m_g[j][i] = m_u[ j][ l1 ] + *(c + j*(m_numNodes-1) + i ) ;
    00478
    00479
    00480                                         }else{
    00481
    00482                                                 m_g[j][i] = m_v[ j][ l1] + *(c + j*(m_numNodes-1) + i ) ;
    00483
    00484                                         }
    00485
    00486                                         // update u and the pointer for p
    00487
    00488                                         if(m_g[j][i] < m_u[i][l2] ){
    00489
    00490                                                 m_u[i][l2] = m_g[j][i];
    00491                                                 m_px[i][l2] =  j;
    00492
    00493                                         }
    00494
    00495
    00496                                 }//end of second for loop on j, j > i
    00497
    00498
    00499                                 //now calculate the second best solution and point
    00500                                 //right now m_px[i][l2] points to the best j node
    00501
    00502                                 for(j =m_numHubs; j < m_numNodes; j++){ //we are coming from node j
    00503
    00504                                         if(j != i){
    00505
    00506                                                 if( (m_g[j][i] < m_v[i][l2] ) && (m_px[i][l2] != j)  ){ // kipp the && gives the second best
    00507
    00508                                                         //if(g(j, i) < v(i, l2),
    00509
    00510                                                         m_v[i][l2] = m_g[j][i];
    00511                                                         m_tx[i][l2] = j;
    00512
    00513
    00514                                                 }
    00515
    00516                                         }
    00517
    00518
    00519                                 }//end second best calculation, anothe for loop on j
    00520
    00521                                 //now if l2 = l  we are done
    00522                                 if(l2 == l ){
    00523
    00524                                         if(  m_u[i][l2] + *(c + i*(m_numNodes-1) + k )  < rcost){
    00525
    00526                                                 rcost = m_u[i][l2] + *(c + i*(m_numNodes-1) + k );
    00527
    00528                                                 bestLastNode = i;
    00529                                         }
    00530
    00531                                 }
    00532
    00533
    00534                         }//end if on demand less than l2
    00535
    00536
    00537                 }//i loop
    00538
    00539
    00540         }//l2 loop
    00541
    00542
    00543         //std::cout << "best Last Node = "  << bestLastNode << std::endl;
    00544
    00545         // now get the path that gives the reduced cost
    00546
    00547         int currentNode;
    00548         int successorNode;
    00549         int lvalue;
    00550
    00551         //initialize
    00552         // we work backwords from bestLastNode
    00553         //in our recursion we recurse on the currentNode and assume
    00554         //it is in the optimal path
    00555
    00556         //node bestLastNode is a higher number than k
    00557         *(m_varIdx + (*kountVar)++) = startPnt + bestLastNode*(m_numNodes - 1)  +  k ;
    00558
    00559
    00560
    00561
    00562         //by successor, I mean node after current node in the path
    00563         successorNode = k;
    00564         currentNode = bestLastNode;
    00565         //the lvalue is the demand through the currentNode
    00566         lvalue = l ;
    00567
    00568
    00569         while(currentNode != k){
    00570
    00571
    00572                 if( m_px[ currentNode][ lvalue ] != successorNode){
    00573
    00574
    00575
    00576                         //update nodes
    00577                         successorNode = currentNode;
    00578                         currentNode = m_px[ currentNode][ lvalue ];
    00579
    00580
    00581                         if(currentNode - successorNode > 0){
    00582                                  //below diagonal
    00583                                 *(m_varIdx  + (*kountVar)++) = startPnt + currentNode*(m_numNodes - 1)  +  successorNode;
    00584
    00585                         }else{
    00586                                  //above diagonal
    00587                                 *(m_varIdx + (*kountVar)++) = startPnt + currentNode*(m_numNodes - 1)  +  successorNode  - 1 ;
    00588
    00589                         }
    00590
    00591
    00592                 }else{ //take second best
    00593
    00594
    00595                         //update nodes
    00596                         successorNode = currentNode;
    00597                         currentNode = m_tx[ currentNode][ lvalue ];
    00598
    00599                         if(currentNode - successorNode > 0){
    00600                                  //below diagonal
    00601                                 *(m_varIdx  + (*kountVar)++) = startPnt + currentNode*(m_numNodes - 1)  +  successorNode;
    00602
    00603                         }else{
    00604                                  //above diagonal
    00605                                 *(m_varIdx + (*kountVar)++) = startPnt + currentNode*(m_numNodes - 1)  +  successorNode  - 1 ;
    00606
    00607                         }
    00608
    00609                 }
    00610
    00611                 //update lvalue
    00612                 lvalue = lvalue - *(m_demand + successorNode);
    00613
    00614
    00615
    00616         }
    00617
    00618
    00619         //go from node k to bestLastNode -- already done in loop above
    00620         //*(m_varIdx + (*kountVar)++) = startPnt + k*(m_numNodes - 1)  +  currentNode - 1;
    00621
    00622
    00623
    00624
    00625         return rcost;
    00626
    00627 }//end qroute
    00628