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 BipartiteGraphPartialOrdering::CheckVertexOrdering(string s_VertexOrderingVariant)
00029 {
00030 if(m_s_VertexOrderingVariant.compare(s_VertexOrderingVariant) == 0)
00031 {
00032 return(_TRUE);
00033 }
00034
00035 if(m_s_VertexOrderingVariant.compare("ALL") != 0)
00036 {
00037 m_s_VertexOrderingVariant = s_VertexOrderingVariant;
00038 }
00039
00040 return(_FALSE);
00041 }
00042
00043
00044
00045 BipartiteGraphPartialOrdering::BipartiteGraphPartialOrdering()
00046 {
00047 Clear();
00048 }
00049
00050
00051
00052 BipartiteGraphPartialOrdering::~BipartiteGraphPartialOrdering()
00053 {
00054 Clear();
00055 }
00056
00057
00058
00059 void BipartiteGraphPartialOrdering::Clear()
00060 {
00061 BipartiteGraphInputOutput::Clear();
00062
00063 m_d_OrderingTime = _UNKNOWN;
00064
00065 m_s_VertexOrderingVariant.clear();
00066
00067 m_vi_OrderedVertices.clear();
00068
00069 return;
00070 }
00071
00072
00073
00074 void BipartiteGraphPartialOrdering::Reset()
00075 {
00076 m_d_OrderingTime = _UNKNOWN;
00077
00078 m_s_VertexOrderingVariant.clear();
00079
00080 m_vi_OrderedVertices.clear();
00081
00082 return;
00083 }
00084
00085
00086
00087 int BipartiteGraphPartialOrdering::RowNaturalOrdering()
00088 {
00089 if(CheckVertexOrdering("ROW_NATURAL"))
00090 {
00091 return(_TRUE);
00092 }
00093
00094 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00095
00096 m_vi_OrderedVertices.clear();
00097 m_vi_OrderedVertices.reserve((unsigned) i_LeftVertexCount);
00098
00099 for(int i = 0; i<i_LeftVertexCount; i++)
00100 {
00101 m_vi_OrderedVertices.push_back(i);
00102 }
00103
00104 return(_TRUE);
00105 }
00106
00107
00108
00109 int BipartiteGraphPartialOrdering::ColumnNaturalOrdering()
00110 {
00111 if(CheckVertexOrdering("COLUMN_NATURAL"))
00112 {
00113 return(_TRUE);
00114 }
00115
00116 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00117 int i_RightVertexCount = STEP_DOWN((signed) m_vi_RightVertices.size());
00118
00119 m_vi_OrderedVertices.clear();
00120 m_vi_OrderedVertices.reserve((unsigned) i_RightVertexCount);
00121
00122 for(int i = 0; i < i_RightVertexCount; i++)
00123 {
00124 m_vi_OrderedVertices.push_back(i + i_LeftVertexCount);
00125 }
00126
00127 return(_TRUE);
00128 }
00129
00130 int BipartiteGraphPartialOrdering::RowRandomOrdering() {
00131 if(CheckVertexOrdering("ROW_RANDOM"))
00132 {
00133 return(_TRUE);
00134 }
00135
00136 m_s_VertexOrderingVariant = "ROW_RANDOM";
00137
00138 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00139
00140 m_vi_OrderedVertices.clear();
00141 m_vi_OrderedVertices.resize((unsigned) i_LeftVertexCount);
00142
00143 for(unsigned int i = 0; i<i_LeftVertexCount; i++) {
00144 m_vi_OrderedVertices[i] = i;
00145 }
00146
00147 randomOrdering(m_vi_OrderedVertices);
00148
00149 return(_TRUE);
00150 }
00151
00152 int BipartiteGraphPartialOrdering::ColumnRandomOrdering() {
00153 if(CheckVertexOrdering("COLUMN_RANDOM"))
00154 {
00155 return(_TRUE);
00156 }
00157
00158 m_s_VertexOrderingVariant = "COLUMN_RANDOM";
00159
00160 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00161 int i_RightVertexCount = STEP_DOWN((signed) m_vi_RightVertices.size());
00162
00163 m_vi_OrderedVertices.clear();
00164 m_vi_OrderedVertices.resize((unsigned) i_RightVertexCount);
00165
00166 for(unsigned int i = 0; i<i_RightVertexCount; i++) {
00167 m_vi_OrderedVertices[i] = i + i_LeftVertexCount;
00168 }
00169
00170 randomOrdering(m_vi_OrderedVertices);
00171
00172 return(_TRUE);
00173 }
00174
00175
00176 int BipartiteGraphPartialOrdering::RowLargestFirstOrdering()
00177 {
00178 if(CheckVertexOrdering("ROW_LARGEST_FIRST"))
00179 {
00180 return(_TRUE);
00181 }
00182
00183 int i, j, k;
00184 int i_DegreeCount = 0;
00185 int i_VertexCount, i_Current;
00186 vector<int> vi_Visited;
00187 vector< vector<int> > vvi_GroupedVertexDegree;
00188
00189 i_VertexCount = (int)m_vi_LeftVertices.size () - 1;
00190 m_i_MaximumVertexDegree = 0;
00191 m_i_MinimumVertexDegree = i_VertexCount;
00192 vvi_GroupedVertexDegree.resize ( i_VertexCount );
00193 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00194
00195
00196 for ( i=0; i<i_VertexCount; ++i )
00197 {
00198
00199
00200
00201 i_DegreeCount = 0;
00202
00203 for ( j=m_vi_LeftVertices [i]; j<m_vi_LeftVertices [i+1]; ++j )
00204 {
00205 i_Current = m_vi_Edges [j];
00206
00207 for ( k=m_vi_RightVertices [i_Current]; k<m_vi_RightVertices [i_Current+1]; ++k )
00208 {
00209
00210
00211 if ( m_vi_Edges [k] != i && vi_Visited [m_vi_Edges [k]] != i )
00212 {
00213 ++i_DegreeCount;
00214
00215 vi_Visited [m_vi_Edges [k]] = i;
00216 }
00217 }
00218 }
00219
00220 vvi_GroupedVertexDegree [i_DegreeCount].push_back ( i );
00221
00222 if ( m_i_MaximumVertexDegree < i_DegreeCount )
00223 {
00224 m_i_MaximumVertexDegree = i_DegreeCount;
00225 }
00226 else if (m_i_MinimumVertexDegree > i_DegreeCount)
00227 {
00228 m_i_MinimumVertexDegree = i_DegreeCount;
00229 }
00230 }
00231
00232 if(i_VertexCount <2) m_i_MinimumVertexDegree = i_DegreeCount;
00233
00234
00235 m_vi_OrderedVertices.clear ();
00236
00237 for ( i=m_i_MaximumVertexDegree; i>=m_i_MinimumVertexDegree; --i )
00238 {
00239
00240
00241 j = (int)vvi_GroupedVertexDegree [i].size ();
00242
00243 for ( k=0; k<j; ++k )
00244 {
00245 m_vi_OrderedVertices.push_back ( vvi_GroupedVertexDegree [i][k] );
00246 }
00247 }
00248
00249 return(_TRUE);
00250 }
00251
00252
00253
00254
00255 int BipartiteGraphPartialOrdering::ColumnLargestFirstOrdering()
00256 {
00257 if(CheckVertexOrdering("COLUMN_LARGEST_FIRST"))
00258 {
00259 return(_TRUE);
00260 }
00261
00262 int i, j, k;
00263 int i_DegreeCount = 0;
00264 int i_VertexCount, i_Current;
00265 vector<int> vi_Visited;
00266 vector< vector<int> > vvi_GroupedVertexDegree;
00267
00268 i_VertexCount = (int)m_vi_RightVertices.size () - 1;
00269 m_i_MaximumVertexDegree = 0;
00270 m_i_MinimumVertexDegree = i_VertexCount;
00271 vvi_GroupedVertexDegree.resize ( i_VertexCount );
00272 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00273
00274 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00275
00276
00277 for ( i=0; i<i_VertexCount; ++i )
00278 {
00279
00280
00281
00282 i_DegreeCount = 0;
00283
00284 for ( j=m_vi_RightVertices [i]; j<m_vi_RightVertices [i+1]; ++j )
00285 {
00286 i_Current = m_vi_Edges [j];
00287
00288 for ( k=m_vi_LeftVertices [i_Current]; k<m_vi_LeftVertices [i_Current+1]; ++k )
00289 {
00290
00291
00292 if ( m_vi_Edges [k] != i && vi_Visited [m_vi_Edges [k]] != i )
00293 {
00294 ++i_DegreeCount;
00295
00296 vi_Visited [m_vi_Edges [k]] = i;
00297 }
00298 }
00299 }
00300 vvi_GroupedVertexDegree [i_DegreeCount].push_back ( i );
00301
00302 if ( m_i_MaximumVertexDegree < i_DegreeCount )
00303 {
00304 m_i_MaximumVertexDegree = i_DegreeCount;
00305 }
00306 else if (m_i_MinimumVertexDegree > i_DegreeCount)
00307 {
00308 m_i_MinimumVertexDegree = i_DegreeCount;
00309 }
00310 }
00311
00312 if(i_VertexCount <2) m_i_MinimumVertexDegree = i_DegreeCount;
00313
00314
00315 m_vi_OrderedVertices.clear ();
00316
00317
00318 for ( i=m_i_MaximumVertexDegree; i>=m_i_MinimumVertexDegree; --i )
00319 {
00320
00321
00322 j = (int)vvi_GroupedVertexDegree [i].size ();
00323
00324 for ( k=0; k<j; ++k )
00325 {
00326 m_vi_OrderedVertices.push_back ( vvi_GroupedVertexDegree [i][k] + i_LeftVertexCount);
00327 }
00328 }
00329
00330 return(_TRUE);
00331 }
00332
00333
00334
00335
00336 int BipartiteGraphPartialOrdering::RowSmallestLastOrdering()
00337 {
00338 if(CheckVertexOrdering("ROW_SMALLEST_LAST"))
00339 {
00340 return(_TRUE);
00341 }
00342
00343 int i, j, k, u, l;
00344 int i_Current;
00345 int i_SelectedVertex, i_SelectedVertexCount;
00346 int i_VertexCount;
00347 int i_VertexCountMinus1;
00348 int i_HighestInducedVertexDegree, i_InducedVertexDegree;
00349 vector<int> vi_InducedVertexDegree;
00350 vector<int> vi_Visited;
00351 vector< vector<int> > vvi_GroupedInducedVertexDegree;
00352 vector< int > vi_VertexLocation;
00353
00354
00355
00356 i_SelectedVertex = _UNKNOWN;
00357 i_VertexCount = (int)m_vi_LeftVertices.size () - 1;
00358 i_VertexCountMinus1 = i_VertexCount - 1;
00359 i_HighestInducedVertexDegree = 0;
00360 vi_Visited.clear();
00361 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00362 m_vi_OrderedVertices.clear();
00363 m_vi_OrderedVertices.resize(i_VertexCount, _UNKNOWN);
00364
00365 vi_InducedVertexDegree.clear();
00366 vi_InducedVertexDegree.reserve((unsigned) i_VertexCount);
00367 vvi_GroupedInducedVertexDegree.clear();
00368 vvi_GroupedInducedVertexDegree.resize((unsigned) i_VertexCount);
00369 vi_VertexLocation.clear();
00370 vi_VertexLocation.reserve((unsigned) i_VertexCount);
00371
00372 for ( i=0; i<i_VertexCount; ++i )
00373 {
00374
00375
00376
00377 i_InducedVertexDegree = 0;
00378
00379 for ( j=m_vi_LeftVertices[i]; j<m_vi_LeftVertices[i+1]; ++j )
00380 {
00381 i_Current = m_vi_Edges [j];
00382
00383 for (k=m_vi_RightVertices[i_Current]; k<m_vi_RightVertices[i_Current+1]; ++k)
00384 {
00385
00386
00387 if ( m_vi_Edges [k] != i && vi_Visited [m_vi_Edges [k]] != i )
00388 {
00389 ++i_InducedVertexDegree;
00390
00391 vi_Visited [m_vi_Edges [k]] = i;
00392 }
00393 }
00394 }
00395
00396
00397 vi_InducedVertexDegree.push_back ( i_InducedVertexDegree );
00398
00399
00400 vvi_GroupedInducedVertexDegree [i_InducedVertexDegree].push_back ( i );
00401
00402 vi_VertexLocation.push_back(vvi_GroupedInducedVertexDegree[i_InducedVertexDegree].size() - 1);
00403
00404
00405 if ( i_HighestInducedVertexDegree < i_InducedVertexDegree )
00406 {
00407 i_HighestInducedVertexDegree = i_InducedVertexDegree;
00408 }
00409 }
00410
00411 i_SelectedVertexCount = 0;
00412
00413 vi_Visited.clear ();
00414 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00415
00416
00417 int iMin = 1;
00418
00419
00420
00421 while ( i_SelectedVertexCount < i_VertexCount )
00422 {
00423 if(iMin != 0 && vvi_GroupedInducedVertexDegree[iMin - 1].size() != _FALSE)
00424 iMin--;
00425
00426
00427 for ( i=iMin; i<(i_HighestInducedVertexDegree+1); ++i )
00428 {
00429
00430 if ( vvi_GroupedInducedVertexDegree[i].size () != 0 )
00431 {
00432 i_SelectedVertex = vvi_GroupedInducedVertexDegree[i].back ();
00433
00434 vvi_GroupedInducedVertexDegree[i].pop_back();
00435 break;
00436 }
00437 else
00438 iMin++;
00439 }
00440
00441
00442
00443 for ( i=m_vi_LeftVertices [i_SelectedVertex]; i<m_vi_LeftVertices [i_SelectedVertex+1]; ++i )
00444 {
00445
00446 i_Current = m_vi_Edges [i];
00447
00448
00449 for ( j=m_vi_RightVertices [i_Current]; j<m_vi_RightVertices [i_Current+1]; ++j )
00450 {
00451
00452
00453 if ( m_vi_Edges [j] == i_SelectedVertex || vi_Visited [m_vi_Edges [j]] == i_SelectedVertex )
00454 {
00455
00456 continue;
00457 }
00458
00459 u = m_vi_Edges[j];
00460
00461
00462 if ( vi_InducedVertexDegree [u] == _UNKNOWN )
00463 {
00464
00465 continue;
00466 }
00467
00468
00469 vi_Visited [u] = i_SelectedVertex;
00470
00471
00472
00473 if(vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].size() > 1)
00474 {
00475 l = vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].back();
00476
00477 vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]][vi_VertexLocation[u]] = l;
00478
00479
00480 vi_VertexLocation[l] = vi_VertexLocation[u];
00481 }
00482
00483 vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].pop_back();
00484
00485
00486 vi_InducedVertexDegree[u]--;
00487
00488
00489 vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].push_back(u);
00490
00491
00492 vi_VertexLocation[u] = vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].size() - 1;
00493
00494
00495
00496 }
00497 }
00498
00499
00500
00501 vi_InducedVertexDegree [i_SelectedVertex] = _UNKNOWN;
00502
00503 m_vi_OrderedVertices[i_VertexCountMinus1 - i_SelectedVertexCount] = i_SelectedVertex;
00504
00505
00506 ++i_SelectedVertexCount;
00507 }
00508
00509
00510 vi_InducedVertexDegree.clear();
00511 vi_VertexLocation.clear();
00512 vvi_GroupedInducedVertexDegree.clear();
00513
00514 return(_TRUE);
00515 }
00516
00517
00518
00519 int BipartiteGraphPartialOrdering::ColumnSmallestLastOrdering()
00520 {
00521 if(CheckVertexOrdering("COLUMN_SMALLEST_LAST"))
00522 {
00523 return(_TRUE);
00524 }
00525
00526 int i, j, k, u, l;
00527 int i_Current;
00528 int i_SelectedVertex, i_SelectedVertexCount;
00529 int i_VertexCount;
00530 int i_VertexCountMinus1;
00531 int i_HighestInducedVertexDegree, i_InducedVertexDegree;
00532 vector<int> vi_InducedVertexDegree;
00533 vector<int> vi_Visited;
00534 vector< vector<int> > vvi_GroupedInducedVertexDegree;
00535 vector< int > vi_VertexLocation;
00536
00537
00538 i_SelectedVertex = _UNKNOWN;
00539 i_VertexCount = (int)m_vi_RightVertices.size () - 1;
00540 i_VertexCountMinus1 = i_VertexCount - 1;
00541 i_HighestInducedVertexDegree = 0;
00542 vi_Visited.clear();
00543 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00544 m_vi_OrderedVertices.clear();
00545 m_vi_OrderedVertices.resize(i_VertexCount, _UNKNOWN);
00546
00547 vi_InducedVertexDegree.clear();
00548 vi_InducedVertexDegree.reserve((unsigned) i_VertexCount);
00549 vvi_GroupedInducedVertexDegree.clear();
00550 vvi_GroupedInducedVertexDegree.resize((unsigned) i_VertexCount);
00551 vi_VertexLocation.clear();
00552 vi_VertexLocation.reserve((unsigned) i_VertexCount);
00553
00554 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00555
00556 for ( i=0; i<i_VertexCount; ++i )
00557 {
00558
00559
00560
00561 i_InducedVertexDegree = 0;
00562
00563 for ( j=m_vi_RightVertices [i]; j<m_vi_RightVertices [i+1]; ++j )
00564 {
00565 i_Current = m_vi_Edges [j];
00566
00567 for ( k=m_vi_LeftVertices [i_Current]; k<m_vi_LeftVertices [i_Current+1]; ++k )
00568 {
00569
00570
00571 if ( m_vi_Edges [k] != i && vi_Visited [m_vi_Edges [k]] != i )
00572 {
00573 ++i_InducedVertexDegree;
00574
00575 vi_Visited [m_vi_Edges [k]] = i;
00576 }
00577 }
00578 }
00579
00580
00581 vi_InducedVertexDegree.push_back ( i_InducedVertexDegree );
00582
00583
00584 vvi_GroupedInducedVertexDegree [i_InducedVertexDegree].push_back ( i );
00585
00586 vi_VertexLocation.push_back(vvi_GroupedInducedVertexDegree[i_InducedVertexDegree].size() - 1);
00587
00588
00589 if ( i_HighestInducedVertexDegree < i_InducedVertexDegree )
00590 {
00591 i_HighestInducedVertexDegree = i_InducedVertexDegree;
00592 }
00593 }
00594
00595 i_SelectedVertexCount = 0;
00596
00597 vi_Visited.clear ();
00598 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00599
00600
00601 int iMin = 1;
00602
00603
00604
00605 while ( i_SelectedVertexCount < i_VertexCount )
00606 {
00607 if(iMin != 0 && vvi_GroupedInducedVertexDegree[iMin - 1].size() != _FALSE)
00608 iMin--;
00609
00610
00611 for ( i=iMin; i<(i_HighestInducedVertexDegree+1); ++i )
00612 {
00613
00614 if ( vvi_GroupedInducedVertexDegree[i].size () != 0 )
00615 {
00616 i_SelectedVertex = vvi_GroupedInducedVertexDegree[i].back ();
00617
00618 vvi_GroupedInducedVertexDegree[i].pop_back();
00619 break;
00620 }
00621 else
00622 iMin++;
00623 }
00624
00625
00626
00627 for ( i=m_vi_RightVertices [i_SelectedVertex]; i<m_vi_RightVertices [i_SelectedVertex+1]; ++i )
00628 {
00629
00630 i_Current = m_vi_Edges [i];
00631
00632
00633 for ( j=m_vi_LeftVertices [i_Current]; j<m_vi_LeftVertices [i_Current+1]; ++j )
00634 {
00635
00636
00637 if ( m_vi_Edges [j] == i_SelectedVertex || vi_Visited [m_vi_Edges [j]] == i_SelectedVertex )
00638 {
00639
00640 continue;
00641 }
00642
00643 u = m_vi_Edges[j];
00644
00645
00646 if ( vi_InducedVertexDegree [u] == _UNKNOWN )
00647 {
00648
00649 continue;
00650 }
00651
00652
00653 vi_Visited [u] = i_SelectedVertex;
00654
00655
00656
00657 if(vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].size() > 1)
00658 {
00659 l = vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].back();
00660
00661 vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]][vi_VertexLocation[u]] = l;
00662
00663
00664 vi_VertexLocation[l] = vi_VertexLocation[u];
00665 }
00666
00667 vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].pop_back();
00668
00669
00670 vi_InducedVertexDegree[u]--;
00671
00672
00673 vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].push_back(u);
00674
00675
00676 vi_VertexLocation[u] = vvi_GroupedInducedVertexDegree[vi_InducedVertexDegree[u]].size() - 1;
00677 }
00678 }
00679
00680
00681
00682 vi_InducedVertexDegree [i_SelectedVertex] = _UNKNOWN;
00683
00684 m_vi_OrderedVertices[i_VertexCountMinus1 - i_SelectedVertexCount] = i_SelectedVertex + i_LeftVertexCount;
00685
00686
00687 ++i_SelectedVertexCount;
00688 }
00689
00690
00691 vi_InducedVertexDegree.clear();
00692 vi_VertexLocation.clear();
00693 vvi_GroupedInducedVertexDegree.clear();
00694
00695 return(_TRUE);
00696 }
00697
00698
00699
00700
00701 int BipartiteGraphPartialOrdering::RowIncidenceDegreeOrdering()
00702 {
00703 if(CheckVertexOrdering("ROW_INCIDENCE_DEGREE"))
00704 {
00705 return(_TRUE);
00706 }
00707
00708 int i, j, k, l, u;
00709 int i_Current;
00710 int i_HighestDegreeVertex, m_i_MaximumVertexDegree;
00711 int i_VertexCount, i_VertexDegree, i_IncidenceVertexDegree;
00712 int i_SelectedVertex, i_SelectedVertexCount;
00713 vector<int> vi_IncidenceVertexDegree;
00714 vector<int> vi_Visited;
00715 vector< vector <int> > vvi_GroupedIncidenceVertexDegree;
00716 vector< int > vi_VertexLocation;
00717
00718
00719 i_SelectedVertex = _UNKNOWN;
00720 i_VertexCount = (int)m_vi_LeftVertices.size () - 1;
00721 vvi_GroupedIncidenceVertexDegree.clear();
00722 vvi_GroupedIncidenceVertexDegree.resize ( i_VertexCount );
00723 i_HighestDegreeVertex = _UNKNOWN;
00724 m_i_MaximumVertexDegree = _UNKNOWN;
00725 i_IncidenceVertexDegree = 0;
00726 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00727 m_vi_OrderedVertices.clear();
00728 m_vi_OrderedVertices.reserve(i_VertexCount);
00729 vi_VertexLocation.clear();
00730 vi_VertexLocation.reserve(i_VertexCount);
00731 vi_IncidenceVertexDegree.clear();
00732 vi_IncidenceVertexDegree.reserve(i_VertexCount);
00733
00734 for ( i=0; i<i_VertexCount; ++i )
00735 {
00736
00737
00738
00739 i_VertexDegree = 0;
00740
00741 for ( j=m_vi_LeftVertices [i]; j<m_vi_LeftVertices [i+1]; ++j )
00742 {
00743 i_Current = m_vi_Edges [j];
00744
00745 for ( k=m_vi_RightVertices [i_Current]; k<m_vi_RightVertices [i_Current+1]; ++k )
00746 {
00747
00748
00749 if ( m_vi_Edges [k] != i && vi_Visited [m_vi_Edges [k]] != i )
00750 {
00751 ++i_VertexDegree;
00752
00753 vi_Visited [m_vi_Edges [k]] = i;
00754 }
00755 }
00756 }
00757 vi_IncidenceVertexDegree.push_back ( i_IncidenceVertexDegree );
00758 vvi_GroupedIncidenceVertexDegree [i_IncidenceVertexDegree].push_back ( i );
00759 vi_VertexLocation.push_back ( vvi_GroupedIncidenceVertexDegree [i_IncidenceVertexDegree].size () -1 );
00760
00761 if ( m_i_MaximumVertexDegree < i_VertexDegree )
00762 {
00763 m_i_MaximumVertexDegree = i_VertexDegree;
00764 i_HighestDegreeVertex = i;
00765 }
00766 }
00767
00768
00769 i_SelectedVertexCount = 0;
00770 vi_Visited.clear ();
00771 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00772
00773 int iMax = 0;
00774
00775 while ( i_SelectedVertexCount < i_VertexCount )
00776 {
00777 if(iMax != m_i_MaximumVertexDegree && vvi_GroupedIncidenceVertexDegree[iMax + 1].size() != _FALSE)
00778 iMax++;
00779
00780 for ( i=iMax; i>=0; i-- )
00781 {
00782 if ( (int)vvi_GroupedIncidenceVertexDegree [i].size () != 0 )
00783 {
00784 i_SelectedVertex = vvi_GroupedIncidenceVertexDegree [i].back ();
00785
00786 vvi_GroupedIncidenceVertexDegree [i].pop_back();
00787 break;
00788 }
00789 }
00790
00791
00792
00793 for ( i=m_vi_LeftVertices [i_SelectedVertex]; i<m_vi_LeftVertices [i_SelectedVertex+1]; ++i )
00794 {
00795 i_Current = m_vi_Edges [i];
00796
00797 for ( j=m_vi_RightVertices [i_Current]; j<m_vi_RightVertices [i_Current+1]; ++j )
00798 {
00799 u = m_vi_Edges[j];
00800
00801
00802 if ( vi_IncidenceVertexDegree [u] == _UNKNOWN )
00803 {
00804
00805 continue;
00806 }
00807
00808
00809
00810 if ( u == i_SelectedVertex || vi_Visited [u] == i_SelectedVertex )
00811 {
00812
00813 continue;
00814 }
00815
00816
00817 if(vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]].size() > 1)
00818 {
00819 l = vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]].back();
00820
00821 vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]][vi_VertexLocation[u]] = l;
00822
00823 vi_VertexLocation[l] = vi_VertexLocation[u];
00824 }
00825
00826
00827 vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]].pop_back();
00828
00829
00830 vi_Visited [u] = i_SelectedVertex;
00831
00832
00833
00834 ++vi_IncidenceVertexDegree [u];
00835
00836 vvi_GroupedIncidenceVertexDegree [vi_IncidenceVertexDegree [u]].push_back ( u );
00837
00838 vi_VertexLocation [u] = vvi_GroupedIncidenceVertexDegree [vi_IncidenceVertexDegree [u]].size () -1;
00839 }
00840 }
00841
00842
00843
00844 vi_IncidenceVertexDegree [i_SelectedVertex] = _UNKNOWN;
00845
00846 m_vi_OrderedVertices.push_back ( i_SelectedVertex );
00847
00848 ++i_SelectedVertexCount;
00849 }
00850
00851
00852 vi_IncidenceVertexDegree.clear();
00853 vi_VertexLocation.clear();
00854 vvi_GroupedIncidenceVertexDegree.clear();
00855
00856 return(_TRUE);
00857
00858 }
00859
00860
00861
00862
00863 int BipartiteGraphPartialOrdering::ColumnIncidenceDegreeOrdering()
00864 {
00865 if(CheckVertexOrdering("COLUMN_INCIDENCE_DEGREE"))
00866 {
00867 return(_TRUE);
00868 }
00869
00870 int i, j, k, l, u;
00871 int i_Current;
00872 int i_HighestDegreeVertex, m_i_MaximumVertexDegree;
00873 int i_VertexCount, i_VertexDegree, i_IncidenceVertexDegree;
00874 int i_SelectedVertex, i_SelectedVertexCount;
00875 vector<int> vi_IncidenceVertexDegree;
00876 vector<int> vi_Visited;
00877 vector< vector<int> > vvi_GroupedIncidenceVertexDegree;
00878 vector< int > vi_VertexLocation;
00879
00880
00881 i_SelectedVertex = _UNKNOWN;
00882 i_VertexCount = (int)m_vi_RightVertices.size () - 1;
00883 vvi_GroupedIncidenceVertexDegree.resize ( i_VertexCount );
00884 i_HighestDegreeVertex = _UNKNOWN;
00885 m_i_MaximumVertexDegree = _UNKNOWN;
00886 i_IncidenceVertexDegree = 0;
00887 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00888 m_vi_OrderedVertices.clear();
00889
00890 int i_LeftVertexCount = STEP_DOWN((signed) m_vi_LeftVertices.size());
00891
00892
00893 for ( i=0; i<i_VertexCount; ++i )
00894 {
00895
00896
00897
00898 i_VertexDegree = 0;
00899
00900 for ( j=m_vi_RightVertices [i]; j<m_vi_RightVertices [i+1]; ++j )
00901 {
00902 i_Current = m_vi_Edges [j];
00903
00904 for ( k=m_vi_LeftVertices [i_Current]; k<m_vi_LeftVertices [i_Current+1]; ++k )
00905 {
00906
00907
00908 if ( m_vi_Edges [k] != i && vi_Visited [m_vi_Edges [k]] != i )
00909 {
00910 ++i_VertexDegree;
00911
00912 vi_Visited [m_vi_Edges [k]] = i;
00913 }
00914 }
00915 }
00916 vi_IncidenceVertexDegree.push_back ( i_IncidenceVertexDegree );
00917 vvi_GroupedIncidenceVertexDegree [i_IncidenceVertexDegree].push_back ( i );
00918 vi_VertexLocation.push_back ( vvi_GroupedIncidenceVertexDegree [i_IncidenceVertexDegree].size () - 1);
00919
00920 if ( m_i_MaximumVertexDegree < i_VertexDegree )
00921 {
00922 m_i_MaximumVertexDegree = i_VertexDegree;
00923 i_HighestDegreeVertex = i;
00924 }
00925 }
00926
00927
00928 i_SelectedVertexCount = 0;
00929 vi_Visited.clear ();
00930 vi_Visited.resize ( i_VertexCount, _UNKNOWN );
00931
00932 int iMax = 0;
00933
00934 while ( i_SelectedVertexCount < i_VertexCount )
00935 {
00936 if(iMax != m_i_MaximumVertexDegree && vvi_GroupedIncidenceVertexDegree[iMax + 1].size() != _FALSE)
00937 iMax++;
00938
00939
00940 for ( i=m_i_MaximumVertexDegree; i>=0; i-- )
00941 {
00942 if ( (int)vvi_GroupedIncidenceVertexDegree [i].size () != 0 )
00943 {
00944 i_SelectedVertex = vvi_GroupedIncidenceVertexDegree [i].back ();
00945
00946 vvi_GroupedIncidenceVertexDegree [i].pop_back();
00947 break;
00948 }
00949 }
00950
00951
00952
00953 for ( i=m_vi_RightVertices [i_SelectedVertex]; i<m_vi_RightVertices [i_SelectedVertex+1]; ++i )
00954 {
00955 i_Current = m_vi_Edges [i];
00956
00957 for ( j=m_vi_LeftVertices [i_Current]; j<m_vi_LeftVertices [i_Current+1]; ++j )
00958 {
00959 u = m_vi_Edges[j];
00960
00961
00962 if ( vi_IncidenceVertexDegree [u] == _UNKNOWN )
00963 {
00964
00965 continue;
00966 }
00967
00968
00969
00970 if ( u == i_SelectedVertex || vi_Visited [u] == i_SelectedVertex )
00971 {
00972
00973 continue;
00974 }
00975
00976
00977 if(vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]].size() > 1)
00978 {
00979 l = vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]].back();
00980
00981 vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]][vi_VertexLocation[u]] = l;
00982
00983 vi_VertexLocation[l] = vi_VertexLocation[u];
00984 }
00985
00986
00987 vvi_GroupedIncidenceVertexDegree[vi_IncidenceVertexDegree[u]].pop_back();
00988
00989
00990 vi_Visited [u] = i_SelectedVertex;
00991
00992
00993
00994 ++vi_IncidenceVertexDegree [u];
00995
00996 vvi_GroupedIncidenceVertexDegree [vi_IncidenceVertexDegree [u]].push_back ( u );
00997
00998 vi_VertexLocation [u] = vvi_GroupedIncidenceVertexDegree [vi_IncidenceVertexDegree [u]].size () -1;
00999 }
01000 }
01001
01002
01003
01004 vi_IncidenceVertexDegree [i_SelectedVertex] = _UNKNOWN;
01005
01006 m_vi_OrderedVertices.push_back ( i_SelectedVertex + i_LeftVertexCount );
01007
01008 ++i_SelectedVertexCount;
01009 }
01010
01011
01012 vi_IncidenceVertexDegree.clear();
01013 vi_VertexLocation.clear();
01014 vvi_GroupedIncidenceVertexDegree.clear();
01015
01016 return(_TRUE);
01017 }
01018
01019
01020
01021 string BipartiteGraphPartialOrdering::GetVertexOrderingVariant()
01022 {
01023 if(m_s_VertexOrderingVariant.compare("ROW_NATURAL") == 0)
01024 {
01025 return("Row Natural");
01026 }
01027 else
01028 if(m_s_VertexOrderingVariant.compare("COLUMN_NATURAL") == 0)
01029 {
01030 return("Column Natural");
01031 }
01032 else
01033 if(m_s_VertexOrderingVariant.compare("ROW_LARGEST_FIRST") == 0)
01034 {
01035 return("Row Largest First");
01036 }
01037 else
01038 if(m_s_VertexOrderingVariant.compare("COLUMN_LARGEST_FIRST") == 0)
01039 {
01040 return("Column Largest First");
01041 }
01042 else
01043 if(m_s_VertexOrderingVariant.compare("ROW_SMALLEST_LAST") == 0)
01044 {
01045 return("Row Smallest Last");
01046 }
01047 else
01048 if(m_s_VertexOrderingVariant.compare("COLUMN_SMALLEST_LAST") == 0)
01049 {
01050 return("Column Smallest Last");
01051 }
01052 else
01053 if(m_s_VertexOrderingVariant.compare("ROW_INCIDENCE_DEGREE") == 0)
01054 {
01055 return("Row Incidence Degree");
01056 }
01057 else
01058 if(m_s_VertexOrderingVariant.compare("COLUMN_INCIDENCE_DEGREE") == 0)
01059 {
01060 return("Column Incidence Degree");
01061 }
01062 else
01063 {
01064 return("Unknown");
01065 }
01066 }
01067
01068
01069 void BipartiteGraphPartialOrdering::GetOrderedVertices(vector<int> &output)
01070 {
01071 output = (m_vi_OrderedVertices);
01072 }
01073
01074 int BipartiteGraphPartialOrdering::OrderVertices(string s_OrderingVariant, string s_ColoringVariant) {
01075 s_ColoringVariant = toUpper(s_ColoringVariant);
01076 s_OrderingVariant = toUpper(s_OrderingVariant);
01077
01078 if(s_ColoringVariant == "ROW_PARTIAL_DISTANCE_TWO")
01079 {
01080 if((s_OrderingVariant.compare("NATURAL") == 0))
01081 {
01082 return(RowNaturalOrdering());
01083 }
01084 else
01085 if((s_OrderingVariant.compare("LARGEST_FIRST") == 0))
01086 {
01087 return(RowLargestFirstOrdering());
01088 }
01089 else
01090 if((s_OrderingVariant.compare("SMALLEST_LAST") == 0))
01091 {
01092 return(RowSmallestLastOrdering());
01093 }
01094 else
01095 if((s_OrderingVariant.compare("INCIDENCE_DEGREE") == 0))
01096 {
01097 return(RowIncidenceDegreeOrdering());
01098 }
01099 else
01100 if((s_OrderingVariant.compare("RANDOM") == 0))
01101 {
01102 return(RowRandomOrdering());
01103 }
01104 else
01105 {
01106 cerr<<endl;
01107 cerr<<"Unknown Ordering Method";
01108 cerr<<endl;
01109 }
01110 }
01111 else
01112 if(s_ColoringVariant == "COLUMN_PARTIAL_DISTANCE_TWO")
01113 {
01114 if((s_OrderingVariant.compare("NATURAL") == 0))
01115 {
01116 return(ColumnNaturalOrdering());
01117 }
01118 else
01119 if((s_OrderingVariant.compare("LARGEST_FIRST") == 0))
01120 {
01121 return(ColumnLargestFirstOrdering());
01122 }
01123 else
01124 if((s_OrderingVariant.compare("SMALLEST_LAST") == 0))
01125 {
01126 return(ColumnSmallestLastOrdering());
01127 }
01128 else
01129 if((s_OrderingVariant.compare("INCIDENCE_DEGREE") == 0))
01130 {
01131 return(ColumnIncidenceDegreeOrdering());
01132 }
01133 else
01134 if((s_OrderingVariant.compare("RANDOM") == 0))
01135 {
01136 return(ColumnRandomOrdering());
01137 }
01138 else
01139 {
01140 cerr<<endl;
01141 cerr<<"Unknown Ordering Method: "<<s_OrderingVariant;
01142 cerr<<endl;
01143 }
01144 }
01145 else
01146 {
01147 cerr<<endl;
01148 cerr<<"Invalid s_ColoringVariant = \""<<s_ColoringVariant<<"\", must be either \"COLUMN_PARTIAL_DISTANCE_TWO\" or \"ROW_PARTIAL_DISTANCE_TWO\".";
01149 cerr<<endl;
01150 }
01151
01152 return(_TRUE);
01153 }
01154
01155
01156 void BipartiteGraphPartialOrdering::PrintVertexOrdering() {
01157 cout<<"PrintVertexOrdering() "<<m_s_VertexOrderingVariant<<endl;
01158 for(unsigned int i=0; i<m_vi_OrderedVertices.size();i++) {
01159
01160 cout<<"\t["<<setw(5)<<i<<"] "<<setw(5)<<m_vi_OrderedVertices[i]<<endl;
01161 }
01162 cout<<endl;
01163 }
01164
01165 double BipartiteGraphPartialOrdering::GetVertexOrderingTime() {
01166 return m_d_OrderingTime;
01167 }
01168
01169 }