00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include "ColPackHeaders.h"
00022
00023 using namespace std;
00024
00025 namespace ColPack
00026 {
00027
00028 int GraphColoring::FindCycle(int i_Vertex, int i_AdjacentVertex, int i_DistanceOneVertex, int i_SetID, vector<int> & vi_CandidateColors, vector<int> & vi_FirstVisitedOne, vector<int> & vi_FirstVisitedTwo)
00029 {
00030 int i_VertexOne, i_VertexTwo;
00031
00032 if(i_SetID != _UNKNOWN)
00033 {
00034 i_VertexOne = vi_FirstVisitedOne[i_SetID];
00035 i_VertexTwo = vi_FirstVisitedTwo[i_SetID];
00036
00037 if(i_VertexOne != i_Vertex)
00038 {
00039 vi_FirstVisitedOne[i_SetID] = i_Vertex;
00040 vi_FirstVisitedTwo[i_SetID] = i_AdjacentVertex;
00041 }
00042 else
00043 if((i_VertexOne == i_Vertex) && (i_VertexTwo != i_AdjacentVertex))
00044 {
00045 vi_CandidateColors[m_vi_VertexColors[i_DistanceOneVertex]] = i_Vertex;
00046
00047 #if DEBUG == 1401
00048
00049 cout<<"DEBUG 1401 | Acyclic Coloring | Found Cycle | Vertex "<<STEP_UP(i_Vertex)<<endl;
00050 #endif
00051
00052 }
00053 }
00054
00055 return(_TRUE);
00056 }
00057
00058
00059
00060
00061 int GraphColoring::UpdateSet(int i_Vertex, int i_AdjacentVertex, int i_DistanceOneVertex, map< int, map<int, int> > & mimi2_VertexEdgeMap, vector<int> & vi_FirstSeenOne, vector<int> & vi_FirstSeenTwo, vector<int> & vi_FirstSeenThree)
00062 {
00063 int i_ColorID;
00064
00065 int i_VertexOne, i_VertexTwo, i_VertexThree;
00066
00067 i_ColorID = m_vi_VertexColors[i_AdjacentVertex];
00068
00069 i_VertexOne = vi_FirstSeenOne[i_ColorID];
00070 i_VertexTwo = vi_FirstSeenTwo[i_ColorID];
00071 i_VertexThree = vi_FirstSeenThree[i_ColorID];
00072
00073 if(i_VertexOne != i_Vertex)
00074 {
00075 vi_FirstSeenOne[i_ColorID] = i_Vertex;
00076 vi_FirstSeenTwo[i_ColorID] = i_AdjacentVertex;
00077 vi_FirstSeenThree[i_ColorID] = i_DistanceOneVertex;
00078 }
00079 else
00080 {
00081 if(i_VertexTwo < i_VertexThree)
00082 {
00083 return(mimi2_VertexEdgeMap[i_VertexTwo][i_VertexThree]);
00084 }
00085 else
00086 {
00087 return(mimi2_VertexEdgeMap[i_VertexThree][i_VertexTwo]);
00088 }
00089 }
00090
00091 return(_UNKNOWN);
00092 }
00093
00094
00095
00096 int GraphColoring::SearchDepthFirst(int i_RootVertex, int i_ParentVertex, int i_Vertex, vector<int> & vi_TouchedVertices)
00097 {
00098 int i;
00099
00100 int i_VertexCount;
00101
00102 int i_ViolationCount;
00103
00104 i_ViolationCount = _FALSE;
00105
00106 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00107
00108 for(i=m_vi_Vertices[i_Vertex]; i<m_vi_Vertices[STEP_UP(i_Vertex)]; i++)
00109 {
00110 if(m_vi_Edges[i] == i_ParentVertex)
00111 {
00112 continue;
00113 }
00114
00115 if(m_vi_Edges[i] == i_RootVertex)
00116 {
00117 i_ViolationCount++;
00118
00119 if(i_ViolationCount == _TRUE)
00120 {
00121 cout<<endl;
00122 cout<<"Acyclic Coloring | Violation Check | "<<m_s_InputFile<<endl;
00123 cout<<endl;
00124 }
00125
00126 cout<<"Violation "<<i_ViolationCount<<"\t : "<<STEP_UP(i_RootVertex)<<" ["<<STEP_UP(m_vi_VertexColors[i_RootVertex])<<"] ... "<<STEP_UP(i_ParentVertex)<<" ["<<STEP_UP(m_vi_VertexColors[i_ParentVertex])<<"] - "<<STEP_UP(i_Vertex)<<" ["<<STEP_UP(m_vi_VertexColors[i_Vertex])<<"] - "<<STEP_UP(m_vi_Edges[i])<<" ["<<STEP_UP(m_vi_VertexColors[m_vi_Edges[i]])<<"]"<<endl;
00127 }
00128
00129 if(m_vi_VertexColors[m_vi_Edges[i]] == m_vi_VertexColors[i_Vertex])
00130 {
00131 i_ViolationCount++;
00132
00133 if(i_ViolationCount == _TRUE)
00134 {
00135 cout<<endl;
00136 cout<<"Acyclic Coloring | Violation Check | "<<m_s_InputFile<<endl;
00137 cout<<endl;
00138 }
00139
00140 cout<<"Violation "<<i_ViolationCount<<"\t : "<<STEP_UP(i_Vertex)<<" ["<<STEP_UP(m_vi_VertexColors[i_Vertex])<<"] - "<<STEP_UP(m_vi_Edges[i])<<" ["<<STEP_UP(m_vi_VertexColors[m_vi_Edges[i]])<<"]"<<endl;
00141
00142 }
00143
00144 if(vi_TouchedVertices[m_vi_Edges[i]] == _TRUE)
00145 {
00146 continue;
00147 }
00148
00149 if(m_vi_VertexColors[m_vi_Edges[i]] != m_vi_VertexColors[i_ParentVertex])
00150 {
00151 continue;;
00152 }
00153
00154 vi_TouchedVertices[m_vi_Edges[i]] = _TRUE;
00155
00156 i_ViolationCount = SearchDepthFirst(i_RootVertex, i_Vertex, m_vi_Edges[i], vi_TouchedVertices);
00157
00158 }
00159
00160 return(i_ViolationCount);
00161
00162 }
00163
00164
00165
00166 int GraphColoring::CheckVertexColoring(string s_GraphColoringVariant)
00167 {
00168 if(m_s_VertexColoringVariant.compare(s_GraphColoringVariant) == 0)
00169 {
00170 return(_TRUE);
00171 }
00172
00173 if(m_s_VertexColoringVariant.compare("ALL") != 0)
00174 {
00175 m_s_VertexColoringVariant = s_GraphColoringVariant;
00176 }
00177
00178 if(m_s_VertexOrderingVariant.empty())
00179 {
00180 NaturalOrdering();
00181 }
00182
00183 return(_FALSE);
00184 }
00185
00186
00187
00188 int GraphColoring::CalculateVertexColorClasses()
00189 {
00190 if(m_s_VertexColoringVariant.empty())
00191 {
00192 return(_FALSE);
00193 }
00194
00195 int i_TotalVertexColors = STEP_UP(m_i_VertexColorCount);
00196
00197 m_vi_VertexColorFrequency.clear();
00198 m_vi_VertexColorFrequency.resize((unsigned) i_TotalVertexColors, _FALSE);
00199
00200 int i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00201
00202 for(int i = 0; i < i_VertexCount; i++)
00203 {
00204 m_vi_VertexColorFrequency[m_vi_VertexColors[i]]++;
00205 }
00206
00207 for(int i = 0; i < i_TotalVertexColors; i++)
00208 {
00209 if(m_i_LargestColorClassSize < m_vi_VertexColorFrequency[i])
00210 {
00211 m_i_LargestColorClass = i;
00212
00213 m_i_LargestColorClassSize = m_vi_VertexColorFrequency[i];
00214
00215 }
00216
00217 if(m_i_SmallestColorClassSize == _UNKNOWN)
00218 {
00219 m_i_SmallestColorClass = i;
00220
00221 m_i_SmallestColorClassSize = m_vi_VertexColorFrequency[i];
00222 }
00223 else
00224 if(m_i_SmallestColorClassSize > m_vi_VertexColorFrequency[i])
00225 {
00226 m_i_SmallestColorClass = i;
00227
00228 m_i_SmallestColorClassSize = m_vi_VertexColorFrequency[i];
00229 }
00230 }
00231
00232 m_d_AverageColorClassSize = i_TotalVertexColors / i_VertexCount;
00233
00234 return(_TRUE);
00235 }
00236
00237
00238
00239 GraphColoring::GraphColoring() : GraphOrdering()
00240 {
00241 Clear();
00242
00243 Seed_init();
00244 }
00245
00246
00247
00248 GraphColoring::~GraphColoring()
00249 {
00250 Clear();
00251
00252 Seed_reset();
00253 }
00254
00255
00256 void GraphColoring::Clear()
00257 {
00258 GraphOrdering::Clear();
00259
00260 m_i_VertexColorCount = _UNKNOWN;
00261
00262 m_i_LargestColorClass = _UNKNOWN;
00263 m_i_SmallestColorClass = _UNKNOWN;
00264
00265 m_i_LargestColorClassSize = _UNKNOWN;
00266 m_i_SmallestColorClassSize = _UNKNOWN;
00267
00268 m_d_AverageColorClassSize = _UNKNOWN;
00269
00270 m_i_ColoringUnits = _UNKNOWN;
00271
00272 m_d_ColoringTime = _UNKNOWN;
00273 m_d_CheckingTime = _UNKNOWN;
00274
00275 m_s_VertexColoringVariant.clear();
00276
00277 m_vi_VertexColors.clear();
00278
00279 m_vi_VertexColorFrequency.clear();
00280
00281
00282 return;
00283 }
00284
00285
00286
00287 int GraphColoring::DistanceOneColoring()
00288 {
00289
00290
00291
00292
00293
00294 int i, j;
00295
00296 int i_PresentVertex;
00297
00298 int i_VertexCount;
00299
00300 vector<int> vi_CandidateColors;
00301
00302 m_i_VertexColorCount = _UNKNOWN;
00303
00304 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00305
00306 m_vi_VertexColors.clear();
00307 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00308
00309 vi_CandidateColors.clear();
00310 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00311
00312 for(i=0; i<i_VertexCount; i++)
00313 {
00314 i_PresentVertex = m_vi_OrderedVertices[i];
00315
00316 #if VERBOSE == _TRUE
00317
00318 cout<<"DEBUG 1454 | Distance One Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
00319
00320 #endif
00321
00322 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
00323 {
00324 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
00325 {
00326 continue;
00327 }
00328
00329 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
00330
00331 }
00332
00333 for(j=0; j<i_VertexCount; j++)
00334 {
00335 if(vi_CandidateColors[j] != i_PresentVertex)
00336 {
00337 m_vi_VertexColors[i_PresentVertex] = j;
00338
00339 if(m_i_VertexColorCount < j)
00340 {
00341 m_i_VertexColorCount = j;
00342 }
00343
00344 break;
00345 }
00346 }
00347 }
00348
00349 return(_TRUE);
00350
00351 }
00352
00353
00354
00355 int GraphColoring::DistanceTwoColoring()
00356 {
00357 if(CheckVertexColoring("DISTANCE TWO"))
00358 {
00359 return(_TRUE);
00360 }
00361
00362 int i, j, k;
00363
00364 int i_PresentVertex;
00365
00366 int i_VertexCount;
00367
00368 vector<int> vi_CandidateColors;
00369
00370 m_i_VertexColorCount = _UNKNOWN;
00371
00372 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00373
00374 m_vi_VertexColors.clear();
00375 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00376
00377 vi_CandidateColors.clear();
00378 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00379
00380 for(i=0; i<i_VertexCount; i++)
00381 {
00382 i_PresentVertex = m_vi_OrderedVertices[i];
00383
00384 #if VERBOSE == _TRUE
00385
00386 cout<<"DEBUG 1455 | Distance Two Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
00387
00388 #endif
00389 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
00390 {
00391
00392
00393
00394
00395
00396
00397
00398 if(m_vi_VertexColors[m_vi_Edges[j]] != _UNKNOWN) vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
00399
00400 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
00401 {
00402
00403
00404
00405
00406
00407
00408
00409
00410
00411 if(m_vi_VertexColors[m_vi_Edges[k]] != _UNKNOWN)
00412 {
00413 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00414 }
00415 }
00416 }
00417
00418 for(j=0; j<i_VertexCount; j++)
00419 {
00420 if(vi_CandidateColors[j] != i_PresentVertex)
00421 {
00422 m_vi_VertexColors[i_PresentVertex] = j;
00423
00424 if(m_i_VertexColorCount < j)
00425 {
00426 m_i_VertexColorCount = j;
00427 }
00428
00429 break;
00430 }
00431 }
00432 }
00433
00434 return(_TRUE);
00435 }
00436
00437
00438
00439 int GraphColoring::NaiveStarColoring()
00440 {
00441 if(CheckVertexColoring("NAIVE STAR"))
00442 {
00443 return(_TRUE);
00444 }
00445
00446 int i, j, k, l;
00447
00448 int i_PresentVertex;
00449
00450 int i_VertexCount;
00451
00452 vector<int> vi_CandidateColors;
00453
00454 m_i_VertexColorCount = _UNKNOWN;
00455
00456 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00457
00458 m_vi_VertexColors.clear();
00459 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00460
00461 vi_CandidateColors.clear();
00462 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00463
00464 for(i=0; i<i_VertexCount; i++)
00465 {
00466 i_PresentVertex = m_vi_OrderedVertices[i];
00467
00468 #if VERBOSE == _TRUE
00469
00470 cout<<"DEBUG 1456 | Naive Star Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
00471
00472 #endif
00473
00474 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
00475 {
00476 if(m_vi_VertexColors[m_vi_Edges[j]] != _UNKNOWN)
00477 {
00478 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
00479 }
00480
00481 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
00482 {
00483 if(m_vi_Edges[k] == i_PresentVertex)
00484 {
00485 continue;
00486 }
00487
00488 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
00489 {
00490 continue;
00491 }
00492
00493 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
00494 {
00495 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00496 }
00497 else
00498 {
00499 for(l=m_vi_Vertices[m_vi_Edges[k]]; l<m_vi_Vertices[STEP_UP(m_vi_Edges[k])]; l++)
00500 {
00501 if(m_vi_Edges[l] == m_vi_Edges[j])
00502 {
00503 continue;
00504 }
00505
00506 if(m_vi_VertexColors[m_vi_Edges[l]] == _UNKNOWN)
00507 {
00508 continue;
00509 }
00510
00511 if(m_vi_VertexColors[m_vi_Edges[l]] == m_vi_VertexColors[m_vi_Edges[j]])
00512 {
00513 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00514
00515 break;
00516 }
00517 }
00518 }
00519 }
00520 }
00521
00522 for(j=0; j<i_VertexCount; j++)
00523 {
00524 if(vi_CandidateColors[j] != i_PresentVertex)
00525 {
00526 m_vi_VertexColors[i_PresentVertex] = j;
00527
00528 if(m_i_VertexColorCount < j)
00529 {
00530 m_i_VertexColorCount = j;
00531 }
00532
00533 break;
00534 }
00535 }
00536 }
00537
00538 return(_TRUE);
00539
00540 }
00541
00542
00543
00544 int GraphColoring::RestrictedStarColoring()
00545 {
00546 if(CheckVertexColoring("RESTRICTED STAR"))
00547 {
00548 return(_TRUE);
00549 }
00550
00551 int i, j, k;
00552
00553 int i_PresentVertex;
00554
00555 int i_VertexCount;
00556
00557 vector<int> vi_CandidateColors;
00558
00559 m_i_VertexColorCount = _UNKNOWN;
00560
00561 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00562
00563 m_vi_VertexColors.clear();
00564 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00565
00566 vi_CandidateColors.clear();
00567 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00568
00569 for(i=0; i<i_VertexCount; i++)
00570 {
00571
00572 i_PresentVertex = m_vi_OrderedVertices[i];
00573
00574 #if VERBOSE == _TRUE
00575
00576 cout<<"DEBUG 1457 | Restricted Star Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
00577
00578 #endif
00579
00580 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
00581 {
00582 if(m_vi_VertexColors[m_vi_Edges[j]] != _UNKNOWN)
00583 {
00584 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
00585 }
00586
00587 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
00588 {
00589 if(m_vi_Edges[k] == i_PresentVertex)
00590 {
00591 continue;
00592 }
00593
00594 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
00595 {
00596 continue;
00597 }
00598
00599 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
00600 {
00601 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00602 }
00603 else
00604 if(m_vi_VertexColors[m_vi_Edges[k]] < m_vi_VertexColors[m_vi_Edges[j]])
00605 {
00606 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00607 }
00608 }
00609 }
00610
00611 for(j=0; j<i_VertexCount; j++)
00612 {
00613 if(vi_CandidateColors[j] != i_PresentVertex)
00614 {
00615 m_vi_VertexColors[i_PresentVertex] = j;
00616
00617 if(m_i_VertexColorCount < j)
00618 {
00619 m_i_VertexColorCount = j;
00620 }
00621
00622 break;
00623 }
00624 }
00625 }
00626
00627 return(_TRUE);
00628
00629 }
00630
00631
00632
00633 int GraphColoring::StarColoring()
00634 {
00635 if(CheckVertexColoring("STAR"))
00636 {
00637 return(_TRUE);
00638 }
00639
00640 int i, j, k;
00641
00642 int _FOUND;
00643
00644 int i_ColorID, i_StarID;
00645
00646 int i_PresentVertex;
00647
00648 int i_VertexCount, i_EdgeCount;
00649
00650 int i_VertexOne, i_VertexTwo;
00651
00652 vector<int> vi_MemberEdges;
00653
00654 vector<int> vi_CandidateColors;
00655
00656 vector<int> vi_EdgeStarMap;
00657 vector<int> vi_StarHubMap;
00658
00659 vector<int> vi_FirstTreated;
00660
00661 vector<int> vi_FirstSeenOne, vi_FirstSeenTwo;
00662
00663 map< int, map<int, int> > mimi2_VertexEdgeMap;
00664
00665 m_i_VertexColorCount = _UNKNOWN;
00666
00667 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
00668
00669 i_EdgeCount = (signed) m_vi_Edges.size();
00670
00671 vi_EdgeStarMap.clear();
00672 vi_EdgeStarMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
00673
00674 vi_StarHubMap.clear();
00675 vi_StarHubMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
00676
00677 m_vi_VertexColors.clear();
00678 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00679
00680 vi_CandidateColors.clear();
00681 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
00682
00683 vi_FirstSeenOne.clear();
00684 vi_FirstSeenOne.resize((unsigned) i_VertexCount, _UNKNOWN);
00685
00686 vi_FirstSeenTwo.clear();
00687 vi_FirstSeenTwo.resize((unsigned) i_VertexCount, _UNKNOWN);
00688
00689
00690
00691
00692 vi_FirstTreated.clear();
00693 vi_FirstTreated.resize((unsigned) i_VertexCount, _UNKNOWN);
00694
00695 k = _FALSE;
00696
00697 for(i=0; i<i_VertexCount; i++)
00698 {
00699 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
00700 {
00701 if(i < m_vi_Edges[j])
00702 {
00703 mimi2_VertexEdgeMap[i][m_vi_Edges[j]] = k;
00704
00705 vi_EdgeStarMap[k] = k;
00706
00707 k++;
00708 }
00709 }
00710 }
00711
00712 #if VERBOSE == _TRUE
00713
00714 cout<<endl;
00715
00716 #endif
00717
00718 for(i=0; i<i_VertexCount; i++)
00719 {
00720 i_PresentVertex = m_vi_OrderedVertices[i];
00721
00722 #if VERBOSE == _TRUE
00723
00724 cout<<"DEBUG 1458 | Star Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
00725
00726 #endif
00727 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
00728 {
00729 i_ColorID = m_vi_VertexColors[m_vi_Edges[j]];
00730
00731 if(i_ColorID == _UNKNOWN)
00732 {
00733 continue;
00734 }
00735
00736 vi_CandidateColors[i_ColorID] = i_PresentVertex;
00737
00738 i_VertexOne = vi_FirstSeenOne[i_ColorID];
00739 i_VertexTwo = vi_FirstSeenTwo[i_ColorID];
00740
00741 if(i_VertexOne == i_PresentVertex)
00742 {
00743 if(vi_FirstTreated[i_VertexTwo] != i_PresentVertex)
00744 {
00745
00746 for(k=m_vi_Vertices[i_VertexTwo]; k<m_vi_Vertices[STEP_UP(i_VertexTwo)]; k++)
00747 {
00748 if(m_vi_Edges[k] == i_PresentVertex)
00749 {
00750 continue;
00751 }
00752
00753 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
00754 {
00755 continue;
00756 }
00757
00758 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00759
00760 }
00761
00762 vi_FirstTreated[i_VertexTwo] = i_PresentVertex;
00763 }
00764
00765 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
00766 {
00767 if(m_vi_Edges[k] == i_PresentVertex)
00768 {
00769 continue;
00770 }
00771
00772 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
00773 {
00774 continue;
00775 }
00776
00777 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00778
00779 }
00780
00781 vi_FirstTreated[m_vi_Edges[j]] = i_PresentVertex;
00782 }
00783 else
00784 {
00785 vi_FirstSeenOne[i_ColorID] = i_PresentVertex;
00786 vi_FirstSeenTwo[i_ColorID] = m_vi_Edges[j];
00787
00788 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
00789 {
00790 if(m_vi_Edges[k] == i_PresentVertex)
00791 {
00792 continue;
00793 }
00794
00795 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
00796 {
00797 continue;
00798 }
00799
00800 if(m_vi_Edges[j] < m_vi_Edges[k])
00801 {
00802 if(vi_StarHubMap[vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]]] == m_vi_Edges[k])
00803 {
00804 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00805 }
00806 }
00807 else
00808 {
00809 if(vi_StarHubMap[vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]]] == m_vi_Edges[k])
00810 {
00811 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
00812 }
00813 }
00814 }
00815 }
00816 }
00817
00818 for(j=0; j<i_VertexCount; j++)
00819 {
00820 if(vi_CandidateColors[j] != i_PresentVertex)
00821 {
00822 m_vi_VertexColors[i_PresentVertex] = j;
00823
00824 if(m_i_VertexColorCount < j)
00825 {
00826 m_i_VertexColorCount = j;
00827 }
00828
00829 break;
00830 }
00831 }
00832
00833 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
00834 {
00835 _FOUND = _FALSE;
00836
00837 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
00838 {
00839 continue;
00840 }
00841
00842 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
00843 {
00844 if(m_vi_Edges[k] == i_PresentVertex)
00845 {
00846 continue;
00847 }
00848
00849 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
00850 {
00851 continue;
00852 }
00853
00854 if(m_vi_VertexColors[m_vi_Edges[k]] == m_vi_VertexColors[i_PresentVertex])
00855 {
00856 _FOUND = _TRUE;
00857
00858 if(m_vi_Edges[j] < m_vi_Edges[k])
00859 {
00860 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]];
00861
00862 vi_StarHubMap[i_StarID] = m_vi_Edges[j];
00863
00864 if(i_PresentVertex < m_vi_Edges[j])
00865 {
00866 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
00867 }
00868 else
00869 {
00870 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
00871 }
00872 }
00873 else
00874 {
00875 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]];
00876
00877 vi_StarHubMap[i_StarID] = m_vi_Edges[j];
00878
00879 if(i_PresentVertex < m_vi_Edges[j])
00880 {
00881 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
00882 }
00883 else
00884 {
00885 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
00886 }
00887 }
00888
00889 break;
00890 }
00891 }
00892
00893 if (!_FOUND)
00894 {
00895 i_VertexOne = vi_FirstSeenOne[m_vi_VertexColors[m_vi_Edges[j]]];
00896 i_VertexTwo = vi_FirstSeenTwo[m_vi_VertexColors[m_vi_Edges[j]]];
00897
00898 if((i_VertexOne == i_PresentVertex) && (i_VertexTwo != m_vi_Edges[j]))
00899 {
00900 if(i_PresentVertex < i_VertexTwo)
00901 {
00902 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][i_VertexTwo]];
00903
00904 vi_StarHubMap[i_StarID] = i_PresentVertex;
00905
00906 if(i_PresentVertex < m_vi_Edges[j])
00907 {
00908 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
00909 }
00910 else
00911 {
00912 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
00913 }
00914 }
00915 else
00916 {
00917 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[i_VertexTwo][i_PresentVertex]];
00918
00919 vi_StarHubMap[i_StarID] = i_PresentVertex;
00920
00921 if(i_PresentVertex < m_vi_Edges[j])
00922 {
00923 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
00924 }
00925 else
00926 {
00927 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
00928 }
00929 }
00930 }
00931 }
00932 }
00933 }
00934
00935 #if VERBOSE == _TRUE
00936
00937 cout<<endl;
00938
00939 #endif
00940
00941 #if STATISTICS == _TRUE
00942
00943
00944
00945
00946
00947
00948
00949
00950
00951
00952
00953
00954
00955
00956
00957
00958
00959
00960
00961
00962
00963
00964
00965
00966
00967
00968
00969
00970 #endif
00971
00972 return(_TRUE);
00973
00974 }
00975
00976
00977
00978 int GraphColoring::StarColoring(vector<int> & vi_StarHubMap, vector<int> & vi_EdgeStarMap, map< int, map<int, int> > & mimi2_VertexEdgeMap)
00979 {
00980 int i, j, k;
00981
00982 int _FOUND;
00983
00984 int i_ColorID, i_StarID;
00985
00986 int i_PresentVertex;
00987
00988 int i_VertexCount, i_EdgeCount;
00989
00990 int i_VertexOne, i_VertexTwo;
00991
00992 vector<int> vi_MemberEdges;
00993
00994 vector<int> vi_CandidateColors;
00995
00996 vector<int> vi_FirstTreated;
00997
00998 vector<int> vi_FirstSeenOne, vi_FirstSeenTwo;
00999
01000 m_i_VertexColorCount = _UNKNOWN;
01001
01002 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
01003
01004 i_EdgeCount = (signed) m_vi_Edges.size();
01005
01006 vi_EdgeStarMap.clear();
01007 vi_EdgeStarMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
01008
01009 vi_StarHubMap.clear();
01010 vi_StarHubMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
01011
01012 m_vi_VertexColors.clear();
01013 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
01014
01015 vi_CandidateColors.clear();
01016 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
01017
01018 vi_FirstSeenOne.clear();
01019 vi_FirstSeenOne.resize((unsigned) i_VertexCount, _UNKNOWN);
01020
01021 vi_FirstSeenTwo.clear();
01022 vi_FirstSeenTwo.resize((unsigned) i_VertexCount, _UNKNOWN);
01023
01024
01025
01026
01027 vi_FirstTreated.clear();
01028 vi_FirstTreated.resize((unsigned) i_VertexCount, _UNKNOWN);
01029
01030 k = _FALSE;
01031
01032 for(i=0; i<i_VertexCount; i++)
01033 {
01034 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
01035 {
01036 if(i < m_vi_Edges[j])
01037 {
01038 mimi2_VertexEdgeMap[i][m_vi_Edges[j]] = k;
01039
01040 vi_EdgeStarMap[k] = k;
01041
01042 k++;
01043 }
01044 }
01045 }
01046
01047 #if VERBOSE == _TRUE
01048
01049 cout<<endl;
01050
01051 #endif
01052
01053 for(i=0; i<i_VertexCount; i++)
01054 {
01055 i_PresentVertex = m_vi_OrderedVertices[i];
01056
01057 #if VERBOSE == _TRUE
01058
01059 cout<<"DEBUG 305 | Star Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
01060
01061 #endif
01062 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
01063 {
01064 i_ColorID = m_vi_VertexColors[m_vi_Edges[j]];
01065
01066 if(i_ColorID == _UNKNOWN)
01067 {
01068 continue;
01069 }
01070
01071 vi_CandidateColors[i_ColorID] = i_PresentVertex;
01072
01073 i_VertexOne = vi_FirstSeenOne[i_ColorID];
01074 i_VertexTwo = vi_FirstSeenTwo[i_ColorID];
01075
01076 if(i_VertexOne == i_PresentVertex)
01077 {
01078 if(vi_FirstTreated[i_VertexTwo] != i_PresentVertex)
01079 {
01080
01081 for(k=m_vi_Vertices[i_VertexTwo]; k<m_vi_Vertices[STEP_UP(i_VertexTwo)]; k++)
01082 {
01083 if(m_vi_Edges[k] == i_PresentVertex)
01084 {
01085 continue;
01086 }
01087
01088 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
01089 {
01090 continue;
01091 }
01092
01093 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
01094
01095 }
01096
01097 vi_FirstTreated[i_VertexTwo] = i_PresentVertex;
01098 }
01099
01100 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
01101 {
01102 if(m_vi_Edges[k] == i_PresentVertex)
01103 {
01104 continue;
01105 }
01106
01107 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
01108 {
01109 continue;
01110 }
01111
01112 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
01113
01114 }
01115
01116 vi_FirstTreated[m_vi_Edges[j]] = i_PresentVertex;
01117 }
01118 else
01119 {
01120 vi_FirstSeenOne[i_ColorID] = i_PresentVertex;
01121 vi_FirstSeenTwo[i_ColorID] = m_vi_Edges[j];
01122
01123 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
01124 {
01125 if(m_vi_Edges[k] == i_PresentVertex)
01126 {
01127 continue;
01128 }
01129
01130 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
01131 {
01132 continue;
01133 }
01134
01135 if(m_vi_Edges[j] < m_vi_Edges[k])
01136 {
01137 if(vi_StarHubMap[vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]]] == m_vi_Edges[k])
01138 {
01139 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
01140 }
01141 }
01142 else
01143 {
01144 if(vi_StarHubMap[vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]]] == m_vi_Edges[k])
01145 {
01146 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
01147 }
01148 }
01149 }
01150 }
01151 }
01152
01153 for(j=0; j<i_VertexCount; j++)
01154 {
01155 if(vi_CandidateColors[j] != i_PresentVertex)
01156 {
01157 m_vi_VertexColors[i_PresentVertex] = j;
01158
01159 if(m_i_VertexColorCount < j)
01160 {
01161 m_i_VertexColorCount = j;
01162 }
01163
01164 break;
01165 }
01166 }
01167
01168 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
01169 {
01170 _FOUND = _FALSE;
01171
01172 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
01173 {
01174 continue;
01175 }
01176
01177 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
01178 {
01179 if(m_vi_Edges[k] == i_PresentVertex)
01180 {
01181 continue;
01182 }
01183
01184 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
01185 {
01186 continue;
01187 }
01188
01189 if(m_vi_VertexColors[m_vi_Edges[k]] == m_vi_VertexColors[i_PresentVertex])
01190 {
01191 _FOUND = _TRUE;
01192
01193 if(m_vi_Edges[j] < m_vi_Edges[k])
01194 {
01195 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]];
01196
01197 vi_StarHubMap[i_StarID] = m_vi_Edges[j];
01198
01199 if(i_PresentVertex < m_vi_Edges[j])
01200 {
01201 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
01202 }
01203 else
01204 {
01205 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
01206 }
01207 }
01208 else
01209 {
01210 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]];
01211
01212 vi_StarHubMap[i_StarID] = m_vi_Edges[j];
01213
01214 if(i_PresentVertex < m_vi_Edges[j])
01215 {
01216 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
01217 }
01218 else
01219 {
01220 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
01221 }
01222 }
01223
01224 break;
01225 }
01226 }
01227
01228 if (!_FOUND)
01229 {
01230 i_VertexOne = vi_FirstSeenOne[m_vi_VertexColors[m_vi_Edges[j]]];
01231 i_VertexTwo = vi_FirstSeenTwo[m_vi_VertexColors[m_vi_Edges[j]]];
01232
01233 if((i_VertexOne == i_PresentVertex) && (i_VertexTwo != m_vi_Edges[j]))
01234 {
01235 if(i_PresentVertex < i_VertexTwo)
01236 {
01237 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][i_VertexTwo]];
01238
01239 vi_StarHubMap[i_StarID] = i_PresentVertex;
01240
01241 if(i_PresentVertex < m_vi_Edges[j])
01242 {
01243 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
01244 }
01245 else
01246 {
01247 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
01248 }
01249 }
01250 else
01251 {
01252 i_StarID = vi_EdgeStarMap[mimi2_VertexEdgeMap[i_VertexTwo][i_PresentVertex]];
01253
01254 vi_StarHubMap[i_StarID] = i_PresentVertex;
01255
01256 if(i_PresentVertex < m_vi_Edges[j])
01257 {
01258 vi_EdgeStarMap[mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]]] = i_StarID;
01259 }
01260 else
01261 {
01262 vi_EdgeStarMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex]] = i_StarID;
01263 }
01264 }
01265 }
01266 }
01267 }
01268 }
01269
01270 #if VERBOSE == _TRUE
01271
01272 cout<<endl;
01273
01274 #endif
01275
01276 #if STATISTICS == _TRUE
01277
01278 vector<int> vi_Hubs;
01279
01280 vi_Hubs.resize((unsigned) i_EdgeCount/2, _FALSE);
01281
01282 m_i_ColoringUnits = _FALSE;
01283
01284 for(i=0; i<i_EdgeCount/2; i++)
01285 {
01286 if(vi_StarHubMap[i] == _UNKNOWN)
01287 {
01288 m_i_ColoringUnits++;
01289
01290 continue;
01291 }
01292
01293 if(vi_Hubs[vi_StarHubMap[i]] == _FALSE)
01294 {
01295 vi_Hubs[vi_StarHubMap[i]] = _TRUE;
01296
01297 m_i_ColoringUnits++;
01298 }
01299 }
01300
01301 #endif
01302
01303 return(_TRUE);
01304
01305 }
01306
01307
01308
01309
01310
01311 int GraphColoring::CheckStarColoring()
01312 {
01313 int i, j, k, l;
01314
01315 int i_VertexCount, i_EdgeCount;
01316
01317 int i_FirstColor, i_SecondColor, i_ThirdColor, i_FourthColor;
01318
01319 int i_ViolationCount;
01320
01321 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
01322
01323 i_EdgeCount = (signed) m_vi_Edges.size();
01324
01325 i_ViolationCount = _FALSE;
01326
01327 for(i=0; i<i_VertexCount; i++)
01328 {
01329 i_FirstColor = m_vi_VertexColors[i];
01330
01331 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
01332 {
01333 i_SecondColor = m_vi_VertexColors[m_vi_Edges[j]];
01334
01335 if(i_SecondColor == i_FirstColor)
01336 {
01337 i_ViolationCount++;
01338
01339 if(i_ViolationCount == _TRUE)
01340 {
01341 cout<<endl;
01342 cout<<"Star Coloring | Violation Check | "<<m_s_InputFile<<endl;
01343 cout<<endl;
01344 }
01345
01346 cout<<"Violation "<<i_ViolationCount<<"\t : "<<STEP_UP(i)<<" ["<<STEP_UP(i_FirstColor)<<"] - "<<STEP_UP(m_vi_Edges[j])<<" ["<<STEP_UP(i_SecondColor)<<"]"<<endl;
01347
01348 continue;
01349 }
01350
01351 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
01352 {
01353 if(m_vi_Edges[k] == i)
01354 {
01355 continue;
01356 }
01357
01358 i_ThirdColor = m_vi_VertexColors[m_vi_Edges[k]];
01359
01360 if(i_ThirdColor == i_SecondColor)
01361 {
01362 i_ViolationCount++;
01363
01364 if(i_ViolationCount == _TRUE)
01365 {
01366 cout<<endl;
01367 cout<<"Star Coloring | Violation Check | "<<m_s_InputFile<<endl;
01368 cout<<endl;
01369 }
01370
01371 cout<<"Violation "<<i_ViolationCount<<"\t : "<<STEP_UP(m_vi_Edges[j])<<" ["<<STEP_UP(i_SecondColor)<<"] - "<<STEP_UP(m_vi_Edges[k])<<" ["<<STEP_UP(i_ThirdColor)<<"]"<<endl;
01372
01373 continue;
01374 }
01375
01376 if(i_ThirdColor != i_FirstColor)
01377 {
01378 continue;
01379 }
01380
01381 if(i_ThirdColor == i_FirstColor)
01382 {
01383 for(l=m_vi_Vertices[m_vi_Edges[k]]; l<m_vi_Vertices[STEP_UP(m_vi_Edges[k])]; l++)
01384 {
01385 if((m_vi_Edges[l] == m_vi_Edges[j]))
01386 {
01387 continue;
01388 }
01389
01390 i_FourthColor = m_vi_VertexColors[m_vi_Edges[l]];
01391
01392 if(i_FourthColor == i_ThirdColor)
01393 {
01394 i_ViolationCount++;
01395
01396 if(i_ViolationCount == _TRUE)
01397 {
01398 cout<<endl;
01399 cout<<"Star Coloring | Violation Check | "<<m_s_InputFile<<endl;
01400 cout<<endl;
01401 }
01402
01403 cout<<"Violation "<<i_ViolationCount<<"\t : "<<STEP_UP(m_vi_Edges[k])<<" ["<<STEP_UP(i_ThirdColor)<<"] - "<<STEP_UP(m_vi_Edges[l])<<" ["<<STEP_UP(i_FourthColor)<<"]"<<endl;
01404
01405 }
01406
01407 if(i_FourthColor == i_SecondColor)
01408 {
01409 i_ViolationCount++;
01410
01411 if(i_ViolationCount == _TRUE)
01412 {
01413 cout<<endl;
01414 cout<<"Star Coloring | Violation Check | "<<m_s_InputFile<<endl;
01415 cout<<endl;
01416 }
01417
01418 cout<<"Violation "<<i_ViolationCount<<"\t : "<<STEP_UP(i)<<" ["<<STEP_UP(i_FirstColor)<<"] - "<<STEP_UP(m_vi_Edges[j])<<" ["<<STEP_UP(i_SecondColor)<<"] - "<<STEP_UP(m_vi_Edges[k])<<" ["<<STEP_UP(i_ThirdColor)<<"] - "<<STEP_UP(m_vi_Edges[l])<<" ["<<STEP_UP(i_FourthColor)<<"]"<<endl;
01419
01420 continue;
01421 }
01422 }
01423 }
01424 }
01425 }
01426 }
01427
01428 if(i_ViolationCount)
01429 {
01430 cout<<endl;
01431 cout<<"[Total Violations = "<<i_ViolationCount<<"]"<<endl;
01432 cout<<endl;
01433 }
01434
01435 return(i_ViolationCount);
01436 }
01437
01438
01439
01440
01441 int GraphColoring::AcyclicColoring()
01442 {
01443 if(CheckVertexColoring("ACYCLIC"))
01444 {
01445 return(_TRUE);
01446 }
01447
01448 int i, j, k;
01449
01450 int i_VertexCount, i_EdgeCount;
01451
01452 int i_AdjacentEdgeID, i_EdgeID, i_SetID;
01453
01454 int i_PresentVertex;
01455
01456 vector<int> vi_CandidateColors;
01457
01458 vector<int> vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree;
01459 vector<int> vi_FirstVisitedOne, vi_FirstVisitedTwo;
01460
01461 #if DISJOINT_SETS == _FALSE
01462
01463 int l;
01464
01465 int i_MemberCount;
01466
01467 int i_SmallerSetID, i_BiggerSetID;
01468
01469 vector<int> vi_DisjointSets;
01470 vector<int> vi_MemberEdges;
01471 vector<int> vi_EdgeSetMap;
01472
01473 vector< vector<int> > v2i_SetEdgeMap;
01474
01475 #endif
01476
01477 #if DISJOINT_SETS == _TRUE
01478
01479 int i_SetOneID, i_SetTwoID;
01480
01481 #endif
01482
01483 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
01484
01485 k = _FALSE;
01486
01487 m_mimi2_VertexEdgeMap.clear();
01488
01489 for(i=0; i<i_VertexCount; i++)
01490 {
01491 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
01492 {
01493 if(i < m_vi_Edges[j])
01494 {
01495 m_mimi2_VertexEdgeMap[i][m_vi_Edges[j]] = k;
01496
01497 k++;
01498 }
01499 }
01500 }
01501
01502 #if DEBUG == 1461
01503
01504 i_EdgeCount = (signed) v2i_EdgeVertexMap.size();
01505
01506 cout<<endl;
01507 cout<<"DEBUG 1461 | Acyclic Coloring | Edge Vertex Map"<<endl;
01508 cout<<endl;
01509
01510 for(i=0; i<i_EdgeCount; i++)
01511 {
01512 cout<<"Edge "<<STEP_UP(i)<<"\t"<<" : "<<STEP_UP(v2i_EdgeVertexMap[i][0])<<" - "<<STEP_UP(v2i_EdgeVertexMap[i][1])<<endl;
01513 }
01514
01515 cout<<endl;
01516 cout<<"DEBUG 1461 | Acyclic Coloring | Vertex Edge Map"<<endl;
01517 cout<<endl;
01518
01519 for(i=0; i<i_EdgeCount; i++)
01520 {
01521 cout<<"Vertex "<<STEP_UP(v2i_EdgeVertexMap[i][0])<<" - Vertex "<<STEP_UP(v2i_EdgeVertexMap[i][1])<<"\t"<<" : "<<STEP_UP(m_mimi2_VertexEdgeMap[v2i_EdgeVertexMap[i][0]][v2i_EdgeVertexMap[i][1]])<<endl;
01522
01523 }
01524
01525 cout<<endl;
01526
01527 #endif
01528
01529 i_EdgeCount = (signed) m_vi_Edges.size();
01530
01531 m_vi_VertexColors.clear();
01532 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
01533
01534 vi_CandidateColors.clear();
01535 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
01536
01537 vi_FirstSeenOne.clear();
01538 vi_FirstSeenOne.resize((unsigned) i_VertexCount, _UNKNOWN);
01539
01540 vi_FirstSeenTwo.clear();
01541 vi_FirstSeenTwo.resize((unsigned) i_VertexCount, _UNKNOWN);
01542
01543 vi_FirstSeenThree.clear();
01544 vi_FirstSeenThree.resize((unsigned) i_VertexCount, _UNKNOWN);
01545
01546 vi_FirstVisitedOne.clear();
01547 vi_FirstVisitedOne.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
01548
01549 vi_FirstVisitedTwo.clear();
01550 vi_FirstVisitedTwo.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
01551
01552 #if DISJOINT_SETS == _FALSE
01553
01554 vi_MemberEdges.clear();
01555
01556 vi_EdgeSetMap.clear();
01557 vi_EdgeSetMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
01558
01559 v2i_SetEdgeMap.clear();
01560 v2i_SetEdgeMap.resize((unsigned) i_EdgeCount/2);
01561
01562 vi_DisjointSets.clear();
01563
01564 #endif
01565
01566 #if DISJOINT_SETS == _TRUE
01567
01568 m_ds_DisjointSets.SetSize(i_EdgeCount/2);
01569
01570 #endif
01571
01572 #if VERBOSE == _TRUE
01573
01574 cout<<endl;
01575
01576 #endif
01577
01578 m_i_VertexColorCount = _UNKNOWN;
01579
01580 for(i=0; i<i_VertexCount; i++)
01581 {
01582 i_PresentVertex = m_vi_OrderedVertices[i];
01583
01584 #if VERBOSE == _TRUE
01585
01586 cout<<"DEBUG 1461 | Acyclic Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
01587
01588 #endif
01589
01590 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
01591 {
01592 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
01593 {
01594 continue;
01595 }
01596
01597 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
01598 }
01599
01600 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
01601 {
01602 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
01603 {
01604 continue;
01605 }
01606
01607 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
01608 {
01609 if(m_vi_Edges[k] == i_PresentVertex)
01610 {
01611 continue;
01612 }
01613
01614 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
01615 {
01616 continue;
01617 }
01618
01619 if(vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] == i_PresentVertex)
01620 {
01621 continue;
01622 }
01623
01624 #if DISJOINT_SETS == _TRUE
01625
01626 if(m_vi_Edges[j] < m_vi_Edges[k])
01627 {
01628 i_SetID = m_ds_DisjointSets.FindAndCompress(m_mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]);
01629 }
01630 else
01631 {
01632 i_SetID = m_ds_DisjointSets.FindAndCompress(m_mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]);
01633 }
01634 #endif
01635
01636 #if DISJOINT_SETS == _FALSE
01637
01638 if(m_vi_Edges[j] < m_vi_Edges[k])
01639 {
01640 i_SetID = vi_EdgeSetMap[m_mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]];
01641 }
01642 else
01643 {
01644 i_SetID = vi_EdgeSetMap[m_mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]];
01645 }
01646 #endif
01647
01648 FindCycle(i_PresentVertex, m_vi_Edges[j], m_vi_Edges[k], i_SetID, vi_CandidateColors, vi_FirstVisitedOne, vi_FirstVisitedTwo);
01649 }
01650 }
01651
01652 for(j=0; j<i_VertexCount; j++)
01653 {
01654 if(vi_CandidateColors[j] != i_PresentVertex)
01655 {
01656 m_vi_VertexColors[i_PresentVertex] = j;
01657
01658 if(m_i_VertexColorCount < j)
01659 {
01660 m_i_VertexColorCount = j;
01661 }
01662
01663 break;
01664 }
01665 }
01666
01667 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
01668 {
01669 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
01670 {
01671 continue;
01672 }
01673
01674 if(i_PresentVertex < m_vi_Edges[j])
01675 {
01676 i_EdgeID = m_mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]];
01677 }
01678 else
01679 {
01680 i_EdgeID = m_mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex];
01681 }
01682
01683 #if DISJOINT_SETS == _FALSE
01684
01685 vi_DisjointSets.push_back(_TRUE);
01686
01687 vi_EdgeSetMap[i_EdgeID] = STEP_DOWN((signed) vi_DisjointSets.size());
01688
01689 v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].push_back(i_EdgeID);
01690 #endif
01691
01692 i_AdjacentEdgeID = UpdateSet(i_PresentVertex, m_vi_Edges[j], i_PresentVertex, m_mimi2_VertexEdgeMap, vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree);
01693
01694 if(i_AdjacentEdgeID != _UNKNOWN)
01695 {
01696
01697 #if DISJOINT_SETS == _TRUE
01698
01699 i_SetOneID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
01700 i_SetTwoID = m_ds_DisjointSets.FindAndCompress(i_AdjacentEdgeID);
01701
01702 #if DEBUG == 1461
01703
01704 cout<<endl;
01705 cout<<"DEBUG 1461 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(i_SetOneID)<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(i_SetTwoID)<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
01706 cout<<endl;
01707
01708 #endif
01709
01710 m_ds_DisjointSets.UnionBySize(i_SetOneID, i_SetTwoID);
01711
01712 #endif
01713
01714 #if DISJOINT_SETS == _FALSE
01715
01716 #if DEBUG == 1461
01717
01718 cout<<endl;
01719 cout<<"DEBUG 1461 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(vi_EdgeSetMap[i_EdgeID])<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(vi_EdgeSetMap[i_AdjacentEdgeID])<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
01720 cout<<endl;
01721
01722 #endif
01723
01724 if(v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].size() < v2i_SetEdgeMap[vi_EdgeSetMap[i_AdjacentEdgeID]].size())
01725 {
01726 i_SmallerSetID = vi_EdgeSetMap[i_EdgeID];
01727 i_BiggerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
01728 }
01729 else
01730 {
01731 i_BiggerSetID = vi_EdgeSetMap[i_EdgeID];
01732 i_SmallerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
01733 }
01734
01735 vi_MemberEdges.clear();
01736 vi_MemberEdges.swap(v2i_SetEdgeMap[i_BiggerSetID]);
01737
01738 vi_DisjointSets[i_BiggerSetID] = _FALSE;
01739
01740 i_MemberCount = (signed) vi_MemberEdges.size();
01741
01742 for(k=0; k<i_MemberCount; k++)
01743 {
01744 vi_EdgeSetMap[vi_MemberEdges[k]] = i_SmallerSetID;
01745
01746 v2i_SetEdgeMap[i_SmallerSetID].push_back(vi_MemberEdges[k]);
01747 }
01748 #endif
01749 }
01750 }
01751
01752
01753 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
01754 {
01755 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
01756 {
01757 continue;
01758 }
01759
01760 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
01761 {
01762 if(m_vi_Edges[k] == i_PresentVertex)
01763 {
01764 continue;
01765 }
01766
01767 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
01768 {
01769 continue;
01770 }
01771
01772 if(m_vi_VertexColors[m_vi_Edges[k]] == m_vi_VertexColors[i_PresentVertex])
01773 {
01774 if(m_vi_Edges[j] < m_vi_Edges[k])
01775 {
01776 i_AdjacentEdgeID = m_mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]];
01777 }
01778 else
01779 {
01780 i_AdjacentEdgeID = m_mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]];
01781 }
01782
01783 i_EdgeID = UpdateSet(i_PresentVertex, m_vi_Edges[j], m_vi_Edges[k], m_mimi2_VertexEdgeMap, vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree);
01784
01785 if(i_EdgeID != _UNKNOWN)
01786 {
01787
01788 #if DISJOINT_SETS == _TRUE
01789
01790 i_SetOneID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
01791 i_SetTwoID = m_ds_DisjointSets.FindAndCompress(i_AdjacentEdgeID);
01792
01793 #if DEBUG == 1461
01794 cout<<endl;
01795 cout<<"DEBUG 1461 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(i_SetOneID)<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(i_SetTwoID)<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
01796 cout<<endl;
01797
01798 #endif
01799
01800 m_ds_DisjointSets.UnionBySize(i_SetOneID, i_SetTwoID);
01801
01802 #endif
01803
01804 #if DISJOINT_SETS == _FALSE
01805
01806 #if DEBUG == 1461
01807 cout<<endl;
01808 cout<<"DEBUG 1461 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(vi_EdgeSetMap[i_EdgeID])<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(vi_EdgeSetMap[i_AdjacentEdgeID])<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
01809 cout<<endl;
01810
01811 #endif
01812
01813 if(v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].size() < v2i_SetEdgeMap[vi_EdgeSetMap[i_AdjacentEdgeID]].size())
01814 {
01815 i_SmallerSetID = vi_EdgeSetMap[i_EdgeID];
01816 i_BiggerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
01817 }
01818 else
01819 {
01820 i_BiggerSetID = vi_EdgeSetMap[i_EdgeID];
01821 i_SmallerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
01822 }
01823
01824 vi_MemberEdges.clear();
01825 vi_MemberEdges.swap(v2i_SetEdgeMap[i_BiggerSetID]);
01826
01827 vi_DisjointSets[i_BiggerSetID] = _FALSE;
01828
01829 i_MemberCount = (signed) vi_MemberEdges.size();
01830
01831 for(l=0; l<i_MemberCount; l++)
01832 {
01833 vi_EdgeSetMap[vi_MemberEdges[l]] = i_SmallerSetID;
01834
01835 v2i_SetEdgeMap[i_SmallerSetID].push_back(vi_MemberEdges[l]);
01836 }
01837
01838 #endif
01839 }
01840 }
01841 }
01842 }
01843
01844 #if DEBUG == 1461
01845
01846 cout<<endl;
01847 cout<<"DEBUG 1461 | Acyclic Coloring | Vertex Colors | Upto Vertex "<<STEP_UP(i)<<endl;
01848 cout<<endl;
01849
01850 for(j=0; j<i_VertexCount; j++)
01851 {
01852 cout<<"Vertex "<<STEP_UP(j)<<" = "<<STEP_UP(m_vi_VertexColors[j])<<endl;
01853 }
01854 #endif
01855
01856 }
01857
01858
01859 #if DEBUG == 1461
01860
01861 #if DISJOINT_SETS == _FALSE
01862
01863 i_EdgeCount = (signed) v2i_EdgeVertexMap.size();
01864
01865 cout<<endl;
01866 cout<<"DEBUG 1461 | Acyclic Coloring | Edge Set Map"<<endl;
01867 cout<<endl;
01868
01869 for(i=0; i<i_EdgeCount; i++)
01870 {
01871 cout<<STEP_UP(i)<<"\t"<<" : "<<STEP_UP(vi_EdgeSetMap[i])<<endl;
01872 }
01873
01874 cout<<endl;
01875 cout<<"DEBUG 1461 | Acyclic Coloring | Set Edge Map"<<endl;
01876 cout<<endl;
01877
01878 for(i=0; i<i_EdgeCount; i++)
01879 {
01880 i_MemberCount = (signed) v2i_SetEdgeMap[i].size();
01881
01882 if(i_MemberCount == _FALSE)
01883 {
01884 continue;
01885 }
01886
01887 cout<<"Set "<<STEP_UP(i)<<"\t"<<" : ";
01888
01889 for(j=0; j<i_MemberCount; j++)
01890 {
01891 if(j == STEP_DOWN(i_MemberCount))
01892 {
01893 cout<<STEP_UP(v2i_SetEdgeMap[i][j])<<" ("<<i_MemberCount<<")"<<endl;
01894 }
01895 else
01896 {
01897 cout<<STEP_UP(v2i_SetEdgeMap[i][j])<<", ";
01898 }
01899 }
01900 }
01901
01902 cout<<endl;
01903
01904 #endif
01905
01906 #endif
01907
01908 #if VERBOSE == _TRUE
01909
01910 cout<<endl;
01911
01912 #endif
01913
01914 #if STATISTICS == _TRUE
01915
01916 #if DISJOINT_SETS == _TRUE
01917
01918 m_i_ColoringUnits = m_ds_DisjointSets.Count();
01919
01920
01921 #elif DISJOINT_SETS == _FALSE
01922
01923 int i_SetSize;
01924
01925 m_i_ColoringUnits = _FALSE;
01926
01927 i_SetSize = (unsigned) v2i_SetEdgeMap.size();
01928
01929 for(i=0; i<i_SetSize; i++)
01930 {
01931 if(v2i_SetEdgeMap[i].empty())
01932 {
01933 continue;
01934 }
01935
01936 m_i_ColoringUnits++;
01937 }
01938
01939 #endif
01940
01941 #endif
01942
01943
01944
01945
01946
01947
01948 return(_TRUE);
01949
01950 }
01951
01952 int GraphColoring::AcyclicColoring_ForIndirectRecovery() {
01953
01954 if(CheckVertexColoring("ACYCLIC"))
01955 {
01956 return(_TRUE);
01957 }
01958
01959 int i, j, k;
01960
01961 int i_VertexCount, i_EdgeCount;
01962
01963 int i_AdjacentEdgeID, i_EdgeID, i_SetID;
01964
01965 int i_PresentVertex;
01966
01967 vector<int> vi_CandidateColors;
01968
01969 vector<int> vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree;
01970 vector<int> vi_FirstVisitedOne, vi_FirstVisitedTwo;
01971
01972
01973
01974 #if DISJOINT_SETS == _FALSE
01975
01976 int l;
01977
01978 int i_MemberCount;
01979
01980 int i_SmallerSetID, i_BiggerSetID;
01981
01982 vector<int> vi_DisjointSets;
01983 vector<int> vi_MemberEdges;
01984 vector<int> vi_EdgeSetMap;
01985
01986 vector< vector<int> > v2i_SetEdgeMap;
01987
01988 #endif
01989
01990 #if DISJOINT_SETS == _TRUE
01991
01992 int i_SetOneID, i_SetTwoID;
01993
01994
01995
01996 #endif
01997
01998 if(m_s_VertexColoringVariant.compare("ALL") != 0)
01999 {
02000 m_s_VertexColoringVariant = "ACYCLIC";
02001 }
02002
02003 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
02004
02005 k=_FALSE;
02006
02007
02008
02009 m_mimi2_VertexEdgeMap.clear();
02010
02011 for(i=0; i<i_VertexCount; i++)
02012 {
02013 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
02014 {
02015 if(i < m_vi_Edges[j])
02016 {
02017 m_mimi2_VertexEdgeMap[i][m_vi_Edges[j]] = k;
02018
02019 k++;
02020 }
02021 }
02022 }
02023
02024
02025 #if DEBUG == 1462
02026
02027 cout<<endl;
02028 cout<<"DEBUG 1462 | Acyclic Coloring | Edge Vertex Map"<<endl;
02029 cout<<endl;
02030
02031 for(i=0; i<i_VertexCount; i++)
02032 {
02033 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
02034 {
02035 if(i < m_vi_Edges[j])
02036 {
02037 cout<<"Edge "<<STEP_UP(m_mimi2_VertexEdgeMap[i][m_vi_Edges[j]])<<"\t"<<" : "<<STEP_UP(i)<<" - "<<STEP_UP(m_vi_Edges[j])<<endl;
02038 }
02039 }
02040 }
02041
02042 cout<<endl;
02043
02044 #endif
02045
02046 i_EdgeCount = (signed) m_vi_Edges.size();
02047
02048 m_vi_VertexColors.clear();
02049 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
02050
02051 vi_CandidateColors.clear();
02052 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
02053
02054 vi_FirstSeenOne.clear();
02055 vi_FirstSeenOne.resize((unsigned) i_VertexCount, _UNKNOWN);
02056
02057 vi_FirstSeenTwo.clear();
02058 vi_FirstSeenTwo.resize((unsigned) i_VertexCount, _UNKNOWN);
02059
02060 vi_FirstSeenThree.clear();
02061 vi_FirstSeenThree.resize((unsigned) i_VertexCount, _UNKNOWN);
02062
02063 vi_FirstVisitedOne.clear();
02064 vi_FirstVisitedOne.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
02065
02066 vi_FirstVisitedTwo.clear();
02067 vi_FirstVisitedTwo.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
02068
02069
02070 #if DISJOINT_SETS == _FALSE
02071
02072 vi_MemberEdges.clear();
02073
02074 vi_EdgeSetMap.clear();
02075 vi_EdgeSetMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
02076
02077 v2i_SetEdgeMap.clear();
02078 v2i_SetEdgeMap.resize((unsigned) i_EdgeCount/2);
02079
02080 vi_DisjointSets.clear();
02081
02082 #endif
02083
02084 #if DISJOINT_SETS == _TRUE
02085
02086 m_ds_DisjointSets.SetSize(i_EdgeCount/2);
02087
02088 #endif
02089
02090 #if VERBOSE == _TRUE
02091
02092 cout<<endl;
02093
02094 #endif
02095
02096 m_i_VertexColorCount = _UNKNOWN;
02097
02098
02099 for(i=0; i<i_VertexCount; i++)
02100 {
02101
02102 i_PresentVertex = m_vi_OrderedVertices[i];
02103
02104
02105 #if VERBOSE == _TRUE
02106
02107
02108 cout<<"DEBUG 1462 | Acyclic Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
02109
02110 #endif
02111
02112 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02113 {
02114 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02115 {
02116 continue;
02117 }
02118
02119 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
02120 }
02121
02122 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02123 {
02124 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02125 {
02126 continue;
02127 }
02128
02129 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
02130 {
02131 if(m_vi_Edges[k] == i_PresentVertex)
02132 {
02133 continue;
02134 }
02135
02136 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
02137 {
02138 continue;
02139 }
02140
02141 if(vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] == i_PresentVertex)
02142 {
02143 continue;
02144 }
02145
02146 #if DISJOINT_SETS == _TRUE
02147
02148 if(m_vi_Edges[j] < m_vi_Edges[k])
02149 {
02150 i_SetID = m_ds_DisjointSets.FindAndCompress(m_mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]);
02151 }
02152 else
02153 {
02154 i_SetID = m_ds_DisjointSets.FindAndCompress(m_mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]);
02155 }
02156 #endif
02157
02158 #if DISJOINT_SETS == _FALSE
02159
02160 if(m_vi_Edges[j] < m_vi_Edges[k])
02161 {
02162 i_SetID = vi_EdgeSetMap[m_mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]];
02163 }
02164 else
02165 {
02166 i_SetID = vi_EdgeSetMap[m_mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]];
02167 }
02168 #endif
02169
02170 FindCycle(i_PresentVertex, m_vi_Edges[j], m_vi_Edges[k], i_SetID, vi_CandidateColors, vi_FirstVisitedOne, vi_FirstVisitedTwo);
02171 }
02172 }
02173
02174 for(j=0; j<i_VertexCount; j++)
02175 {
02176 if(vi_CandidateColors[j] != i_PresentVertex)
02177 {
02178 m_vi_VertexColors[i_PresentVertex] = j;
02179
02180 if(m_i_VertexColorCount < j)
02181 {
02182 m_i_VertexColorCount = j;
02183 }
02184
02185 break;
02186 }
02187 }
02188
02189 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02190 {
02191 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02192 {
02193 continue;
02194 }
02195
02196 if(i_PresentVertex < m_vi_Edges[j])
02197 {
02198 i_EdgeID = m_mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]];
02199 }
02200 else
02201 {
02202 i_EdgeID = m_mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex];
02203 }
02204
02205 #if DISJOINT_SETS == _FALSE
02206
02207 vi_DisjointSets.push_back(_TRUE);
02208
02209 vi_EdgeSetMap[i_EdgeID] = STEP_DOWN((signed) vi_DisjointSets.size());
02210
02211 v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].push_back(i_EdgeID);
02212 #endif
02213
02214
02215 i_AdjacentEdgeID = UpdateSet(i_PresentVertex, m_vi_Edges[j], i_PresentVertex, m_mimi2_VertexEdgeMap, vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree);
02216
02217 if(i_AdjacentEdgeID != _UNKNOWN)
02218 {
02219
02220 #if DISJOINT_SETS == _TRUE
02221
02222 i_SetOneID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
02223 i_SetTwoID = m_ds_DisjointSets.FindAndCompress(i_AdjacentEdgeID);
02224
02225 #if DEBUG == 1462
02226
02227 cout<<endl;
02228 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(i_SetOneID)<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(i_SetTwoID)<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02229 cout<<endl;
02230
02231 #endif
02232
02233
02234
02235
02236
02237 m_ds_DisjointSets.UnionBySize(i_SetOneID, i_SetTwoID);
02238
02239
02240
02241
02242
02243 #endif
02244
02245 #if DISJOINT_SETS == _FALSE
02246
02247 #if DEBUG == 1462
02248
02249 cout<<endl;
02250 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(vi_EdgeSetMap[i_EdgeID])<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(vi_EdgeSetMap[i_AdjacentEdgeID])<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02251 cout<<endl;
02252
02253 #endif
02254
02255 if(v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].size() < v2i_SetEdgeMap[vi_EdgeSetMap[i_AdjacentEdgeID]].size())
02256 {
02257 i_SmallerSetID = vi_EdgeSetMap[i_EdgeID];
02258 i_BiggerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02259 }
02260 else
02261 {
02262 i_BiggerSetID = vi_EdgeSetMap[i_EdgeID];
02263 i_SmallerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02264 }
02265
02266 vi_MemberEdges.clear();
02267 vi_MemberEdges.swap(v2i_SetEdgeMap[i_BiggerSetID]);
02268
02269 vi_DisjointSets[i_BiggerSetID] = _FALSE;
02270
02271 i_MemberCount = (signed) vi_MemberEdges.size();
02272
02273 for(k=0; k<i_MemberCount; k++)
02274 {
02275 vi_EdgeSetMap[vi_MemberEdges[k]] = i_SmallerSetID;
02276
02277 v2i_SetEdgeMap[i_SmallerSetID].push_back(vi_MemberEdges[k]);
02278 }
02279 #endif
02280 }
02281 }
02282
02283
02284
02285 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02286 {
02287 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02288 {
02289 continue;
02290 }
02291
02292 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
02293 {
02294 if(m_vi_Edges[k] == i_PresentVertex)
02295 {
02296 continue;
02297 }
02298
02299 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
02300 {
02301 continue;
02302 }
02303
02304 if(m_vi_VertexColors[m_vi_Edges[k]] == m_vi_VertexColors[i_PresentVertex])
02305 {
02306 if(m_vi_Edges[j] < m_vi_Edges[k])
02307 {
02308 i_AdjacentEdgeID = m_mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]];
02309 }
02310 else
02311 {
02312 i_AdjacentEdgeID = m_mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]];
02313 }
02314
02315 i_EdgeID = UpdateSet(i_PresentVertex, m_vi_Edges[j], m_vi_Edges[k], m_mimi2_VertexEdgeMap, vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree);
02316
02317 if(i_EdgeID != _UNKNOWN)
02318 {
02319
02320 #if DISJOINT_SETS == _TRUE
02321
02322 i_SetOneID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
02323 i_SetTwoID = m_ds_DisjointSets.FindAndCompress(i_AdjacentEdgeID);
02324
02325 #if DEBUG == 1462
02326 cout<<endl;
02327 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(i_SetOneID)<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(i_SetTwoID)<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02328 cout<<endl;
02329
02330 #endif
02331
02332
02333
02334
02335
02336 m_ds_DisjointSets.UnionBySize(i_SetOneID, i_SetTwoID);
02337
02338
02339
02340
02341
02342 #endif
02343
02344 #if DISJOINT_SETS == _FALSE
02345
02346 #if DEBUG == 1462
02347
02348 cout<<endl;
02349 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(vi_EdgeSetMap[i_EdgeID])<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(vi_EdgeSetMap[i_AdjacentEdgeID])<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02350 cout<<endl;
02351
02352 #endif
02353
02354 if(v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].size() < v2i_SetEdgeMap[vi_EdgeSetMap[i_AdjacentEdgeID]].size())
02355 {
02356 i_SmallerSetID = vi_EdgeSetMap[i_EdgeID];
02357 i_BiggerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02358 }
02359 else
02360 {
02361 i_BiggerSetID = vi_EdgeSetMap[i_EdgeID];
02362 i_SmallerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02363 }
02364
02365 vi_MemberEdges.clear();
02366 vi_MemberEdges.swap(v2i_SetEdgeMap[i_BiggerSetID]);
02367
02368 vi_DisjointSets[i_BiggerSetID] = _FALSE;
02369
02370 i_MemberCount = (signed) vi_MemberEdges.size();
02371
02372 for(l=0; l<i_MemberCount; l++)
02373 {
02374 vi_EdgeSetMap[vi_MemberEdges[l]] = i_SmallerSetID;
02375
02376 v2i_SetEdgeMap[i_SmallerSetID].push_back(vi_MemberEdges[l]);
02377 }
02378
02379 #endif
02380 }
02381 }
02382 }
02383 }
02384
02385 #if DEBUG == 1462
02386
02387 cout<<endl;
02388 cout<<"DEBUG 1462 | Acyclic Coloring | Vertex Colors | Upto Vertex "<<STEP_UP(i_PresentVertex)<<endl;
02389 cout<<endl;
02390
02391 for(j=0; j<i_VertexCount; j++)
02392 {
02393 cout<<"Vertex "<<STEP_UP(j)<<" = "<<STEP_UP(m_vi_VertexColors[j])<<endl;
02394 }
02395 #endif
02396
02397 }
02398
02399
02400 #if DEBUG == 1462
02401
02402 #if DISJOINT_SETS == _FALSE
02403
02404 i_EdgeCount = (signed) v2i_EdgeVertexMap.size();
02405
02406 cout<<endl;
02407 cout<<"DEBUG 1462 | Acyclic Coloring | Edge Set Map"<<endl;
02408 cout<<endl;
02409
02410 for(i=0; i<i_EdgeCount; i++)
02411 {
02412 cout<<STEP_UP(i)<<"\t"<<" : "<<STEP_UP(vi_EdgeSetMap[i])<<endl;
02413 }
02414
02415 cout<<endl;
02416 cout<<"DEBUG 1462 | Acyclic Coloring | Set Edge Map"<<endl;
02417 cout<<endl;
02418
02419 for(i=0; i<i_EdgeCount; i++)
02420 {
02421 i_MemberCount = (signed) v2i_SetEdgeMap[i].size();
02422
02423 if(i_MemberCount == _FALSE)
02424 {
02425 continue;
02426 }
02427
02428 cout<<"Set "<<STEP_UP(i)<<"\t"<<" : ";
02429
02430 for(j=0; j<i_MemberCount; j++)
02431 {
02432 if(j == STEP_DOWN(i_MemberCount))
02433 {
02434 cout<<STEP_UP(v2i_SetEdgeMap[i][j])<<" ("<<i_MemberCount<<")"<<endl;
02435 }
02436 else
02437 {
02438 cout<<STEP_UP(v2i_SetEdgeMap[i][j])<<", ";
02439 }
02440 }
02441 }
02442
02443 cout<<endl;
02444
02445 #endif
02446
02447 #endif
02448
02449 #if VERBOSE == _TRUE
02450
02451 cout<<endl;
02452
02453 #endif
02454
02455
02456
02457
02458
02459
02460 return(_TRUE);
02461 }
02462
02463
02464 int GraphColoring::AcyclicColoring(vector<int> & vi_Sets, map< int, vector<int> > & mivi_VertexSets)
02465 {
02466
02467 if(CheckVertexColoring("ACYCLIC"))
02468 {
02469 return(_TRUE);
02470 }
02471
02472 int i, j, k;
02473
02474 int i_VertexCount, i_EdgeCount;
02475
02476 int i_AdjacentEdgeID, i_EdgeID, i_SetID;
02477
02478 int i_PresentVertex;
02479
02480 vector<int> vi_CandidateColors;
02481
02482 vector<int> vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree;
02483 vector<int> vi_FirstVisitedOne, vi_FirstVisitedTwo;
02484
02485 map< int, map<int, int> > mimi2_VertexEdgeMap;
02486
02487 #if DISJOINT_SETS == _FALSE
02488
02489 int l;
02490
02491 int i_MemberCount;
02492
02493 int i_SmallerSetID, i_BiggerSetID;
02494
02495 vector<int> vi_DisjointSets;
02496 vector<int> vi_MemberEdges;
02497 vector<int> vi_EdgeSetMap;
02498
02499 vector< vector<int> > v2i_SetEdgeMap;
02500
02501 #endif
02502
02503 #if DISJOINT_SETS == _TRUE
02504
02505 int i_SetOneID, i_SetTwoID;
02506
02507
02508
02509 #endif
02510
02511 if(m_s_VertexColoringVariant.compare("ALL") != 0)
02512 {
02513 m_s_VertexColoringVariant = "ACYCLIC";
02514 }
02515
02516 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
02517
02518 k=_FALSE;
02519
02520
02521
02522 mimi2_VertexEdgeMap.clear();
02523
02524 for(i=0; i<i_VertexCount; i++)
02525 {
02526 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
02527 {
02528 if(i < m_vi_Edges[j])
02529 {
02530 mimi2_VertexEdgeMap[i][m_vi_Edges[j]] = k;
02531
02532 k++;
02533 }
02534 }
02535 }
02536
02537 #if DEBUG == 1462
02538
02539 cout<<endl;
02540 cout<<"DEBUG 1462 | Acyclic Coloring | Edge Vertex Map"<<endl;
02541 cout<<endl;
02542
02543 for(i=0; i<i_VertexCount; i++)
02544 {
02545 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
02546 {
02547 if(i < m_vi_Edges[j])
02548 {
02549 cout<<"Edge "<<STEP_UP(mimi2_VertexEdgeMap[i][m_vi_Edges[j]])<<"\t"<<" : "<<STEP_UP(i)<<" - "<<STEP_UP(m_vi_Edges[j])<<endl;
02550 }
02551 }
02552 }
02553
02554 cout<<endl;
02555
02556 #endif
02557
02558 i_EdgeCount = (signed) m_vi_Edges.size();
02559
02560 m_vi_VertexColors.clear();
02561 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
02562
02563 vi_CandidateColors.clear();
02564 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
02565
02566 vi_FirstSeenOne.clear();
02567 vi_FirstSeenOne.resize((unsigned) i_VertexCount, _UNKNOWN);
02568
02569 vi_FirstSeenTwo.clear();
02570 vi_FirstSeenTwo.resize((unsigned) i_VertexCount, _UNKNOWN);
02571
02572 vi_FirstSeenThree.clear();
02573 vi_FirstSeenThree.resize((unsigned) i_VertexCount, _UNKNOWN);
02574
02575 vi_FirstVisitedOne.clear();
02576 vi_FirstVisitedOne.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
02577
02578 vi_FirstVisitedTwo.clear();
02579 vi_FirstVisitedTwo.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
02580
02581
02582 #if DISJOINT_SETS == _FALSE
02583
02584 vi_MemberEdges.clear();
02585
02586 vi_EdgeSetMap.clear();
02587 vi_EdgeSetMap.resize((unsigned) i_EdgeCount/2, _UNKNOWN);
02588
02589 v2i_SetEdgeMap.clear();
02590 v2i_SetEdgeMap.resize((unsigned) i_EdgeCount/2);
02591
02592 vi_DisjointSets.clear();
02593
02594 #endif
02595
02596 #if DISJOINT_SETS == _TRUE
02597
02598 m_ds_DisjointSets.SetSize(i_EdgeCount/2);
02599
02600 #endif
02601
02602 #if VERBOSE == _TRUE
02603
02604 cout<<endl;
02605
02606 #endif
02607
02608 m_i_VertexColorCount = _UNKNOWN;
02609
02610
02611 for(i=0; i<i_VertexCount; i++)
02612 {
02613
02614 i_PresentVertex = m_vi_OrderedVertices[i];
02615
02616
02617 #if VERBOSE == _TRUE
02618
02619
02620 cout<<"DEBUG 1462 | Acyclic Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
02621
02622 #endif
02623
02624 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02625 {
02626 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02627 {
02628 continue;
02629 }
02630
02631 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
02632 }
02633
02634 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02635 {
02636 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02637 {
02638 continue;
02639 }
02640
02641 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
02642 {
02643 if(m_vi_Edges[k] == i_PresentVertex)
02644 {
02645 continue;
02646 }
02647
02648 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
02649 {
02650 continue;
02651 }
02652
02653 if(vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] == i_PresentVertex)
02654 {
02655 continue;
02656 }
02657
02658 #if DISJOINT_SETS == _TRUE
02659
02660 if(m_vi_Edges[j] < m_vi_Edges[k])
02661 {
02662 i_SetID = m_ds_DisjointSets.FindAndCompress(mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]);
02663 }
02664 else
02665 {
02666 i_SetID = m_ds_DisjointSets.FindAndCompress(mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]);
02667 }
02668 #endif
02669
02670 #if DISJOINT_SETS == _FALSE
02671
02672 if(m_vi_Edges[j] < m_vi_Edges[k])
02673 {
02674 i_SetID = vi_EdgeSetMap[mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]]];
02675 }
02676 else
02677 {
02678 i_SetID = vi_EdgeSetMap[mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]]];
02679 }
02680 #endif
02681
02682 FindCycle(i_PresentVertex, m_vi_Edges[j], m_vi_Edges[k], i_SetID, vi_CandidateColors, vi_FirstVisitedOne, vi_FirstVisitedTwo);
02683 }
02684 }
02685
02686 for(j=0; j<i_VertexCount; j++)
02687 {
02688 if(vi_CandidateColors[j] != i_PresentVertex)
02689 {
02690 m_vi_VertexColors[i_PresentVertex] = j;
02691
02692 if(m_i_VertexColorCount < j)
02693 {
02694 m_i_VertexColorCount = j;
02695 }
02696
02697 break;
02698 }
02699 }
02700
02701 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02702 {
02703 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02704 {
02705 continue;
02706 }
02707
02708 if(i_PresentVertex < m_vi_Edges[j])
02709 {
02710 i_EdgeID = mimi2_VertexEdgeMap[i_PresentVertex][m_vi_Edges[j]];
02711 }
02712 else
02713 {
02714 i_EdgeID = mimi2_VertexEdgeMap[m_vi_Edges[j]][i_PresentVertex];
02715 }
02716
02717 #if DISJOINT_SETS == _FALSE
02718
02719 vi_DisjointSets.push_back(_TRUE);
02720
02721 vi_EdgeSetMap[i_EdgeID] = STEP_DOWN((signed) vi_DisjointSets.size());
02722
02723 v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].push_back(i_EdgeID);
02724 #endif
02725
02726
02727 i_AdjacentEdgeID = UpdateSet(i_PresentVertex, m_vi_Edges[j], i_PresentVertex, mimi2_VertexEdgeMap, vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree);
02728
02729 if(i_AdjacentEdgeID != _UNKNOWN)
02730 {
02731
02732 #if DISJOINT_SETS == _TRUE
02733
02734 i_SetOneID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
02735 i_SetTwoID = m_ds_DisjointSets.FindAndCompress(i_AdjacentEdgeID);
02736
02737 #if DEBUG == 1462
02738
02739 cout<<endl;
02740 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(i_SetOneID)<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(i_SetTwoID)<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02741 cout<<endl;
02742
02743 #endif
02744
02745 m_ds_DisjointSets.UnionBySize(i_SetOneID, i_SetTwoID);
02746
02747 #endif
02748
02749 #if DISJOINT_SETS == _FALSE
02750
02751 #if DEBUG == 1462
02752
02753 cout<<endl;
02754 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(vi_EdgeSetMap[i_EdgeID])<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(vi_EdgeSetMap[i_AdjacentEdgeID])<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02755 cout<<endl;
02756
02757 #endif
02758
02759 if(v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].size() < v2i_SetEdgeMap[vi_EdgeSetMap[i_AdjacentEdgeID]].size())
02760 {
02761 i_SmallerSetID = vi_EdgeSetMap[i_EdgeID];
02762 i_BiggerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02763 }
02764 else
02765 {
02766 i_BiggerSetID = vi_EdgeSetMap[i_EdgeID];
02767 i_SmallerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02768 }
02769
02770 vi_MemberEdges.clear();
02771 vi_MemberEdges.swap(v2i_SetEdgeMap[i_BiggerSetID]);
02772
02773 vi_DisjointSets[i_BiggerSetID] = _FALSE;
02774
02775 i_MemberCount = (signed) vi_MemberEdges.size();
02776
02777 for(k=0; k<i_MemberCount; k++)
02778 {
02779 vi_EdgeSetMap[vi_MemberEdges[k]] = i_SmallerSetID;
02780
02781 v2i_SetEdgeMap[i_SmallerSetID].push_back(vi_MemberEdges[k]);
02782 }
02783 #endif
02784 }
02785 }
02786
02787
02788
02789 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
02790 {
02791 if(m_vi_VertexColors[m_vi_Edges[j]] == _UNKNOWN)
02792 {
02793 continue;
02794 }
02795
02796 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
02797 {
02798 if(m_vi_Edges[k] == i_PresentVertex)
02799 {
02800 continue;
02801 }
02802
02803 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
02804 {
02805 continue;
02806 }
02807
02808 if(m_vi_VertexColors[m_vi_Edges[k]] == m_vi_VertexColors[i_PresentVertex])
02809 {
02810 if(m_vi_Edges[j] < m_vi_Edges[k])
02811 {
02812 i_AdjacentEdgeID = mimi2_VertexEdgeMap[m_vi_Edges[j]][m_vi_Edges[k]];
02813 }
02814 else
02815 {
02816 i_AdjacentEdgeID = mimi2_VertexEdgeMap[m_vi_Edges[k]][m_vi_Edges[j]];
02817 }
02818
02819 i_EdgeID = UpdateSet(i_PresentVertex, m_vi_Edges[j], m_vi_Edges[k], mimi2_VertexEdgeMap, vi_FirstSeenOne, vi_FirstSeenTwo, vi_FirstSeenThree);
02820
02821 if(i_EdgeID != _UNKNOWN)
02822 {
02823
02824 #if DISJOINT_SETS == _TRUE
02825
02826 i_SetOneID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
02827 i_SetTwoID = m_ds_DisjointSets.FindAndCompress(i_AdjacentEdgeID);
02828
02829 #if DEBUG == 1462
02830 cout<<endl;
02831 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(i_SetOneID)<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(i_SetTwoID)<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02832 cout<<endl;
02833
02834 #endif
02835
02836 m_ds_DisjointSets.UnionBySize(i_SetOneID, i_SetTwoID);
02837
02838 #endif
02839
02840 #if DISJOINT_SETS == _FALSE
02841
02842 #if DEBUG == 1462
02843
02844 cout<<endl;
02845 cout<<"DEBUG 1462 | Acyclic Coloring | Unify Tree | Tree "<<STEP_UP(vi_EdgeSetMap[i_EdgeID])<<" (Edge "<<STEP_UP(i_EdgeID)<<") and Tree "<<STEP_UP(vi_EdgeSetMap[i_AdjacentEdgeID])<<" (Edge "<<STEP_UP(i_AdjacentEdgeID)<<")"<<endl;
02846 cout<<endl;
02847
02848 #endif
02849
02850 if(v2i_SetEdgeMap[vi_EdgeSetMap[i_EdgeID]].size() < v2i_SetEdgeMap[vi_EdgeSetMap[i_AdjacentEdgeID]].size())
02851 {
02852 i_SmallerSetID = vi_EdgeSetMap[i_EdgeID];
02853 i_BiggerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02854 }
02855 else
02856 {
02857 i_BiggerSetID = vi_EdgeSetMap[i_EdgeID];
02858 i_SmallerSetID = vi_EdgeSetMap[i_AdjacentEdgeID];
02859 }
02860
02861 vi_MemberEdges.clear();
02862 vi_MemberEdges.swap(v2i_SetEdgeMap[i_BiggerSetID]);
02863
02864 vi_DisjointSets[i_BiggerSetID] = _FALSE;
02865
02866 i_MemberCount = (signed) vi_MemberEdges.size();
02867
02868 for(l=0; l<i_MemberCount; l++)
02869 {
02870 vi_EdgeSetMap[vi_MemberEdges[l]] = i_SmallerSetID;
02871
02872 v2i_SetEdgeMap[i_SmallerSetID].push_back(vi_MemberEdges[l]);
02873 }
02874
02875 #endif
02876 }
02877 }
02878 }
02879 }
02880
02881 #if DEBUG == 1462
02882
02883 cout<<endl;
02884 cout<<"DEBUG 1462 | Acyclic Coloring | Vertex Colors | Upto Vertex "<<STEP_UP(i_PresentVertex)<<endl;
02885 cout<<endl;
02886
02887 for(j=0; j<i_VertexCount; j++)
02888 {
02889 cout<<"Vertex "<<STEP_UP(j)<<" = "<<STEP_UP(m_vi_VertexColors[j])<<endl;
02890 }
02891 #endif
02892
02893 }
02894
02895
02896 #if DEBUG == 1462
02897
02898 #if DISJOINT_SETS == _FALSE
02899
02900 i_EdgeCount = (signed) v2i_EdgeVertexMap.size();
02901
02902 cout<<endl;
02903 cout<<"DEBUG 1462 | Acyclic Coloring | Edge Set Map"<<endl;
02904 cout<<endl;
02905
02906 for(i=0; i<i_EdgeCount; i++)
02907 {
02908 cout<<STEP_UP(i)<<"\t"<<" : "<<STEP_UP(vi_EdgeSetMap[i])<<endl;
02909 }
02910
02911 cout<<endl;
02912 cout<<"DEBUG 1462 | Acyclic Coloring | Set Edge Map"<<endl;
02913 cout<<endl;
02914
02915 for(i=0; i<i_EdgeCount; i++)
02916 {
02917 i_MemberCount = (signed) v2i_SetEdgeMap[i].size();
02918
02919 if(i_MemberCount == _FALSE)
02920 {
02921 continue;
02922 }
02923
02924 cout<<"Set "<<STEP_UP(i)<<"\t"<<" : ";
02925
02926 for(j=0; j<i_MemberCount; j++)
02927 {
02928 if(j == STEP_DOWN(i_MemberCount))
02929 {
02930 cout<<STEP_UP(v2i_SetEdgeMap[i][j])<<" ("<<i_MemberCount<<")"<<endl;
02931 }
02932 else
02933 {
02934 cout<<STEP_UP(v2i_SetEdgeMap[i][j])<<", ";
02935 }
02936 }
02937 }
02938
02939 cout<<endl;
02940
02941 #endif
02942
02943 #endif
02944
02945 #if VERBOSE == _TRUE
02946
02947 cout<<endl;
02948
02949 #endif
02950
02951 #if DISJOINT_SETS == _TRUE
02952
02953
02954 vi_Sets.clear();
02955 mivi_VertexSets.clear();
02956
02957 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
02958
02959 for(i=0; i<i_VertexCount; i++)
02960 {
02961 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
02962 {
02963 if(i < m_vi_Edges[j])
02964
02965 {
02966 i_EdgeID = mimi2_VertexEdgeMap[i][m_vi_Edges[j]];
02967
02968 i_SetID = m_ds_DisjointSets.FindAndCompress(i_EdgeID);
02969
02970 if(i_SetID == i_EdgeID)
02971 {
02972 vi_Sets.push_back(i_SetID);
02973 }
02974
02975 mivi_VertexSets[i_SetID].push_back(i);
02976 mivi_VertexSets[i_SetID].push_back(m_vi_Edges[j]);
02977 }
02978 }
02979 }
02980
02981
02982 #endif
02983
02984 #if DISJOINT_SETS == _FALSE
02985
02986 vi_Sets.clear();
02987 mivi_VertexSets.clear();
02988
02989 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
02990
02991 for(i=0; i<i_VertexCount; i++)
02992 {
02993 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
02994 {
02995 if(i < m_vi_Edges[j])
02996 {
02997 i_EdgeID = mimi2_VertexEdgeMap[i][m_vi_Edges[j]];
02998
02999 i_SetID = vi_EdgeSetMap[i_EdgeID];
03000
03001 if(mivi_VertexSets[i_SetID].empty())
03002 {
03003 vi_Sets.push_back(i_SetID);
03004 }
03005
03006 mivi_VertexSets[i_SetID].push_back(i);
03007 mivi_VertexSets[i_SetID].push_back(m_vi_Edges[j]);
03008 }
03009 }
03010 }
03011
03012 #endif
03013
03014 return(_TRUE);
03015 }
03016
03017
03018
03019 int GraphColoring::CheckAcyclicColoring()
03020 {
03021 int i;
03022
03023 int i_VertexCount;
03024
03025 int i_ViolationCount;
03026
03027 vector<int> vi_TouchedVertices;
03028
03029 i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
03030
03031 i_ViolationCount = _FALSE;
03032
03033 for(i=0; i<i_VertexCount; i++)
03034 {
03035 vi_TouchedVertices.clear();
03036 vi_TouchedVertices.resize((unsigned) i_VertexCount, _FALSE);
03037
03038 vi_TouchedVertices[i] = _TRUE;
03039
03040 i_ViolationCount = SearchDepthFirst(i, i, i, vi_TouchedVertices);
03041 }
03042
03043 if(i_ViolationCount)
03044 {
03045 cout<<endl;
03046 cout<<"[Total Violations = "<<i_ViolationCount<<"]"<<endl;
03047 cout<<endl;
03048 }
03049
03050 return(i_ViolationCount);
03051 }
03052
03053
03054
03055 int GraphColoring::TriangularColoring()
03056 {
03057 if(CheckVertexColoring("TRIANGULAR"))
03058 {
03059 return(_TRUE);
03060 }
03061
03062 int i, j, k, l;
03063
03064 int _FOUND;
03065
03066 int i_VertexCount, i_VertexDegree;
03067
03068 int i_HighestVertexDegree;
03069
03070 int i_PresentVertex;
03071
03072 vector<int> vi_VertexHierarchy;
03073
03074 vector< vector<int> > v2i_VertexAdjacency;
03075
03076 i_VertexCount = (signed) m_vi_OrderedVertices.size();
03077
03078 vi_VertexHierarchy.clear();
03079 vi_VertexHierarchy.resize((unsigned) i_VertexCount);
03080
03081 v2i_VertexAdjacency.clear();
03082 v2i_VertexAdjacency.resize((unsigned) i_VertexCount);
03083
03084 for(i=0; i<i_VertexCount; i++)
03085 {
03086 vi_VertexHierarchy[m_vi_OrderedVertices[i]] = i;
03087 }
03088
03089
03090
03091
03092 m_i_VertexColorCount = _UNKNOWN;
03093
03094 for(i=0; i<i_VertexCount; i++)
03095 {
03096 i_PresentVertex = m_vi_OrderedVertices[i];
03097
03098 #if VERBOSE == _TRUE
03099
03100 cout<<"DEBUG 1464 | Triangular Coloring | Processing Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
03101
03102 #endif
03103
03104 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
03105 {
03106 v2i_VertexAdjacency[i_PresentVertex].push_back(m_vi_Edges[j]);
03107
03108 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
03109 {
03110 if(m_vi_Edges[k] == i_PresentVertex)
03111 {
03112 continue;
03113 }
03114
03115 if((vi_VertexHierarchy[m_vi_Edges[j]] > vi_VertexHierarchy[i_PresentVertex]) && (vi_VertexHierarchy[m_vi_Edges[j]] > vi_VertexHierarchy[m_vi_Edges[k]]))
03116 {
03117 _FOUND = _FALSE;
03118
03119 for(l=m_vi_Vertices[m_vi_Edges[k]]; l<m_vi_Vertices[STEP_UP(m_vi_Edges[k])]; l++)
03120 {
03121 if(m_vi_Edges[l] == i_PresentVertex)
03122 {
03123 _FOUND = TRUE;
03124
03125 break;
03126 }
03127 }
03128
03129 if(_FOUND == _FALSE)
03130 {
03131 v2i_VertexAdjacency[i_PresentVertex].push_back(m_vi_Edges[k]);
03132 }
03133 }
03134 }
03135 }
03136 }
03137
03138 m_vi_Vertices.clear();
03139 m_vi_Edges.clear();
03140
03141 i_HighestVertexDegree = _UNKNOWN;
03142
03143 for(i=0; i<i_VertexCount; i++)
03144 {
03145 m_vi_Vertices.push_back((signed) m_vi_Edges.size());
03146
03147 i_VertexDegree = (signed) v2i_VertexAdjacency[i].size();
03148
03149 if(i_HighestVertexDegree < i_VertexDegree)
03150 {
03151 i_HighestVertexDegree = i_VertexDegree;
03152 }
03153
03154 for(j=0; j<i_VertexDegree; j++)
03155 {
03156 m_vi_Edges.push_back(v2i_VertexAdjacency[i][j]);
03157 }
03158
03159 v2i_VertexAdjacency[i].clear();
03160 }
03161
03162 m_vi_Vertices.push_back((signed) m_vi_Edges.size());
03163
03164 #if DEBUG == 1464
03165
03166 int i_EdgeCount;
03167
03168 cout<<endl;
03169 cout<<"DEBUG 1464 | Graph Coloring | Induced Matrix"<<endl;
03170 cout<<endl;
03171
03172 i_VertexCount = (signed) m_vi_Vertices.size();
03173 i_EdgeCount = (signed) m_vi_Edges.size();
03174
03175 for(i=0; i<i_VertexCount; i++)
03176 {
03177 for(j=m_vi_Vertices[i]; j<m_vi_Vertices[STEP_UP(i)]; j++)
03178 {
03179 cout<<"Element["<<STEP_UP(i)<<"]["<<STEP_UP(m_vi_Edges[j])<<"] = 1"<<endl;
03180 }
03181 }
03182
03183 cout<<endl;
03184 cout<<"[Induced Vertices = "<<STEP_DOWN(i_VertexCount)<<"; Induced Edges = "<<i_EdgeCount<<"]"<<endl;
03185 cout<<endl;
03186
03187 #endif
03188
03189 SmallestLastOrdering();
03190
03191 return(DistanceOneColoring());
03192 }
03193
03194
03195
03196
03197 int GraphColoring::ModifiedTriangularColoring()
03198 {
03199 if(CheckVertexColoring("MODIFIED TRIANGULAR"))
03200 {
03201 return(_TRUE);
03202 }
03203
03204 int i, j, k;
03205
03206 int i_VertexCount;
03207
03208 int i_HighestColor;
03209
03210 int i_PresentVertex;
03211
03212 vector<int> vi_CandidateColors;
03213
03214 vector<int> vi_VertexHierarchy;
03215
03216 i_VertexCount = (signed) m_vi_OrderedVertices.size();
03217
03218 vi_VertexHierarchy.clear();
03219 vi_VertexHierarchy.resize((unsigned) i_VertexCount);
03220
03221 for(i=0; i<i_VertexCount; i++)
03222 {
03223 vi_VertexHierarchy[m_vi_OrderedVertices[i]] = i;
03224 }
03225
03226 m_vi_VertexColors.clear();
03227 m_vi_VertexColors.resize((unsigned) i_VertexCount, _UNKNOWN);
03228
03229 vi_CandidateColors.clear();
03230 vi_CandidateColors.resize((unsigned) i_VertexCount, _UNKNOWN);
03231
03232 i_HighestColor = _UNKNOWN;
03233
03234 for(i=0; i<i_VertexCount; i++)
03235 {
03236 i_PresentVertex = m_vi_OrderedVertices[i];
03237
03238 #if VERBOSE == _TRUE
03239
03240 cout<<"DEBUG 1465 | Triangular Coloring | Coloring Vertex "<<STEP_UP(i_PresentVertex)<<"/"<<i_VertexCount<<endl;
03241
03242 #endif
03243
03244 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
03245 {
03246 if(m_vi_VertexColors[m_vi_Edges[j]] != _UNKNOWN)
03247 {
03248 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[j]]] = i_PresentVertex;
03249 }
03250
03251 for(k=m_vi_Vertices[m_vi_Edges[j]]; k<m_vi_Vertices[STEP_UP(m_vi_Edges[j])]; k++)
03252 {
03253 if(m_vi_Edges[k] == i_PresentVertex)
03254 {
03255 continue;
03256 }
03257
03258 if(m_vi_VertexColors[m_vi_Edges[k]] == _UNKNOWN)
03259 {
03260 continue;
03261 }
03262
03263 if((vi_VertexHierarchy[m_vi_Edges[j]] > vi_VertexHierarchy[i_PresentVertex]) && (vi_VertexHierarchy[m_vi_Edges[j]] > vi_VertexHierarchy[m_vi_Edges[k]]))
03264 {
03265 vi_CandidateColors[m_vi_VertexColors[m_vi_Edges[k]]] = i_PresentVertex;
03266 }
03267 }
03268 }
03269
03270 for(j=0; j<i_VertexCount; j++)
03271 {
03272 if(vi_CandidateColors[j] != i_PresentVertex)
03273 {
03274 m_vi_VertexColors[i_PresentVertex] = j;
03275
03276 if(i_HighestColor < j)
03277 {
03278 i_HighestColor = j;
03279 }
03280
03281 break;
03282 }
03283 }
03284 }
03285
03286 return(_TRUE);
03287 }
03288
03289
03290 int GraphColoring::CheckTriangularColoring()
03291 {
03292 return(CheckAcyclicColoring());
03293 }
03294
03295
03296
03297 string GraphColoring::GetVertexColoringVariant()
03298 {
03299 return(m_s_VertexColoringVariant);
03300 }
03301
03302
03303
03304 int GraphColoring::GetVertexColorCount()
03305 {
03306 return(STEP_UP(m_i_VertexColorCount));
03307 }
03308
03309
03310
03311 void GraphColoring::GetVertexColors(vector<int> &output)
03312 {
03313 output = (m_vi_VertexColors);
03314 }
03315
03316
03317
03318 int GraphColoring::GetHubCount()
03319 {
03320 if(CheckVertexColoring("STAR"))
03321 {
03322 return(m_i_ColoringUnits);
03323 }
03324 else
03325 {
03326 return(_UNKNOWN);
03327 }
03328 }
03329
03330
03331
03332 int GraphColoring::GetSetCount()
03333 {
03334 if(CheckVertexColoring("ACYCLIC"))
03335 {
03336 return(m_i_ColoringUnits);
03337 }
03338 else
03339 {
03340 return(_UNKNOWN);
03341 }
03342 }
03343
03344
03345 double GraphColoring::GetVertexColoringTime()
03346 {
03347 return(m_d_ColoringTime);
03348 }
03349
03350
03351 double GraphColoring::GetVertexColoringCheckingTime()
03352 {
03353 return(m_d_CheckingTime);
03354 }
03355
03356
03357 int GraphColoring::PrintVertexColors()
03358 {
03359 int i;
03360
03361 int i_VertexCount;
03362
03363 string _SLASH("/");
03364
03365 StringTokenizer SlashTokenizer(m_s_InputFile, _SLASH);
03366
03367 m_s_InputFile = SlashTokenizer.GetLastToken();
03368
03369 i_VertexCount = (signed) m_vi_VertexColors.size();
03370
03371 cout<<endl;
03372 cout<<m_s_VertexColoringVariant<<" Coloring | "<<m_s_VertexOrderingVariant<<" Ordering | Vertex Colors | "<<m_s_InputFile<<endl;
03373 cout<<endl;
03374
03375 for(i=0; i<i_VertexCount; i++)
03376 {
03377 cout<<"Vertex "<<STEP_UP(i)<<"\t"<<" : "<<STEP_UP(m_vi_VertexColors[i])<<endl;
03378 }
03379
03380 #if STATISTICS == _TRUE
03381
03382 if(m_s_VertexColoringVariant.compare("STAR") == 0)
03383 {
03384 cout<<endl;
03385 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Stars = "<<m_i_ColoringUnits<<"]"<<endl;
03386 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03387 cout<<endl;
03388 }
03389 else
03390 if(m_s_VertexColoringVariant.compare("ACYCLIC") == 0)
03391 {
03392 cout<<endl;
03393 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Sets = "<<m_i_ColoringUnits<<"]"<<endl;
03394 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03395 cout<<endl;
03396 }
03397 else
03398 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03399 {
03400 cout<<endl;
03401 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03402 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03403 cout<<endl;
03404 }
03405 else
03406 {
03407 cout<<endl;
03408 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03409 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03410 cout<<endl;
03411 }
03412
03413 #endif
03414
03415 #if STATISTICS == _FALSE
03416
03417
03418 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03419 {
03420 cout<<endl;
03421 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03422 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Sequencing Time = "<<m_d_SequencingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03423 cout<<endl;
03424 }
03425 else
03426 {
03427 cout<<endl;
03428 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03429 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03430 cout<<endl;
03431
03432 }
03433
03434 #endif
03435
03436 return(_TRUE);
03437 }
03438
03439
03440
03441 int GraphColoring::FileVertexColors()
03442 {
03443 int i;
03444
03445 int i_VertexCount;
03446
03447 string s_InputFile, s_OutputFile;
03448
03449 string s_ColoringExtension, s_OrderingExtension;
03450
03451 string _SLASH("/");
03452
03453 ofstream OutputStream;
03454
03455
03456
03457 if(m_s_VertexOrderingVariant.compare("NATURAL") == 0)
03458 {
03459 s_OrderingExtension = ".N.";
03460 }
03461 else
03462 if(m_s_VertexOrderingVariant.compare("LARGEST_FIRST") == 0)
03463 {
03464 s_OrderingExtension = ".LF.";
03465 }
03466 else
03467 if(m_s_VertexOrderingVariant.compare("DISTANCE_TWO_LARGEST_FIRST") == 0)
03468 {
03469 s_OrderingExtension = ".D2LF.";
03470 }
03471 else
03472 if(m_s_VertexOrderingVariant.compare("SMALLEST_LAST") == 0)
03473 {
03474 s_OrderingExtension = ".SL.";
03475 }
03476 else
03477 if(m_s_VertexOrderingVariant.compare("DISTANCE_TWO_SMALLEST_LAST") == 0)
03478 {
03479 s_OrderingExtension = ".D2SL.";
03480 }
03481 else
03482 if(m_s_VertexOrderingVariant.compare("INCIDENCE_DEGREE") == 0)
03483 {
03484 s_OrderingExtension = ".ID.";
03485 }
03486 else
03487 if(m_s_VertexOrderingVariant.compare("DISTANCE_TWO_INCIDENCE_DEGREE") == 0)
03488 {
03489 s_OrderingExtension = ".D2ID.";
03490 }
03491 else
03492 {
03493 s_OrderingExtension = ".NONE.";
03494 }
03495
03496
03497
03498 if(m_s_VertexColoringVariant.compare("DISTANCE_ONE") == 0)
03499 {
03500 s_ColoringExtension = ".D1.";
03501 }
03502 else
03503 if(m_s_VertexColoringVariant.compare("DISTANCE_TWO") == 0)
03504 {
03505 s_ColoringExtension = ".D2.";
03506 }
03507 else
03508 if(m_s_VertexColoringVariant.compare("NAIVE_STAR") == 0)
03509 {
03510 s_ColoringExtension = ".NS.";
03511 }
03512 else
03513 if(m_s_VertexColoringVariant.compare("RESTRICTED_STAR") == 0)
03514 {
03515 s_ColoringExtension = ".RS.";
03516 }
03517 else
03518 if(m_s_VertexColoringVariant.compare("STAR") == 0)
03519 {
03520 s_ColoringExtension = ".S.";
03521 }
03522 else
03523 if(m_s_VertexColoringVariant.compare("ACYCLIC") == 0)
03524 {
03525 s_ColoringExtension = ".A.";
03526 }
03527 else
03528 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03529 {
03530 s_ColoringExtension = ".T.";
03531 }
03532 else
03533 {
03534 s_ColoringExtension = ".NONE.";
03535 }
03536
03537 StringTokenizer SlashTokenizer(m_s_InputFile, _SLASH);
03538
03539 s_InputFile = SlashTokenizer.GetLastToken();
03540
03541 s_OutputFile = s_InputFile;
03542 s_OutputFile += s_OrderingExtension;
03543 s_OutputFile += s_ColoringExtension;
03544 s_OutputFile += ".out";
03545
03546 OutputStream.open(s_OutputFile.c_str());
03547
03548 i_VertexCount = (signed) m_vi_VertexColors.size();
03549
03550 OutputStream<<endl;
03551 OutputStream<<m_s_VertexColoringVariant<<" Coloring | "<<m_s_VertexOrderingVariant<<" Ordering | Vertex Colors | "<<m_s_InputFile<<endl;
03552 OutputStream<<endl;
03553
03554 for(i=0; i<i_VertexCount; i++)
03555 {
03556 OutputStream<<"Vertex "<<STEP_UP(i)<<"\t"<<" : "<<STEP_UP(m_vi_VertexColors[i])<<endl;
03557 }
03558
03559 #if STATISTICS == _TRUE
03560
03561 if(m_s_VertexColoringVariant.compare("STAR") == 0)
03562 {
03563 OutputStream<<endl;
03564 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Stars = "<<m_i_ColoringUnits<<"]"<<endl;
03565 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03566 OutputStream<<endl;
03567 }
03568 else
03569 if(m_s_VertexColoringVariant.compare("ACYCLIC") == 0)
03570 {
03571 OutputStream<<endl;
03572 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Sets = "<<m_i_ColoringUnits<<"]"<<endl;
03573 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03574 OutputStream<<endl;
03575 }
03576 else
03577 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03578 {
03579 OutputStream<<endl;
03580 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03581 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03582 OutputStream<<endl;
03583 }
03584 else
03585 {
03586 OutputStream<<endl;
03587 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03588 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03589 OutputStream<<endl;
03590 }
03591
03592 #endif
03593
03594 #if STATISTICS == _FALSE
03595
03596 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03597 {
03598 OutputStream<<endl;
03599 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03600 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Sequencing Time = "<<m_d_SequencingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03601 OutputStream<<endl;
03602 }
03603 else
03604 {
03605 OutputStream<<endl;
03606 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03607 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03608 OutputStream<<endl;
03609 }
03610
03611 #endif
03612
03613 OutputStream.close();
03614
03615 return(_TRUE);
03616 }
03617
03618
03619
03620
03621 int GraphColoring::PrintVertexColoringMetrics()
03622 {
03623 cout<<endl;
03624 cout<<m_s_VertexColoringVariant<<" Coloring | "<<m_s_VertexOrderingVariant<<" Ordering | "<<m_s_InputFile<<endl;
03625 cout<<endl;
03626
03627 #if STATISTICS == _TRUE
03628
03629 if(m_s_VertexColoringVariant.compare("STAR") == 0)
03630 {
03631 cout<<endl;
03632 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Stars = "<<m_i_ColoringUnits<<"]"<<endl;
03633 cout<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()/2<<"]"<<endl;
03634 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03635 cout<<endl;
03636 }
03637 else
03638 if(m_s_VertexColoringVariant.compare("ACYCLIC") == 0)
03639 {
03640 cout<<endl;
03641 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Sets = "<<m_i_ColoringUnits<<"]"<<endl;
03642 cout<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()/2<<"]"<<endl;
03643 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03644 cout<<endl;
03645 }
03646 else
03647 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03648 {
03649 cout<<endl;
03650 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03651 cout<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03652 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03653 cout<<endl;
03654 }
03655 else
03656 {
03657 cout<<endl;
03658 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03659 cout<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()/2<<"]"<<endl;
03660 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03661 cout<<endl;
03662 }
03663
03664 #endif
03665
03666 #if STATISTICS == _FALSE
03667
03668 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03669 {
03670 cout<<endl;
03671 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03672 cout<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()/2<<"]"<<endl;
03673 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Sequencing Time = "<<m_d_SequencingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03674 cout<<endl;
03675 }
03676 else
03677 {
03678 cout<<endl;
03679 cout<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03680 cout<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()/2<<"]"<<endl;
03681 cout<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03682 cout<<endl;
03683
03684 }
03685
03686 #endif
03687
03688 return(_TRUE);
03689
03690 }
03691
03692
03693 int GraphColoring::FileVertexColoringMetrics()
03694 {
03695 string s_InputFile, s_OutputFile;
03696
03697 string s_ColoringExtension, s_OrderingExtension;
03698
03699 string _SLASH("/");
03700
03701 ofstream OutputStream;
03702
03703
03704
03705 if(m_s_VertexOrderingVariant.compare("ALL") == 0)
03706 {
03707 s_OrderingExtension = ".ALL.";
03708 }
03709 else
03710 if(m_s_VertexOrderingVariant.compare("NATURAL") == 0)
03711 {
03712 s_OrderingExtension = ".N.";
03713 }
03714 else
03715 if(m_s_VertexOrderingVariant.compare("LARGEST FIRST") == 0)
03716 {
03717 s_OrderingExtension = ".LF.";
03718 }
03719 else
03720 if(m_s_VertexOrderingVariant.compare("DISTANCE TWO LARGEST FIRST") == 0)
03721 {
03722 s_OrderingExtension = ".D2LF.";
03723 }
03724 else
03725 if(m_s_VertexOrderingVariant.compare("SMALLEST LAST") == 0)
03726 {
03727 s_OrderingExtension = ".SL.";
03728 }
03729 else
03730 if(m_s_VertexOrderingVariant.compare("DISTANCE TWO SMALLEST LAST") == 0)
03731 {
03732 s_OrderingExtension = ".D2SL.";
03733 }
03734 else
03735 if(m_s_VertexOrderingVariant.compare("INCIDENCE DEGREE") == 0)
03736 {
03737 s_OrderingExtension = ".ID.";
03738 }
03739 else
03740 if(m_s_VertexOrderingVariant.compare("DISTANCE TWO INCIDENCE DEGREE") == 0)
03741 {
03742 s_OrderingExtension = ".D2ID.";
03743 }
03744 else
03745 {
03746 s_OrderingExtension = ".NONE.";
03747 }
03748
03749
03750
03751 if(m_s_VertexColoringVariant.compare("ALL") == 0)
03752 {
03753 s_ColoringExtension = ".ALL.";
03754 }
03755 else
03756 if(m_s_VertexColoringVariant.compare("DISTANCE ONE") == 0)
03757 {
03758 s_ColoringExtension = ".D1.";
03759 }
03760 else
03761 if(m_s_VertexColoringVariant.compare("DISTANCE TWO") == 0)
03762 {
03763 s_ColoringExtension = ".D2.";
03764 }
03765 else
03766 if(m_s_VertexColoringVariant.compare("NAIVE STAR") == 0)
03767 {
03768 s_ColoringExtension = ".NS.";
03769 }
03770 else
03771 if(m_s_VertexColoringVariant.compare("RESTRICTED STAR") == 0)
03772 {
03773 s_ColoringExtension = ".RS.";
03774 }
03775 else
03776 if(m_s_VertexColoringVariant.compare("STAR") == 0)
03777 {
03778 s_ColoringExtension = ".S.";
03779 }
03780 else
03781 if(m_s_VertexColoringVariant.compare("ACYCLIC") == 0)
03782 {
03783 s_ColoringExtension = ".A.";
03784 }
03785 else
03786 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03787 {
03788 s_ColoringExtension = ".T.";
03789 }
03790 else
03791 {
03792 s_ColoringExtension = ".NONE.";
03793 }
03794
03795 StringTokenizer SlashTokenizer(m_s_InputFile, _SLASH);
03796
03797 s_InputFile = SlashTokenizer.GetLastToken();
03798
03799 s_OutputFile = s_InputFile;
03800 s_OutputFile += s_OrderingExtension;
03801 s_OutputFile += s_ColoringExtension;
03802 s_OutputFile += ".out";
03803
03804 OutputStream.open(s_OutputFile.c_str(), ios::app);
03805
03806 OutputStream<<endl;
03807 OutputStream<<m_s_VertexColoringVariant<<" Coloring | "<<m_s_VertexOrderingVariant<<" Ordering | "<<m_s_InputFile<<endl;
03808 OutputStream<<endl;
03809
03810 #if STATISTICS == _TRUE
03811
03812 if(m_s_VertexColoringVariant.compare("STAR") == 0)
03813 {
03814 OutputStream<<endl;
03815 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Stars = "<<m_i_ColoringUnits<<"]"<<endl;
03816 OutputStream<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03817 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03818 OutputStream<<endl;
03819 }
03820 else
03821 if(m_s_VertexColoringVariant.compare("ACYCLIC") == 0)
03822 {
03823 OutputStream<<endl;
03824 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"; Total Sets = "<<m_i_ColoringUnits<<"]"<<endl;
03825 OutputStream<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03826 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03827 OutputStream<<endl;
03828 }
03829 else
03830 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03831 {
03832 OutputStream<<endl;
03833 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03834 OutputStream<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03835 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03836 OutputStream<<endl;
03837 }
03838 else
03839 {
03840 OutputStream<<endl;
03841 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03842 OutputStream<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03843 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03844 OutputStream<<endl;
03845 }
03846
03847 #endif
03848
03849 #if STATISTICS == _FALSE
03850
03851 if(m_s_VertexColoringVariant.compare("TRIANGULAR") == 0)
03852 {
03853 OutputStream<<endl;
03854 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03855 OutputStream<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03856 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Sequencing Time = "<<m_d_SequencingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03857 OutputStream<<endl;
03858 }
03859 else
03860 {
03861 OutputStream<<endl;
03862 OutputStream<<"[Total Colors = "<<STEP_UP(m_i_VertexColorCount)<<"]"<<endl;
03863 OutputStream<<"[Vertex Count = "<<STEP_DOWN(m_vi_Vertices.size())<<"; Edge Count = "<<m_vi_Edges.size()<<"]"<<endl;
03864 OutputStream<<"[Ordering Time = "<<m_d_OrderingTime<<"; Coloring Time = "<<m_d_ColoringTime<<"]"<<endl;
03865 OutputStream<<endl;
03866 }
03867
03868 #endif
03869
03870 OutputStream.close();
03871
03872 return(_TRUE);
03873
03874 }
03875
03876
03877
03878 void GraphColoring::PrintVertexColorClasses()
03879 {
03880 if(CalculateVertexColorClasses() != _TRUE)
03881 {
03882 cout<<endl;
03883 cout<<"Vertex Color Classes | "<<m_s_VertexColoringVariant<<" Coloring | "<<m_s_VertexOrderingVariant<<" Ordering | "<<m_s_InputFile<<" | Vertex Colors Not Set"<<endl;
03884 cout<<endl;
03885
03886 return;
03887 }
03888
03889 cout<<endl;
03890 cout<<"Vertex Color Classes | "<<m_s_VertexColoringVariant<<" Coloring | "<<m_s_VertexOrderingVariant<<" Ordering | "<<m_s_InputFile<<endl;
03891 cout<<endl;
03892
03893 int i_TotalVertexColors = STEP_UP(m_i_VertexColorCount);
03894
03895 for(int i = 0; i < i_TotalVertexColors; i++)
03896 {
03897 if(m_vi_VertexColorFrequency[i] <= 0)
03898 {
03899 continue;
03900 }
03901
03902 cout<<"Color "<<STEP_UP(i)<<" : "<<m_vi_VertexColorFrequency[i]<<endl;
03903 }
03904
03905 cout<<endl;
03906 cout<<"[Largest Color Class : "<<STEP_UP(m_i_LargestColorClass)<<"; Largest Color Class Size : "<<m_i_LargestColorClassSize<<"]"<<endl;
03907 cout<<"[Smallest Color Class : "<<STEP_UP(m_i_SmallestColorClass)<<"; Smallest Color Class Size : "<<m_i_SmallestColorClassSize<<"]"<<endl;
03908 cout<<"[Average Color Class Size : "<<m_d_AverageColorClassSize<<"]"<<endl;
03909 cout<<endl;
03910
03911 return;
03912 }
03913
03914 double** GraphColoring::GetSeedMatrix(int* i_SeedRowCount, int* i_SeedColumnCount) {
03915
03916 if(seed_available) Seed_reset();
03917
03918 int i_size = m_vi_VertexColors.size();
03919 int i_num_of_colors = m_i_VertexColorCount + 1;
03920 (*i_SeedRowCount) = i_size;
03921 (*i_SeedColumnCount) = i_num_of_colors;
03922 if(i_num_of_colors == 0 || i_size == 0) {return NULL;}
03923 double** Seed = new double*[i_size];
03924
03925
03926 for (int i=0; i<i_size; i++) {
03927 Seed[i] = new double[i_num_of_colors];
03928 for(int j=0; j<i_num_of_colors; j++) Seed[i][j]=0.;
03929 }
03930
03931
03932 for (int i=0; i < i_size; i++) {
03933 Seed[i][m_vi_VertexColors[i]] = 1.;
03934 }
03935
03936 seed_available = true;
03937 i_seed_rowCount = *i_SeedRowCount;
03938 dp2_Seed = Seed;
03939
03940 return Seed;
03941 }
03942
03943 double** GraphColoring::GetSeedMatrix_unmanaged(int* i_SeedRowCount, int* i_SeedColumnCount) {
03944
03945 int i_size = m_vi_VertexColors.size();
03946 int i_num_of_colors = m_i_VertexColorCount + 1;
03947 (*i_SeedRowCount) = i_size;
03948 (*i_SeedColumnCount) = i_num_of_colors;
03949 if(i_num_of_colors == 0 || i_size == 0) {return NULL;}
03950 double** Seed = new double*[i_size];
03951
03952
03953 for (int i=0; i<i_size; i++) {
03954 Seed[i] = new double[i_num_of_colors];
03955 for(int j=0; j<i_num_of_colors; j++) Seed[i][j]=0.;
03956 }
03957
03958
03959 for (int i=0; i < i_size; i++) {
03960 Seed[i][m_vi_VertexColors[i]] = 1.;
03961 }
03962
03963 return Seed;
03964 }
03965
03966 void GraphColoring::Seed_init() {
03967 seed_available = false;
03968
03969 i_seed_rowCount = 0;
03970 dp2_Seed = NULL;
03971 }
03972
03973 void GraphColoring::Seed_reset() {
03974 if(seed_available) {
03975 seed_available = false;
03976
03977 free_2DMatrix(dp2_Seed, i_seed_rowCount);
03978 dp2_Seed = NULL;
03979 i_seed_rowCount = 0;
03980 }
03981 }
03982
03983
03984 int GraphColoring::CheckQuickDistanceTwoColoring(int Verbose)
03985 {
03986 if (m_i_MaximumVertexDegree <= STEP_UP(m_i_VertexColorCount)) return 0;
03987
03988
03989
03990 if (Verbose < 1) return 1;
03991
03992
03993 int i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
03994 int i_VertexWithMaxDegree = -1;
03995 int i_MaximumVertexDegree = -1;
03996 int i_VertexDegree = -1;
03997
03998 for (int i = 0; i < i_VertexCount; i++)
03999 {
04000 i_VertexDegree = m_vi_Vertices[i + 1] - m_vi_Vertices[i];
04001
04002 if(i_MaximumVertexDegree < i_VertexDegree)
04003 {
04004 i_MaximumVertexDegree = i_VertexDegree;
04005 i_VertexWithMaxDegree = i;
04006 }
04007 }
04008
04009 cout<<"VertexWithMaxDegree = "<< i_VertexWithMaxDegree <<"; MaximumVertexDegree = "<< i_MaximumVertexDegree <<endl;
04010 if (Verbose < 2) return 1;
04011
04012 for (int i = m_vi_Vertices[i_VertexWithMaxDegree]; i < m_vi_Vertices[i_VertexWithMaxDegree + 1] - 1; i++) {
04013
04014 for (int j = i + 1; j < m_vi_Vertices[i_VertexWithMaxDegree + 1]; j++) {
04015 if (m_vi_VertexColors[m_vi_Edges[i]] == m_vi_VertexColors[m_vi_Edges[j]])
04016 printf("\t m_vi_VertexColors[m_vi_Edges[i(%d)](%d)](%d) == m_vi_VertexColors[m_vi_Edges[j(%d)](%d)](%d)\n", i, m_vi_Edges[i], m_vi_VertexColors[m_vi_Edges[i]], j, m_vi_Edges[j], m_vi_VertexColors[m_vi_Edges[j]]);
04017 }
04018 }
04019
04020 return 1;
04021 }
04022
04023 int GraphColoring::CheckDistanceTwoColoring(int Verbose) {
04024
04025 int j = 0, k = 0;
04026 int i_PresentVertex, i_DistanceOneVertex, i_DistanceTwoVertex;
04027 int i_VertexCount = STEP_DOWN((signed) m_vi_Vertices.size());
04028
04029 for(i_PresentVertex=0; i_PresentVertex<i_VertexCount; i_PresentVertex++)
04030 {
04031
04032
04033 for(j=m_vi_Vertices[i_PresentVertex]; j<m_vi_Vertices[STEP_UP(i_PresentVertex)]; j++)
04034 {
04035 i_DistanceOneVertex = m_vi_Edges[j];
04036 if(m_vi_VertexColors[i_PresentVertex] == m_vi_VertexColors[i_DistanceOneVertex]) {
04037
04038 if (Verbose < 1) return 1;
04039 printf("D1 VIOLATION. m_vi_VertexColors[i_PresentVertex(%d)](%d) == m_vi_VertexColors[i_DistanceOneVertex(%d)](%d) \n", i_PresentVertex, m_vi_VertexColors[i_PresentVertex], i_DistanceOneVertex, m_vi_VertexColors[i_DistanceOneVertex]);
04040 if (Verbose < 2) return 1;
04041 }
04042
04043
04044 for(k=m_vi_Vertices[i_DistanceOneVertex]; k<m_vi_Vertices[STEP_UP(i_DistanceOneVertex)]; k++)
04045 {
04046 i_DistanceTwoVertex = m_vi_Edges[k];
04047
04048 if(i_DistanceTwoVertex == i_PresentVertex) continue;
04049
04050 if(m_vi_VertexColors[i_PresentVertex] == m_vi_VertexColors[i_DistanceTwoVertex]) {
04051
04052 if (Verbose < 1) return 1;
04053 printf("D2 VIOLATION. m_vi_VertexColors[i_PresentVertex(%d)](%d) == m_vi_VertexColors[i_DistanceTwoVertex(%d)](%d) \n", i_PresentVertex, m_vi_VertexColors[i_PresentVertex], i_DistanceTwoVertex, m_vi_VertexColors[i_DistanceTwoVertex]);
04054 printf("\t i_PresentVertex(%d) and i_DistanceTwoVertex(%d) connected through i_DistanceOneVertex(%d) \n", i_PresentVertex, i_DistanceTwoVertex, i_DistanceOneVertex);
04055 if (Verbose < 2) return 1;
04056 }
04057 }
04058 }
04059 }
04060 return 0;
04061 }
04062
04063 void GraphColoring::SetStringVertexColoringVariant(string s)
04064 {
04065 m_s_VertexColoringVariant = s;
04066 }
04067 }