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 HessianRecovery::DirectRecover_RowCompressedFormat(GraphColoringInterface* g, double** dp2_CompressedMatrix, unsigned int ** uip2_HessianSparsityPattern, double*** dp3_HessianValue) {
00029 if(g==NULL) {
00030 cerr<<"g==NULL"<<endl;
00031 return _FALSE;
00032 }
00033
00034 if(AF_available) reset();
00035
00036 int rowCount = g->GetVertexCount();
00037 int colorCount = g->GetVertexColorCount();
00038
00039 vector<int> vi_VertexColors;
00040 g->GetVertexColors(vi_VertexColors);
00041
00042
00043 int** colorStatistic = new int*[rowCount];
00044
00045
00046 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00047 colorStatistic[i] = new int[colorCount];
00048 for(unsigned int j=0; j < (unsigned int)colorCount; j++) colorStatistic[i][j] = 0;
00049 }
00050
00051
00052 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00053 int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00054 for(unsigned int j=1; j <= (unsigned int)numOfNonZeros; j++) {
00055
00056
00057 colorStatistic[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]++;
00058 }
00059 }
00060
00061
00062 *dp3_HessianValue = new double*[rowCount];
00063 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00064 unsigned int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00065 (*dp3_HessianValue)[i] = new double[numOfNonZeros+1];
00066 (*dp3_HessianValue)[i][0] = numOfNonZeros;
00067 for(unsigned int j=1; j <= numOfNonZeros; j++) (*dp3_HessianValue)[i][j] = 0.;
00068 }
00069
00070
00071 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00072 unsigned int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00073 for(unsigned int j=1; j <= numOfNonZeros; j++) {
00074 if(i == uip2_HessianSparsityPattern[i][j]) {
00075 (*dp3_HessianValue)[i][j] = dp2_CompressedMatrix[i][vi_VertexColors[i]];
00076
00077 }
00078 else {
00079 if(colorStatistic[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]==1) {
00080 (*dp3_HessianValue)[i][j] = dp2_CompressedMatrix[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]];
00081
00082 }
00083 else {
00084 (*dp3_HessianValue)[i][j] = dp2_CompressedMatrix[uip2_HessianSparsityPattern[i][j]][vi_VertexColors[i]];
00085
00086 }
00087 }
00088 }
00089 }
00090
00091 free_2DMatrix(colorStatistic, rowCount);
00092 colorStatistic = NULL;
00093
00094 AF_available = true;
00095 i_AF_rowCount = rowCount;
00096 dp2_AF_Value = *dp3_HessianValue;
00097
00098 return (_TRUE);
00099 }
00100
00101
00102 int HessianRecovery::DirectRecover_SparseSolversFormat(GraphColoringInterface* g, double** dp2_CompressedMatrix, unsigned int ** uip2_HessianSparsityPattern, unsigned int** ip2_RowIndex, unsigned int** ip2_ColumnIndex, double** dp2_HessianValue) {
00103 if(g==NULL) {
00104 cerr<<"g==NULL"<<endl;
00105 return _FALSE;
00106 }
00107
00108 if(SSF_available) {
00109 cout<<"SSF_available="<<SSF_available<<endl; Pause();
00110 reset();
00111 }
00112
00113 int rowCount = g->GetVertexCount();
00114 int colorCount = g->GetVertexColorCount();
00115
00116 vector<int> vi_VertexColors;
00117 g->GetVertexColors(vi_VertexColors);
00118
00119
00120
00121 unsigned int numOfNonZerosInHessianValue = RowCompressedFormat_2_SparseSolversFormat_StructureOnly(uip2_HessianSparsityPattern, rowCount, ip2_RowIndex, ip2_ColumnIndex);
00122
00123
00124 int** colorStatistic = new int*[rowCount];
00125
00126
00127 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00128 colorStatistic[i] = new int[colorCount];
00129 for(unsigned int j=0; j < (unsigned int)colorCount; j++) colorStatistic[i][j] = 0;
00130 }
00131
00132
00133 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00134 int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00135 for(unsigned int j=1; j <= (unsigned int)numOfNonZeros; j++) {
00136
00137
00138 colorStatistic[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]++;
00139 }
00140 }
00141
00142
00143
00144 (*dp2_HessianValue) = new double[numOfNonZerosInHessianValue];
00145 for(unsigned int i=0; i < numOfNonZerosInHessianValue; i++) (*dp2_HessianValue)[i] = 0.;
00146
00147
00148
00149 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00150 unsigned int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00151 unsigned int offset = 0;
00152
00153 for(unsigned int j=1; j <= numOfNonZeros; j++) {
00154 if (i > uip2_HessianSparsityPattern[i][j]) {
00155 offset++;
00156 continue;
00157 }
00158 else if(i == uip2_HessianSparsityPattern[i][j]) {
00159 (*dp2_HessianValue)[(*ip2_RowIndex)[i] + j - offset - 1] = dp2_CompressedMatrix[i][vi_VertexColors[i]];
00160
00161
00162 }
00163 else {
00164 if(colorStatistic[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]==1) {
00165 (*dp2_HessianValue)[(*ip2_RowIndex)[i] + j - offset - 1] = dp2_CompressedMatrix[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]];
00166
00167 }
00168 else {
00169 (*dp2_HessianValue)[(*ip2_RowIndex)[i] + j - offset - 1] = dp2_CompressedMatrix[uip2_HessianSparsityPattern[i][j]][vi_VertexColors[i]];
00170
00171 }
00172 }
00173 }
00174 }
00175
00176 free_2DMatrix(colorStatistic, rowCount);
00177 colorStatistic = NULL;
00178
00179
00180
00181 for(unsigned int i=0; i <= (unsigned int) rowCount ; i++) {
00182 (*ip2_RowIndex)[i]++;
00183 }
00184 for(unsigned int i=0; i < numOfNonZerosInHessianValue; i++) {
00185 (*ip2_ColumnIndex)[i]++;
00186 }
00187
00188
00189 SSF_available = true;
00190 i_SSF_rowCount = rowCount;
00191 ip_SSF_RowIndex = *ip2_RowIndex;
00192 ip_SSF_ColumnIndex = *ip2_ColumnIndex;
00193 dp_SSF_Value = *dp2_HessianValue;
00194
00195 return (_TRUE);
00196 }
00197
00198
00199 int HessianRecovery::DirectRecover_CoordinateFormat(GraphColoringInterface* g, double** dp2_CompressedMatrix, unsigned int ** uip2_HessianSparsityPattern, unsigned int** ip2_RowIndex, unsigned int** ip2_ColumnIndex, double** dp2_HessianValue) {
00200 if(g==NULL) {
00201 cerr<<"g==NULL"<<endl;
00202 return _FALSE;
00203 }
00204
00205 if(CF_available) reset();
00206
00207 int rowCount = g->GetVertexCount();
00208 int colorCount = g->GetVertexColorCount();
00209
00210 vector<int> vi_VertexColors;
00211 g->GetVertexColors(vi_VertexColors);
00212
00213
00214 int** colorStatistic = new int*[rowCount];
00215
00216
00217 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00218 colorStatistic[i] = new int[colorCount];
00219 for(unsigned int j=0; j < (unsigned int)colorCount; j++) colorStatistic[i][j] = 0;
00220 }
00221
00222
00223 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00224 int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00225 for(unsigned int j=1; j <= (unsigned int)numOfNonZeros; j++) {
00226
00227
00228 colorStatistic[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]++;
00229 }
00230 }
00231
00232 vector<unsigned int> RowIndex;
00233 vector<unsigned int> ColumnIndex;
00234 vector<double> HessianValue;
00235
00236
00237 for(unsigned int i=0; i < (unsigned int)rowCount; i++) {
00238 unsigned int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00239 for(unsigned int j=1; j <= numOfNonZeros; j++) {
00240 if(uip2_HessianSparsityPattern[i][j]<i) continue;
00241
00242 if(i == uip2_HessianSparsityPattern[i][j]) {
00243 HessianValue.push_back(dp2_CompressedMatrix[i][vi_VertexColors[i]]);
00244 }
00245 else {
00246 if(colorStatistic[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]==1) {
00247 HessianValue.push_back(dp2_CompressedMatrix[i][vi_VertexColors[uip2_HessianSparsityPattern[i][j]]]);
00248 }
00249 else {
00250 HessianValue.push_back(dp2_CompressedMatrix[uip2_HessianSparsityPattern[i][j]][vi_VertexColors[i]]);
00251 }
00252 }
00253 RowIndex.push_back(i);
00254 ColumnIndex.push_back(uip2_HessianSparsityPattern[i][j]);
00255 }
00256 }
00257
00258 unsigned int numOfNonZeros = RowIndex.size();
00259 (*ip2_RowIndex) = new unsigned int[numOfNonZeros];
00260 (*ip2_ColumnIndex) = new unsigned int[numOfNonZeros];
00261 (*dp2_HessianValue) = new double[numOfNonZeros];
00262
00263 for(int i=0; i < numOfNonZeros; i++) {
00264 (*ip2_RowIndex)[i] = RowIndex[i];
00265 (*ip2_ColumnIndex)[i] = ColumnIndex[i];
00266 (*dp2_HessianValue)[i] = HessianValue[i];
00267 }
00268
00269 free_2DMatrix(colorStatistic, rowCount);
00270 colorStatistic = NULL;
00271
00272
00273 CF_available = true;
00274 i_CF_rowCount = rowCount;
00275 ip_CF_RowIndex = *ip2_RowIndex;
00276 ip_CF_ColumnIndex = *ip2_ColumnIndex;
00277 dp_CF_Value = *dp2_HessianValue;
00278
00279 return (numOfNonZeros);
00280 }
00281
00282 int HessianRecovery::IndirectRecover_RowCompressedFormat(GraphColoringInterface* g, double** dp2_CompressedMatrix, unsigned int ** uip2_HessianSparsityPattern, double*** dp3_HessianValue) {
00283 if(g==NULL) {
00284 cerr<<"g==NULL"<<endl;
00285 return _FALSE;
00286 }
00287
00288 if(AF_available) reset();
00289
00290 int i=0,j=0;
00291 int i_VertexCount = _UNKNOWN;
00292 int i_EdgeID, i_SetID;
00293 vector<int> vi_Sets;
00294 map< int, vector<int> > mivi_VertexSets;
00295
00296 vector<int> vi_Vertices;
00297 g->GetVertices(vi_Vertices);
00298
00299 vector<int> vi_Edges;
00300 g->GetEdges(vi_Edges);
00301
00302 vector<int> vi_VertexColors;
00303 g->GetVertexColors(vi_VertexColors);
00304
00305 map< int, map< int, int> > mimi2_VertexEdgeMap;
00306 g->GetVertexEdgeMap(mimi2_VertexEdgeMap);
00307
00308 DisjointSets ds_DisjointSets;
00309 g->GetDisjointSets(ds_DisjointSets);
00310
00311
00312 vi_Sets.clear();
00313 mivi_VertexSets.clear();
00314
00315 i_VertexCount = g->GetVertexCount();
00316
00317 for(i=0; i<i_VertexCount; i++)
00318 {
00319 for(j=vi_Vertices[i]; j<vi_Vertices[STEP_UP(i)]; j++)
00320 {
00321 if(i < vi_Edges[j])
00322
00323 {
00324 i_EdgeID = mimi2_VertexEdgeMap[i][vi_Edges[j]];
00325
00326 i_SetID = ds_DisjointSets.FindAndCompress(i_EdgeID);
00327
00328 if(i_SetID == i_EdgeID)
00329 {
00330 vi_Sets.push_back(i_SetID);
00331 }
00332
00333 mivi_VertexSets[i_SetID].push_back(i);
00334 mivi_VertexSets[i_SetID].push_back(vi_Edges[j]);
00335 }
00336 }
00337 }
00338
00339 int i_MaximumVertexDegree;
00340
00341 int i_HighestInducedVertexDegree;
00342
00343 int i_LeafVertex, i_ParentVertex, i_PresentVertex;
00344
00345 int i_VertexDegree;
00346
00347 int i_SetCount, i_SetSize;
00348
00349
00350
00351 double d_Value;
00352
00353 vector<int> vi_EvaluatedDiagonals;
00354
00355 vector<int> vi_InducedVertexDegrees;
00356
00357 vector<double> vd_IncludedVertices;
00358
00359 vector< vector<int> > v2i_VertexAdjacency;
00360
00361 vector< vector<double> > v2d_NonzeroAdjacency;
00362
00363 vector< list<int> > vli_GroupedInducedVertexDegrees;
00364
00365 vector< list<int>::iterator > vlit_VertexLocations;
00366
00367 i_MaximumVertexDegree = g->GetMaximumVertexDegree();
00368
00369 #if DEBUG == 5103
00370
00371 cout<<endl;
00372 cout<<"DEBUG 5103 | Hessian Evaluation | Bicolored Sets"<<endl;
00373 cout<<endl;
00374
00375 i_SetCount = (signed) vi_Sets.size();
00376
00377 for(i=0; i<i_SetCount; i++)
00378 {
00379 cout<<STEP_UP(vi_Sets[i])<<"\t"<<" : ";
00380
00381 i_SetSize = (signed) mivi_VertexSets[vi_Sets[i]].size();
00382
00383 for(j=0; j<i_SetSize; j++)
00384 {
00385 if(j == STEP_DOWN(i_SetSize))
00386 {
00387 cout<<STEP_UP(mivi_VertexSets[vi_Sets[i]][j])<<" ("<<i_SetSize<<")"<<endl;
00388 }
00389 else
00390 {
00391 cout<<STEP_UP(mivi_VertexSets[vi_Sets[i]][j])<<", ";
00392 }
00393 }
00394 }
00395
00396 cout<<endl;
00397 cout<<"[Set Count = "<<i_SetCount<<"]"<<endl;
00398 cout<<endl;
00399
00400 #endif
00401
00402
00403 i_VertexCount = g->GetVertexCount();
00404
00405 v2i_VertexAdjacency.clear();
00406 v2i_VertexAdjacency.resize((unsigned) i_VertexCount);
00407
00408 v2d_NonzeroAdjacency.clear();
00409 v2d_NonzeroAdjacency.resize((unsigned) i_VertexCount);
00410
00411 vi_EvaluatedDiagonals.clear();
00412 vi_EvaluatedDiagonals.resize((unsigned) i_VertexCount, _FALSE);
00413
00414 vi_InducedVertexDegrees.clear();
00415 vi_InducedVertexDegrees.resize((unsigned) i_VertexCount, _FALSE);
00416
00417 vd_IncludedVertices.clear();
00418 vd_IncludedVertices.resize((unsigned) i_VertexCount, _UNKNOWN);
00419
00420 i_ParentVertex = _UNKNOWN;
00421
00422 i_SetCount = (signed) vi_Sets.size();
00423
00424 for(i=0; i<i_SetCount; i++)
00425 {
00426 vli_GroupedInducedVertexDegrees.clear();
00427 vli_GroupedInducedVertexDegrees.resize((unsigned) STEP_UP(i_MaximumVertexDegree));
00428
00429 vlit_VertexLocations.clear();
00430 vlit_VertexLocations.resize((unsigned) i_VertexCount);
00431
00432 i_HighestInducedVertexDegree = _UNKNOWN;
00433
00434 i_SetSize = (signed) mivi_VertexSets[vi_Sets[i]].size();
00435
00436 for(j=0; j<i_SetSize; j++)
00437 {
00438 i_PresentVertex = mivi_VertexSets[vi_Sets[i]][j];
00439
00440 vd_IncludedVertices[i_PresentVertex] = _FALSE;
00441
00442 if(vi_InducedVertexDegrees[i_PresentVertex] != _FALSE)
00443 {
00444 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].erase(vlit_VertexLocations[i_PresentVertex]);
00445 }
00446
00447 vi_InducedVertexDegrees[i_PresentVertex]++;
00448
00449 if(i_HighestInducedVertexDegree < vi_InducedVertexDegrees[i_PresentVertex])
00450 {
00451 i_HighestInducedVertexDegree = vi_InducedVertexDegrees[i_PresentVertex];
00452 }
00453
00454 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].push_front(i_PresentVertex);
00455
00456 vlit_VertexLocations[i_PresentVertex] = vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].begin();
00457 }
00458
00459 #if DEBUG == 5103
00460
00461 int k;
00462
00463 list<int>::iterator lit_ListIterator;
00464
00465 cout<<endl;
00466 cout<<"DEBUG 5103 | Hessian Evaluation | Induced Vertex Degrees | Set "<<STEP_UP(i)<<endl;
00467 cout<<endl;
00468
00469 for(j=0; j<STEP_UP(i_HighestInducedVertexDegree); j++)
00470 {
00471 i_SetSize = (signed) vli_GroupedInducedVertexDegrees[j].size();
00472
00473 if(i_SetSize == _FALSE)
00474 {
00475 continue;
00476 }
00477
00478 k = _FALSE;
00479
00480 cout<<"Degree "<<j<<"\t"<<" : ";
00481
00482 for(lit_ListIterator=vli_GroupedInducedVertexDegrees[j].begin(); lit_ListIterator!=vli_GroupedInducedVertexDegrees[j].end(); lit_ListIterator++)
00483 {
00484 if(k == STEP_DOWN(i_SetSize))
00485 {
00486 cout<<STEP_UP(*lit_ListIterator)<<" ("<<i_SetSize<<")"<<endl;
00487 }
00488 else
00489 {
00490 cout<<STEP_UP(*lit_ListIterator)<<", ";
00491 }
00492
00493 k++;
00494 }
00495 }
00496
00497 #endif
00498
00499 #if DEBUG == 5103
00500
00501 cout<<endl;
00502 cout<<"DEBUG 5103 | Hessian Evaluation | Retrieved Elements"<<"| Set "<<STEP_UP(i)<<endl;
00503 cout<<endl;
00504
00505 #endif
00506
00507
00508 for (int index = 0; index < i_VertexCount; index++) {
00509 if(vi_EvaluatedDiagonals[index] == _FALSE)
00510 {
00511 d_Value = dp2_CompressedMatrix[index][vi_VertexColors[index]];
00512
00513 #if DEBUG == 5103
00514
00515 cout<<"Element["<<STEP_UP(index)<<"]["<<STEP_UP(index)<<"] = "<<d_Value<<endl;
00516
00517 #endif
00518 v2i_VertexAdjacency[index].push_back(index);
00519 v2d_NonzeroAdjacency[index].push_back(d_Value);
00520
00521 vi_EvaluatedDiagonals[index] = _TRUE;
00522
00523 }
00524 }
00525
00526 for ( ; ; )
00527 {
00528 if(vli_GroupedInducedVertexDegrees[_TRUE].empty())
00529 {
00530 i_LeafVertex = vli_GroupedInducedVertexDegrees[_FALSE].front();
00531
00532 vi_InducedVertexDegrees[i_LeafVertex] = _FALSE;
00533
00534 vd_IncludedVertices[i_LeafVertex] = _UNKNOWN;
00535
00536 break;
00537 }
00538
00539 i_LeafVertex = vli_GroupedInducedVertexDegrees[_TRUE].front();
00540
00541 vli_GroupedInducedVertexDegrees[_TRUE].pop_front();
00542
00543
00544
00545 for(j=vi_Vertices[i_LeafVertex]; j<vi_Vertices[STEP_UP(i_LeafVertex)]; j++)
00546 {
00547 if(vd_IncludedVertices[vi_Edges[j]] != _UNKNOWN)
00548 {
00549 i_ParentVertex = vi_Edges[j];
00550
00551 break;
00552 }
00553 }
00554
00555 d_Value = dp2_CompressedMatrix[i_LeafVertex][vi_VertexColors[i_ParentVertex]] - vd_IncludedVertices[i_LeafVertex];
00556
00557 vd_IncludedVertices[i_ParentVertex] += d_Value;
00558
00559 vi_InducedVertexDegrees[i_LeafVertex] = _FALSE;
00560 vd_IncludedVertices[i_LeafVertex] = _UNKNOWN;
00561 if(vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].size()>1) {
00562 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].erase(vlit_VertexLocations[i_ParentVertex]);
00563 }
00564 else {
00565 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].pop_back();
00566 }
00567
00568 vi_InducedVertexDegrees[i_ParentVertex]--;
00569 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].push_back(i_ParentVertex);
00570
00571
00572 vlit_VertexLocations[i_ParentVertex] = vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].end();
00573 --vlit_VertexLocations[i_ParentVertex];
00574
00575 v2i_VertexAdjacency[i_LeafVertex].push_back(i_ParentVertex);
00576 v2d_NonzeroAdjacency[i_LeafVertex].push_back(d_Value);
00577
00578 v2i_VertexAdjacency[i_ParentVertex].push_back(i_LeafVertex);
00579 v2d_NonzeroAdjacency[i_ParentVertex].push_back(d_Value);
00580
00581
00582 #if DEBUG == 5103
00583
00584 cout<<"Element["<<STEP_UP(i_LeafVertex)<<"]["<<STEP_UP(i_ParentVertex)<<"] = "<<d_Value<<endl;
00585 #endif
00586
00587 }
00588 }
00589
00590
00591
00592 *dp3_HessianValue = new double*[i_VertexCount];
00593 for(unsigned int i=0; i < (unsigned int)i_VertexCount; i++) {
00594 unsigned int numOfNonZeros = uip2_HessianSparsityPattern[i][0];
00595 (*dp3_HessianValue)[i] = new double[numOfNonZeros+1];
00596 (*dp3_HessianValue)[i][0] = numOfNonZeros;
00597 for(unsigned int j=1; j <= numOfNonZeros; j++) (*dp3_HessianValue)[i][j] = 0.;
00598 }
00599
00600
00601 for(i=0; i<i_VertexCount; i++) {
00602 int NumOfNonzeros = uip2_HessianSparsityPattern[i][0];
00603 i_VertexDegree = (signed) v2i_VertexAdjacency[i].size();
00604 for(j=1; j<=NumOfNonzeros; j++) {
00605 int targetColumnID = uip2_HessianSparsityPattern[i][j];
00606 for (int k=0; k<i_VertexDegree; k++) {
00607 if(targetColumnID == v2i_VertexAdjacency[i][k]) {
00608 (*dp3_HessianValue)[i][j] = v2d_NonzeroAdjacency[i][k];
00609 break;
00610 }
00611 }
00612 }
00613 }
00614
00615 #undef DEBUG
00616
00617 AF_available = true;
00618 i_AF_rowCount = i_VertexCount;
00619 dp2_AF_Value = *dp3_HessianValue;
00620
00621 return (_TRUE);
00622 }
00623
00624
00625 int HessianRecovery::IndirectRecover_SparseSolversFormat(GraphColoringInterface* g, double** dp2_CompressedMatrix, unsigned int ** uip2_HessianSparsityPattern, unsigned int** ip2_RowIndex, unsigned int** ip2_ColumnIndex, double** dp2_HessianValue) {
00626 if(g==NULL) {
00627 cerr<<"g==NULL"<<endl;
00628 return _FALSE;
00629 }
00630
00631
00632 if(SSF_available) {
00633 cout<<"SSF_available="<<SSF_available<<endl; Pause();
00634 reset();
00635 }
00636
00637 unsigned int numOfNonZerosInHessianValue = RowCompressedFormat_2_SparseSolversFormat_StructureOnly(uip2_HessianSparsityPattern, g->GetVertexCount(), ip2_RowIndex, ip2_ColumnIndex);
00638
00639 int i=0,j=0;
00640 int i_VertexCount = _UNKNOWN;
00641 int i_EdgeID, i_SetID;
00642 vector<int> vi_Sets;
00643 map< int, vector<int> > mivi_VertexSets;
00644
00645 vector<int> vi_Vertices;
00646 g->GetVertices(vi_Vertices);
00647
00648 vector<int> vi_Edges;
00649 g->GetEdges(vi_Edges);
00650
00651 vector<int> vi_VertexColors;
00652 g->GetVertexColors(vi_VertexColors);
00653
00654 map< int, map< int, int> > mimi2_VertexEdgeMap;
00655 g->GetVertexEdgeMap(mimi2_VertexEdgeMap);
00656
00657 DisjointSets ds_DisjointSets;
00658 g->GetDisjointSets(ds_DisjointSets);
00659
00660
00661 vi_Sets.clear();
00662 mivi_VertexSets.clear();
00663
00664 i_VertexCount = g->GetVertexCount();
00665
00666 for(i=0; i<i_VertexCount; i++)
00667 {
00668 for(j=vi_Vertices[i]; j<vi_Vertices[STEP_UP(i)]; j++)
00669 {
00670 if(i < vi_Edges[j])
00671
00672 {
00673 i_EdgeID = mimi2_VertexEdgeMap[i][vi_Edges[j]];
00674
00675 i_SetID = ds_DisjointSets.FindAndCompress(i_EdgeID);
00676
00677 if(i_SetID == i_EdgeID)
00678 {
00679 vi_Sets.push_back(i_SetID);
00680 }
00681
00682 mivi_VertexSets[i_SetID].push_back(i);
00683 mivi_VertexSets[i_SetID].push_back(vi_Edges[j]);
00684 }
00685 }
00686 }
00687
00688 int i_MaximumVertexDegree;
00689
00690 int i_HighestInducedVertexDegree;
00691
00692 int i_LeafVertex, i_ParentVertex, i_PresentVertex;
00693
00694 int i_VertexDegree;
00695
00696 int i_SetCount, i_SetSize;
00697
00698 double d_Value;
00699
00700 vector<int> vi_EvaluatedDiagonals;
00701
00702 vector<int> vi_InducedVertexDegrees;
00703
00704 vector<double> vd_IncludedVertices;
00705
00706 vector< vector<int> > v2i_VertexAdjacency;
00707
00708 vector< vector<double> > v2d_NonzeroAdjacency;
00709
00710 vector< list<int> > vli_GroupedInducedVertexDegrees;
00711
00712 vector< list<int>::iterator > vlit_VertexLocations;
00713
00714 i_MaximumVertexDegree = g->GetMaximumVertexDegree();
00715
00716 #if DEBUG == 5103
00717
00718 cout<<endl;
00719 cout<<"DEBUG 5103 | Hessian Evaluation | Bicolored Sets"<<endl;
00720 cout<<endl;
00721
00722 i_SetCount = (signed) vi_Sets.size();
00723
00724 for(i=0; i<i_SetCount; i++)
00725 {
00726 cout<<STEP_UP(vi_Sets[i])<<"\t"<<" : ";
00727
00728 i_SetSize = (signed) mivi_VertexSets[vi_Sets[i]].size();
00729
00730 for(j=0; j<i_SetSize; j++)
00731 {
00732 if(j == STEP_DOWN(i_SetSize))
00733 {
00734 cout<<STEP_UP(mivi_VertexSets[vi_Sets[i]][j])<<" ("<<i_SetSize<<")"<<endl;
00735 }
00736 else
00737 {
00738 cout<<STEP_UP(mivi_VertexSets[vi_Sets[i]][j])<<", ";
00739 }
00740 }
00741 }
00742
00743 cout<<endl;
00744 cout<<"[Set Count = "<<i_SetCount<<"]"<<endl;
00745 cout<<endl;
00746
00747 #endif
00748
00749
00750 i_VertexCount = g->GetVertexCount();
00751
00752 v2i_VertexAdjacency.clear();
00753 v2i_VertexAdjacency.resize((unsigned) i_VertexCount);
00754
00755 v2d_NonzeroAdjacency.clear();
00756 v2d_NonzeroAdjacency.resize((unsigned) i_VertexCount);
00757
00758 vi_EvaluatedDiagonals.clear();
00759 vi_EvaluatedDiagonals.resize((unsigned) i_VertexCount, _FALSE);
00760
00761 vi_InducedVertexDegrees.clear();
00762 vi_InducedVertexDegrees.resize((unsigned) i_VertexCount, _FALSE);
00763
00764 vd_IncludedVertices.clear();
00765 vd_IncludedVertices.resize((unsigned) i_VertexCount, _UNKNOWN);
00766
00767 i_ParentVertex = _UNKNOWN;
00768
00769 i_SetCount = (signed) vi_Sets.size();
00770
00771 for(i=0; i<i_SetCount; i++)
00772 {
00773 vli_GroupedInducedVertexDegrees.clear();
00774 vli_GroupedInducedVertexDegrees.resize((unsigned) STEP_UP(i_MaximumVertexDegree));
00775
00776 vlit_VertexLocations.clear();
00777 vlit_VertexLocations.resize((unsigned) i_VertexCount);
00778
00779 i_HighestInducedVertexDegree = _UNKNOWN;
00780
00781 i_SetSize = (signed) mivi_VertexSets[vi_Sets[i]].size();
00782
00783 for(j=0; j<i_SetSize; j++)
00784 {
00785 i_PresentVertex = mivi_VertexSets[vi_Sets[i]][j];
00786
00787 vd_IncludedVertices[i_PresentVertex] = _FALSE;
00788
00789 if(vi_InducedVertexDegrees[i_PresentVertex] != _FALSE)
00790 {
00791 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].erase(vlit_VertexLocations[i_PresentVertex]);
00792 }
00793
00794 vi_InducedVertexDegrees[i_PresentVertex]++;
00795
00796 if(i_HighestInducedVertexDegree < vi_InducedVertexDegrees[i_PresentVertex])
00797 {
00798 i_HighestInducedVertexDegree = vi_InducedVertexDegrees[i_PresentVertex];
00799 }
00800
00801 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].push_front(i_PresentVertex);
00802
00803 vlit_VertexLocations[i_PresentVertex] = vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].begin();
00804 }
00805
00806 #if DEBUG == 5103
00807
00808 int k;
00809
00810 list<int>::iterator lit_ListIterator;
00811
00812 cout<<endl;
00813 cout<<"DEBUG 5103 | Hessian Evaluation | Induced Vertex Degrees | Set "<<STEP_UP(i)<<endl;
00814 cout<<endl;
00815
00816 for(j=0; j<STEP_UP(i_HighestInducedVertexDegree); j++)
00817 {
00818 i_SetSize = (signed) vli_GroupedInducedVertexDegrees[j].size();
00819
00820 if(i_SetSize == _FALSE)
00821 {
00822 continue;
00823 }
00824
00825 k = _FALSE;
00826
00827 cout<<"Degree "<<j<<"\t"<<" : ";
00828
00829 for(lit_ListIterator=vli_GroupedInducedVertexDegrees[j].begin(); lit_ListIterator!=vli_GroupedInducedVertexDegrees[j].end(); lit_ListIterator++)
00830 {
00831 if(k == STEP_DOWN(i_SetSize))
00832 {
00833 cout<<STEP_UP(*lit_ListIterator)<<" ("<<i_SetSize<<")"<<endl;
00834 }
00835 else
00836 {
00837 cout<<STEP_UP(*lit_ListIterator)<<", ";
00838 }
00839
00840 k++;
00841 }
00842 }
00843
00844 #endif
00845
00846 #if DEBUG == 5103
00847
00848 cout<<endl;
00849 cout<<"DEBUG 5103 | Hessian Evaluation | Retrieved Elements"<<"| Set "<<STEP_UP(i)<<endl;
00850 cout<<endl;
00851
00852 #endif
00853
00854
00855 for (int index = 0; index < i_VertexCount; index++) {
00856 if(vi_EvaluatedDiagonals[index] == _FALSE)
00857 {
00858 d_Value = dp2_CompressedMatrix[index][vi_VertexColors[index]];
00859
00860 #if DEBUG == 5103
00861
00862 cout<<"Element["<<STEP_UP(index)<<"]["<<STEP_UP(index)<<"] = "<<d_Value<<endl;
00863
00864 #endif
00865 v2i_VertexAdjacency[index].push_back(index);
00866 v2d_NonzeroAdjacency[index].push_back(d_Value);
00867
00868 vi_EvaluatedDiagonals[index] = _TRUE;
00869
00870 }
00871 }
00872
00873 for ( ; ; )
00874 {
00875 if(vli_GroupedInducedVertexDegrees[_TRUE].empty())
00876 {
00877 i_LeafVertex = vli_GroupedInducedVertexDegrees[_FALSE].front();
00878
00879 vi_InducedVertexDegrees[i_LeafVertex] = _FALSE;
00880
00881 vd_IncludedVertices[i_LeafVertex] = _UNKNOWN;
00882
00883 break;
00884 }
00885
00886 i_LeafVertex = vli_GroupedInducedVertexDegrees[_TRUE].front();
00887
00888 vli_GroupedInducedVertexDegrees[_TRUE].pop_front();
00889
00890
00891 for(j=vi_Vertices[i_LeafVertex]; j<vi_Vertices[STEP_UP(i_LeafVertex)]; j++)
00892 {
00893 if(vd_IncludedVertices[vi_Edges[j]] != _UNKNOWN)
00894 {
00895 i_ParentVertex = vi_Edges[j];
00896
00897 break;
00898 }
00899 }
00900
00901 d_Value = dp2_CompressedMatrix[i_LeafVertex][vi_VertexColors[i_ParentVertex]] - vd_IncludedVertices[i_LeafVertex];
00902
00903 vd_IncludedVertices[i_ParentVertex] += d_Value;
00904
00905 vi_InducedVertexDegrees[i_LeafVertex] = _FALSE;
00906 vd_IncludedVertices[i_LeafVertex] = _UNKNOWN;
00907 if(vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].size()>1) {
00908 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].erase(vlit_VertexLocations[i_ParentVertex]);
00909 }
00910 else {
00911 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].pop_back();
00912 }
00913
00914 vi_InducedVertexDegrees[i_ParentVertex]--;
00915 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].push_back(i_ParentVertex);
00916
00917
00918 vlit_VertexLocations[i_ParentVertex] = vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].end();
00919 --vlit_VertexLocations[i_ParentVertex];
00920
00921 v2i_VertexAdjacency[i_LeafVertex].push_back(i_ParentVertex);
00922 v2d_NonzeroAdjacency[i_LeafVertex].push_back(d_Value);
00923
00924 v2i_VertexAdjacency[i_ParentVertex].push_back(i_LeafVertex);
00925 v2d_NonzeroAdjacency[i_ParentVertex].push_back(d_Value);
00926
00927 #if DEBUG == 5103
00928
00929 cout<<"Element["<<STEP_UP(i_LeafVertex)<<"]["<<STEP_UP(i_ParentVertex)<<"] = "<<d_Value<<endl;
00930 #endif
00931
00932 }
00933 }
00934
00935
00936
00937
00938
00939 (*dp2_HessianValue) = new double[numOfNonZerosInHessianValue];
00940 for(unsigned int i=0; i < numOfNonZerosInHessianValue; i++) (*dp2_HessianValue)[i] = 0.;
00941
00942
00943 for(i=0; i<i_VertexCount; i++) {
00944 int NumOfNonzeros = uip2_HessianSparsityPattern[i][0];
00945 int offset = 0;
00946 i_VertexDegree = (signed) v2i_VertexAdjacency[i].size();
00947 for(j=1; j<=NumOfNonzeros; j++) {
00948 if( i > uip2_HessianSparsityPattern[i][j] ) {
00949 offset++;
00950 continue;
00951 }
00952 int targetColumnID = uip2_HessianSparsityPattern[i][j];
00953 for (int k=0; k<i_VertexDegree; k++) {
00954 if(targetColumnID == v2i_VertexAdjacency[i][k]) {
00955 (*dp2_HessianValue)[(*ip2_RowIndex)[i] + j - offset - 1] = v2d_NonzeroAdjacency[i][k];
00956 break;
00957 }
00958 }
00959 }
00960 }
00961
00962 #undef DEBUG
00963
00964
00965
00966 for(unsigned int i=0; i <= (unsigned int) i_VertexCount ; i++) {
00967 (*ip2_RowIndex)[i]++;
00968 }
00969 for(unsigned int i=0; i < numOfNonZerosInHessianValue; i++) {
00970 (*ip2_ColumnIndex)[i]++;
00971 }
00972
00973 SSF_available = true;
00974 i_SSF_rowCount = i_VertexCount;
00975 ip_SSF_RowIndex = *ip2_RowIndex;
00976 ip_SSF_ColumnIndex = *ip2_ColumnIndex;
00977 dp_SSF_Value = *dp2_HessianValue;
00978
00979 return (_TRUE);
00980 }
00981
00982 int HessianRecovery::IndirectRecover_CoordinateFormat(GraphColoringInterface* g, double** dp2_CompressedMatrix, unsigned int ** uip2_HessianSparsityPattern, unsigned int** ip2_RowIndex, unsigned int** ip2_ColumnIndex, double** dp2_HessianValue) {
00983
00984 if(g==NULL) {
00985 cerr<<"g==NULL"<<endl;
00986 return _FALSE;
00987 }
00988
00989 if(CF_available) reset();
00990
00991 int i=0,j=0;
00992 int i_VertexCount = _UNKNOWN;
00993 int i_EdgeID, i_SetID;
00994 vector<int> vi_Sets;
00995 map< int, vector<int> > mivi_VertexSets;
00996 vector<int> vi_Vertices;
00997 g->GetVertices(vi_Vertices);
00998 vector<int> vi_Edges;
00999 g->GetEdges(vi_Edges);
01000 vector<int> vi_VertexColors;
01001 g->GetVertexColors(vi_VertexColors);
01002 map< int, map< int, int> > mimi2_VertexEdgeMap;
01003 g->GetVertexEdgeMap(mimi2_VertexEdgeMap);
01004 DisjointSets ds_DisjointSets;
01005 g->GetDisjointSets(ds_DisjointSets);
01006
01007
01008 vi_Sets.clear();
01009 mivi_VertexSets.clear();
01010
01011 i_VertexCount = g->GetVertexCount();
01012
01013 for(i=0; i<i_VertexCount; i++)
01014 {
01015 for(j=vi_Vertices[i]; j<vi_Vertices[STEP_UP(i)]; j++)
01016 {
01017 if(i < vi_Edges[j])
01018
01019 {
01020 i_EdgeID = mimi2_VertexEdgeMap[i][vi_Edges[j]];
01021
01022 i_SetID = ds_DisjointSets.FindAndCompress(i_EdgeID);
01023
01024 if(i_SetID == i_EdgeID)
01025 {
01026 vi_Sets.push_back(i_SetID);
01027 }
01028
01029 mivi_VertexSets[i_SetID].push_back(i);
01030 mivi_VertexSets[i_SetID].push_back(vi_Edges[j]);
01031 }
01032 }
01033 }
01034
01035 int i_MaximumVertexDegree;
01036
01037 int i_HighestInducedVertexDegree;
01038
01039 int i_LeafVertex, i_ParentVertex, i_PresentVertex;
01040
01041 int i_VertexDegree;
01042
01043 int i_SetCount, i_SetSize;
01044
01045 double d_Value;
01046
01047 vector<int> vi_EvaluatedDiagonals;
01048
01049 vector<int> vi_InducedVertexDegrees;
01050
01051 vector<double> vd_IncludedVertices;
01052
01053 vector< vector<int> > v2i_VertexAdjacency;
01054
01055 vector< vector<double> > v2d_NonzeroAdjacency;
01056
01057 vector< list<int> > vli_GroupedInducedVertexDegrees;
01058
01059 vector< list<int>::iterator > vlit_VertexLocations;
01060
01061 i_MaximumVertexDegree = g->GetMaximumVertexDegree();
01062
01063 #if DEBUG == 5103
01064
01065 cout<<endl;
01066 cout<<"DEBUG 5103 | Hessian Evaluation | Bicolored Sets"<<endl;
01067 cout<<endl;
01068
01069 i_SetCount = (signed) vi_Sets.size();
01070
01071 for(i=0; i<i_SetCount; i++)
01072 {
01073 cout<<STEP_UP(vi_Sets[i])<<"\t"<<" : ";
01074
01075 i_SetSize = (signed) mivi_VertexSets[vi_Sets[i]].size();
01076
01077 for(j=0; j<i_SetSize; j++)
01078 {
01079 if(j == STEP_DOWN(i_SetSize))
01080 {
01081 cout<<STEP_UP(mivi_VertexSets[vi_Sets[i]][j])<<" ("<<i_SetSize<<")"<<endl;
01082 }
01083 else
01084 {
01085 cout<<STEP_UP(mivi_VertexSets[vi_Sets[i]][j])<<", ";
01086 }
01087 }
01088 }
01089
01090 cout<<endl;
01091 cout<<"[Set Count = "<<i_SetCount<<"]"<<endl;
01092 cout<<endl;
01093
01094 #endif
01095
01096
01097 i_VertexCount = g->GetVertexCount();
01098
01099 v2i_VertexAdjacency.clear();
01100 v2i_VertexAdjacency.resize((unsigned) i_VertexCount);
01101
01102 v2d_NonzeroAdjacency.clear();
01103 v2d_NonzeroAdjacency.resize((unsigned) i_VertexCount);
01104
01105 vi_EvaluatedDiagonals.clear();
01106 vi_EvaluatedDiagonals.resize((unsigned) i_VertexCount, _FALSE);
01107
01108 vi_InducedVertexDegrees.clear();
01109 vi_InducedVertexDegrees.resize((unsigned) i_VertexCount, _FALSE);
01110
01111 vd_IncludedVertices.clear();
01112 vd_IncludedVertices.resize((unsigned) i_VertexCount, _UNKNOWN);
01113
01114 i_ParentVertex = _UNKNOWN;
01115
01116 i_SetCount = (signed) vi_Sets.size();
01117
01118 for(i=0; i<i_SetCount; i++)
01119 {
01120 vli_GroupedInducedVertexDegrees.clear();
01121 vli_GroupedInducedVertexDegrees.resize((unsigned) STEP_UP(i_MaximumVertexDegree));
01122
01123 vlit_VertexLocations.clear();
01124 vlit_VertexLocations.resize((unsigned) i_VertexCount);
01125
01126 i_HighestInducedVertexDegree = _UNKNOWN;
01127
01128 i_SetSize = (signed) mivi_VertexSets[vi_Sets[i]].size();
01129
01130 for(j=0; j<i_SetSize; j++)
01131 {
01132 i_PresentVertex = mivi_VertexSets[vi_Sets[i]][j];
01133
01134 vd_IncludedVertices[i_PresentVertex] = _FALSE;
01135
01136 if(vi_InducedVertexDegrees[i_PresentVertex] != _FALSE)
01137 {
01138 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].erase(vlit_VertexLocations[i_PresentVertex]);
01139 }
01140
01141 vi_InducedVertexDegrees[i_PresentVertex]++;
01142
01143 if(i_HighestInducedVertexDegree < vi_InducedVertexDegrees[i_PresentVertex])
01144 {
01145 i_HighestInducedVertexDegree = vi_InducedVertexDegrees[i_PresentVertex];
01146 }
01147 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].push_front(i_PresentVertex);
01148
01149 vlit_VertexLocations[i_PresentVertex] = vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_PresentVertex]].begin();
01150 }
01151
01152 #if DEBUG == 5103
01153
01154 int k;
01155
01156 list<int>::iterator lit_ListIterator;
01157
01158 cout<<endl;
01159 cout<<"DEBUG 5103 | Hessian Evaluation | Induced Vertex Degrees | Set "<<STEP_UP(i)<<endl;
01160 cout<<endl;
01161
01162 for(j=0; j<STEP_UP(i_HighestInducedVertexDegree); j++)
01163 {
01164 i_SetSize = (signed) vli_GroupedInducedVertexDegrees[j].size();
01165
01166 if(i_SetSize == _FALSE)
01167 {
01168 continue;
01169 }
01170
01171 k = _FALSE;
01172
01173 cout<<"Degree "<<j<<"\t"<<" : ";
01174
01175 for(lit_ListIterator=vli_GroupedInducedVertexDegrees[j].begin(); lit_ListIterator!=vli_GroupedInducedVertexDegrees[j].end(); lit_ListIterator++)
01176 {
01177 if(k == STEP_DOWN(i_SetSize))
01178 {
01179 cout<<STEP_UP(*lit_ListIterator)<<" ("<<i_SetSize<<")"<<endl;
01180 }
01181 else
01182 {
01183 cout<<STEP_UP(*lit_ListIterator)<<", ";
01184 }
01185
01186 k++;
01187 }
01188 }
01189
01190 #endif
01191
01192 #if DEBUG == 5103
01193
01194 cout<<endl;
01195 cout<<"DEBUG 5103 | Hessian Evaluation | Retrieved Elements"<<"| Set "<<STEP_UP(i)<<endl;
01196 cout<<endl;
01197
01198 #endif
01199
01200 for (int index = 0; index < i_VertexCount; index++) {
01201 if(vi_EvaluatedDiagonals[index] == _FALSE)
01202 {
01203 d_Value = dp2_CompressedMatrix[index][vi_VertexColors[index]];
01204
01205 #if DEBUG == 5103
01206
01207 cout<<"Element["<<STEP_UP(index)<<"]["<<STEP_UP(index)<<"] = "<<d_Value<<endl;
01208
01209 #endif
01210 v2i_VertexAdjacency[index].push_back(index);
01211 v2d_NonzeroAdjacency[index].push_back(d_Value);
01212
01213 vi_EvaluatedDiagonals[index] = _TRUE;
01214
01215 }
01216 }
01217
01218 for ( ; ; )
01219 {
01220 if(vli_GroupedInducedVertexDegrees[_TRUE].empty())
01221 {
01222 i_LeafVertex = vli_GroupedInducedVertexDegrees[_FALSE].front();
01223
01224 vi_InducedVertexDegrees[i_LeafVertex] = _FALSE;
01225
01226 vd_IncludedVertices[i_LeafVertex] = _UNKNOWN;
01227
01228 break;
01229 }
01230
01231 i_LeafVertex = vli_GroupedInducedVertexDegrees[_TRUE].front();
01232
01233 vli_GroupedInducedVertexDegrees[_TRUE].pop_front();
01234
01235
01236 for(j=vi_Vertices[i_LeafVertex]; j<vi_Vertices[STEP_UP(i_LeafVertex)]; j++)
01237 {
01238 if(vd_IncludedVertices[vi_Edges[j]] != _UNKNOWN)
01239 {
01240 i_ParentVertex = vi_Edges[j];
01241
01242 break;
01243 }
01244 }
01245
01246 d_Value = dp2_CompressedMatrix[i_LeafVertex][vi_VertexColors[i_ParentVertex]] - vd_IncludedVertices[i_LeafVertex];
01247
01248 vd_IncludedVertices[i_ParentVertex] += d_Value;
01249
01250 vi_InducedVertexDegrees[i_LeafVertex] = _FALSE;
01251 vd_IncludedVertices[i_LeafVertex] = _UNKNOWN;
01252 if(vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].size()>1) {
01253 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].erase(vlit_VertexLocations[i_ParentVertex]);
01254 }
01255 else {
01256 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].pop_back();
01257 }
01258
01259 vi_InducedVertexDegrees[i_ParentVertex]--;
01260 vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].push_back(i_ParentVertex);
01261
01262
01263 vlit_VertexLocations[i_ParentVertex] = vli_GroupedInducedVertexDegrees[vi_InducedVertexDegrees[i_ParentVertex]].end();
01264 --vlit_VertexLocations[i_ParentVertex];
01265
01266 v2i_VertexAdjacency[i_LeafVertex].push_back(i_ParentVertex);
01267 v2d_NonzeroAdjacency[i_LeafVertex].push_back(d_Value);
01268
01269 v2i_VertexAdjacency[i_ParentVertex].push_back(i_LeafVertex);
01270 v2d_NonzeroAdjacency[i_ParentVertex].push_back(d_Value);
01271
01272
01273 #if DEBUG == 5103
01274
01275 cout<<"Element["<<STEP_UP(i_LeafVertex)<<"]["<<STEP_UP(i_ParentVertex)<<"] = "<<d_Value<<endl;
01276 #endif
01277
01278 }
01279 }
01280
01281 vector<unsigned int> RowIndex;
01282 vector<unsigned int> ColumnIndex;
01283 vector<double> HessianValue;
01284
01285
01286 for(i=0; i<i_VertexCount; i++) {
01287 int NumOfNonzeros = uip2_HessianSparsityPattern[i][0];
01288 i_VertexDegree = (signed) v2i_VertexAdjacency[i].size();
01289 for(j=1; j<=NumOfNonzeros; j++) {
01290 int targetColumnID = uip2_HessianSparsityPattern[i][j];
01291 if(targetColumnID<i) continue;
01292 for (int k=0; k<i_VertexDegree; k++) {
01293 if(targetColumnID == v2i_VertexAdjacency[i][k]) {
01294 HessianValue.push_back(v2d_NonzeroAdjacency[i][k]);
01295 break;
01296 }
01297 }
01298 RowIndex.push_back(i);
01299 ColumnIndex.push_back(uip2_HessianSparsityPattern[i][j]);
01300 }
01301 }
01302
01303 unsigned int numOfNonZeros = RowIndex.size();
01304 (*ip2_RowIndex) = new unsigned int[numOfNonZeros];
01305 (*ip2_ColumnIndex) = new unsigned int[numOfNonZeros];
01306 (*dp2_HessianValue) = new double[numOfNonZeros];
01307
01308 for(int i=0; i < numOfNonZeros; i++) {
01309 (*ip2_RowIndex)[i] = RowIndex[i];
01310 (*ip2_ColumnIndex)[i] = ColumnIndex[i];
01311 (*dp2_HessianValue)[i] = HessianValue[i];
01312 }
01313
01314
01315 CF_available = true;
01316 i_CF_rowCount = numOfNonZeros;
01317 ip_CF_RowIndex = *ip2_RowIndex;
01318 ip_CF_ColumnIndex = *ip2_ColumnIndex;
01319 dp_CF_Value = *dp2_HessianValue;
01320
01321 return (numOfNonZeros);
01322 }
01323 }