SNAP Library 2.0, User Reference  2013-05-13 16:33:57
SNAP, a general purpose, high performance system for analysis and manipulation of large networks
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
ggen.cpp
Go to the documentation of this file.
00001 
00002 // Graph Generators
00003 namespace TSnap {
00004 
00005 PBPGraph GenRndBipart(const int& LeftNodes, const int& RightNodes, const int& Edges, TRnd& Rnd) {
00006   PBPGraph G = TBPGraph::New();
00007   for (int i = 0; i < LeftNodes; i++) { G->AddNode(i, true); }
00008   for (int i = 0; i < RightNodes; i++) { G->AddNode(LeftNodes+i, false); }
00009   IAssertR(Edges <= LeftNodes*RightNodes, "Too many edges in the bipartite graph!");
00010   for (int edges = 0; edges < Edges; ) {
00011     const int LNId = Rnd.GetUniDevInt(LeftNodes);
00012     const int RNId = LeftNodes + Rnd.GetUniDevInt(RightNodes);
00013     if (G->AddEdge(LNId, RNId) != -2) { edges++; } // is new edge
00014   }
00015   return G;
00016 }
00017 
00018 PUNGraph GenRndDegK(const int& Nodes, const int& NodeDeg, const int& NSwitch, TRnd& Rnd) {
00019   // create degree sequence
00020   TIntV DegV(Nodes, 0);
00021   int DegSum=0;
00022   for (int i = 0; i < Nodes; i++) {
00023     DegV.Add(NodeDeg);
00024     DegSum += NodeDeg;
00025   }
00026   IAssert(DegSum % 2 == 0);
00027   PUNGraph G = GenDegSeq(DegV, Rnd); // get some graph that obeys the degree sequnce
00028   return GenRewire(G, NSwitch, Rnd);  // make it random
00029 }
00030 
00034 PUNGraph GenRndPowerLaw(const int& Nodes, const double& PowerExp, const bool& ConfModel, TRnd& Rnd) {
00035   TIntV DegSeqV;
00036   uint DegSum=0;
00037   for (int n = 0; n < Nodes; n++) {
00038     const int Val = (int) TMath::Round(Rnd.GetPowerDev(PowerExp));
00039     if (! (Val >= 1 && Val < Nodes/2)) { n--; continue; } // skip nodes with too large degree
00040     DegSeqV.Add(Val);
00041     DegSum += Val;
00042   }
00043   printf("%d nodes, %u edges\n", Nodes, DegSum);
00044   if (DegSum % 2 == 1) { DegSeqV[0] += 1; }
00045   if (ConfModel) {
00046     // use configuration model -- fast but does not exactly obey the degree sequence
00047     return GenConfModel(DegSeqV, Rnd);
00048   } else {
00049     PUNGraph G = TSnap::GenDegSeq(DegSeqV, Rnd);
00050     return TSnap::GenRewire(G, 10, Rnd);
00051   }
00052 }
00053 
00058 PUNGraph GenDegSeq(const TIntV& DegSeqV, TRnd& Rnd) {
00059   const int Nodes = DegSeqV.Len();
00060   PUNGraph GraphPt = TUNGraph::New();
00061   TUNGraph& Graph = *GraphPt;
00062   Graph.Reserve(Nodes, -1);
00063   TIntH DegH(DegSeqV.Len(), true);
00064   
00065   IAssertR(DegSeqV.IsSorted(false), "DegSeqV must be sorted in descending order.");
00066   int DegSum=0, edge=0;
00067   for (int node = 0; node < Nodes; node++) {
00068     IAssert(Graph.AddNode(node) == node);
00069     DegH.AddDat(node, DegSeqV[node]);
00070     DegSum += DegSeqV[node];
00071   }
00072   IAssert(DegSum % 2 == 0);
00073   while (! DegH.Empty()) {
00074     // pick random nodes and connect
00075     const int NId1 = DegH.GetKey(DegH.GetRndKeyId(TInt::Rnd, 0.5));
00076     const int NId2 = DegH.GetKey(DegH.GetRndKeyId(TInt::Rnd, 0.5));
00077     IAssert(DegH.IsKey(NId1) && DegH.IsKey(NId2));
00078     if (NId1 == NId2) {
00079       if (DegH.GetDat(NId1) == 1) { continue; }
00080       // find rnd edge, break it, and connect the endpoints to the nodes
00081       const TIntPr Edge = TSnapDetail::GetRndEdgeNonAdjNode(GraphPt, NId1, -1);
00082       if (Edge.Val1==-1) { continue; }
00083       Graph.DelEdge(Edge.Val1, Edge.Val2);
00084       Graph.AddEdge(Edge.Val1, NId1);
00085       Graph.AddEdge(NId1, Edge.Val2);
00086       if (DegH.GetDat(NId1) == 2) { DegH.DelKey(NId1); }
00087       else { DegH.GetDat(NId1) -= 2; }
00088     } else {
00089       if (! Graph.IsEdge(NId1, NId2)) {
00090         Graph.AddEdge(NId1, NId2); }  // good edge
00091       else {
00092         // find rnd edge, break and cross-connect
00093         const TIntPr Edge = TSnapDetail::GetRndEdgeNonAdjNode(GraphPt, NId1, NId2);
00094         if (Edge.Val1==-1) {continue; }
00095         Graph.DelEdge(Edge.Val1, Edge.Val2);
00096         Graph.AddEdge(NId1, Edge.Val1);
00097         Graph.AddEdge(NId2, Edge.Val2);
00098       }
00099       if (DegH.GetDat(NId1)==1) { DegH.DelKey(NId1); }
00100       else { DegH.GetDat(NId1) -= 1; }
00101       if (DegH.GetDat(NId2)==1) { DegH.DelKey(NId2); }
00102       else { DegH.GetDat(NId2) -= 1; }
00103     }
00104     if (++edge % 1000 == 0) {
00105       printf("\r %dk / %dk", edge/1000, DegSum/2000); }
00106   }
00107   return GraphPt;
00108 }
00109 
00110 
00119 PUNGraph GenConfModel(const TIntV& DegSeqV, TRnd& Rnd) {
00120   const int Nodes = DegSeqV.Len();
00121   PUNGraph GraphPt = TUNGraph::New();
00122   TUNGraph& Graph = *GraphPt;
00123   Graph.Reserve(Nodes, -1);
00124   TIntV NIdDegV(DegSeqV.Len(), 0);
00125   int DegSum=0, edges=0;
00126   for (int node = 0; node < Nodes; node++) {
00127     Graph.AddNode(node);
00128     for (int d = 0; d < DegSeqV[node]; d++) { NIdDegV.Add(node); }
00129     DegSum += DegSeqV[node];
00130   }
00131   NIdDegV.Shuffle(Rnd);
00132   TIntPrSet EdgeH(DegSum/2); // set of all edges, is faster than graph edge lookup
00133   if (DegSum % 2 != 0) {
00134     printf("Seg seq is odd [%d]: ", DegSeqV.Len());
00135     for (int d = 0; d < TMath::Mn(100, DegSeqV.Len()); d++) { printf("  %d", (int)DegSeqV[d]); }
00136     printf("\n");
00137   }
00138   int u=0, v=0;
00139   for (int c = 0; NIdDegV.Len() > 1; c++) {
00140     u = Rnd.GetUniDevInt(NIdDegV.Len());
00141     while ((v = Rnd.GetUniDevInt(NIdDegV.Len())) == u) { }
00142     if (u > v) { Swap(u, v); }
00143     const int E1 = NIdDegV[u];
00144     const int E2 = NIdDegV[v];
00145     if (v == NIdDegV.Len()-1) { NIdDegV.DelLast(); } 
00146     else { NIdDegV[v] = NIdDegV.Last();  NIdDegV.DelLast(); }
00147     if (u == NIdDegV.Len()-1) { NIdDegV.DelLast(); } 
00148     else { NIdDegV[u] = NIdDegV.Last();  NIdDegV.DelLast(); }
00149     if (E1 == E2 || EdgeH.IsKey(TIntPr(E1, E2))) { continue; }
00150     EdgeH.AddKey(TIntPr(E1, E2));
00151     Graph.AddEdge(E1, E2);
00152     edges++;
00153     if (c % (DegSum/100+1) == 0) { printf("\r configuration model: iter %d: edges: %d, left: %d", c, edges, NIdDegV.Len()/2); }
00154   }
00155   printf("\n");
00156   return GraphPt;
00157 }
00158 
00165 PUNGraph GenRewire(const PUNGraph& OrigGraph, const int& NSwitch, TRnd& Rnd) {
00166   const int Nodes = OrigGraph->GetNodes();
00167   const int Edges = OrigGraph->GetEdges();
00168   PUNGraph GraphPt = TUNGraph::New();
00169   TUNGraph& Graph = *GraphPt;
00170   Graph.Reserve(Nodes, -1);
00171   TExeTm ExeTm;
00172   // generate a graph that satisfies the constraints
00173   printf("Randomizing edges (%d, %d)...\n", Nodes, Edges);
00174   TIntPrSet EdgeSet(Edges);
00175   for (TUNGraph::TNodeI NI = OrigGraph->BegNI(); NI < OrigGraph->EndNI(); NI++) {
00176     const int NId = NI.GetId();
00177     for (int e = 0; e < NI.GetOutDeg(); e++) {
00178       if (NId <= NI.GetOutNId(e)) { continue; }
00179       EdgeSet.AddKey(TIntPr(NId, NI.GetOutNId(e)));
00180     }
00181     Graph.AddNode(NI.GetId());
00182   }
00183   // edge switching
00184   uint skip=0;
00185   for (uint swps = 0; swps < 2*uint(Edges)*uint(NSwitch); swps++) {
00186     const int keyId1 = EdgeSet.GetRndKeyId(Rnd);
00187     const int keyId2 = EdgeSet.GetRndKeyId(Rnd);
00188     if (keyId1 == keyId2) { skip++; continue; }
00189     const TIntPr& E1 = EdgeSet[keyId1];
00190     const TIntPr& E2 = EdgeSet[keyId2];
00191     TIntPr NewE1(E1.Val1, E2.Val1), NewE2(E1.Val2, E2.Val2);
00192     if (NewE1.Val1 > NewE1.Val2) { Swap(NewE1.Val1, NewE1.Val2); }
00193     if (NewE2.Val1 > NewE2.Val2) { Swap(NewE2.Val1, NewE2.Val2); }
00194     if (NewE1!=NewE2 && NewE1.Val1!=NewE1.Val2 && NewE2.Val1!=NewE2.Val2 && ! EdgeSet.IsKey(NewE1) && ! EdgeSet.IsKey(NewE2)) {
00195       EdgeSet.DelKeyId(keyId1);  EdgeSet.DelKeyId(keyId2);
00196       EdgeSet.AddKey(TIntPr(NewE1));
00197       EdgeSet.AddKey(TIntPr(NewE2));
00198     } else { skip++; }
00199     if (swps % Edges == 0) {
00200       printf("\r  %uk/%uk: %uk skip [%s]", swps/1000u, 2*uint(Edges)*uint(NSwitch)/1000u, skip/1000u, ExeTm.GetStr());
00201       if (ExeTm.GetSecs() > 2*3600) { printf(" *** Time limit!\n"); break; } // time limit 2 hours
00202     }
00203   }
00204   printf("\r  total %uk switchings attempted, %uk skiped  [%s]\n", 2*uint(Edges)*uint(NSwitch)/1000u, skip/1000u, ExeTm.GetStr());
00205   for (int e = 0; e < EdgeSet.Len(); e++) {
00206     Graph.AddEdge(EdgeSet[e].Val1, EdgeSet[e].Val2); }
00207   return GraphPt;
00208 }
00209 
00216 PNGraph GenRewire(const PNGraph& OrigGraph, const int& NSwitch, TRnd& Rnd) {
00217   const int Nodes = OrigGraph->GetNodes();
00218   const int Edges = OrigGraph->GetEdges();
00219   PNGraph GraphPt = TNGraph::New();
00220   TNGraph& Graph = *GraphPt;
00221   Graph.Reserve(Nodes, -1);
00222   TExeTm ExeTm;
00223   // generate a graph that satisfies the constraints
00224   printf("Randomizing edges (%d, %d)...\n", Nodes, Edges);
00225   TIntPrSet EdgeSet(Edges);
00226   for (TNGraph::TNodeI NI = OrigGraph->BegNI(); NI < OrigGraph->EndNI(); NI++) {
00227     const int NId = NI.GetId();
00228     for (int e = 0; e < NI.GetOutDeg(); e++) {
00229       EdgeSet.AddKey(TIntPr(NId, NI.GetOutNId(e))); }
00230     Graph.AddNode(NI);
00231   }
00232   // edge switching
00233   uint skip=0;
00234   for (uint swps = 0; swps < 2*uint(Edges)*uint(NSwitch); swps++) {
00235     const int keyId1 = EdgeSet.GetRndKeyId(Rnd);
00236     const int keyId2 = EdgeSet.GetRndKeyId(Rnd);
00237     if (keyId1 == keyId2) { skip++; continue; }
00238     const TIntPr& E1 = EdgeSet[keyId1];
00239     const TIntPr& E2 = EdgeSet[keyId2];
00240     TIntPr NewE1(E1.Val1, E2.Val1), NewE2(E1.Val2, E2.Val2);
00241     if (NewE1.Val1!=NewE2.Val1 && NewE1.Val2!=NewE2.Val1 && NewE1.Val2!=NewE2.Val1 && NewE1.Val2!=NewE2.Val2 && ! EdgeSet.IsKey(NewE1) && ! EdgeSet.IsKey(NewE2)) {
00242       EdgeSet.DelKeyId(keyId1);  EdgeSet.DelKeyId(keyId2);
00243       EdgeSet.AddKey(TIntPr(NewE1));
00244       EdgeSet.AddKey(TIntPr(NewE2));
00245     } else { skip++; }
00246     if (swps % Edges == 0) {
00247       printf("\r  %uk/%uk: %uk skip [%s]", swps/1000u, 2*uint(Edges)*uint(NSwitch)/1000u, skip/1000u, ExeTm.GetStr());
00248       if (ExeTm.GetSecs() > 2*3600) { printf(" *** Time limit!\n"); break; } // time limit 2 hours
00249     }
00250   }
00251   printf("\r  total %uk switchings attempted, %uk skiped  [%s]\n", 2*uint(Edges)*uint(NSwitch)/1000u, skip/1000u, ExeTm.GetStr());
00252   for (int e = 0; e < EdgeSet.Len(); e++) {
00253     Graph.AddEdge(EdgeSet[e].Val1, EdgeSet[e].Val2); }
00254   return GraphPt;
00255 }
00256 
00263 PBPGraph GenRewire(const PBPGraph& OrigGraph, const int& NSwitch, TRnd& Rnd) {
00264   const int Nodes = OrigGraph->GetNodes();
00265   const int Edges = OrigGraph->GetEdges();
00266   PBPGraph GraphPt = TBPGraph::New();
00267   TBPGraph& Graph = *GraphPt;
00268   Graph.Reserve(Nodes, -1);
00269   TExeTm ExeTm;
00270   // generate a graph that satisfies the constraints
00271   printf("Randomizing edges (%d, %d)...\n", Nodes, Edges);
00272   TIntPrSet EdgeSet(Edges);
00273   for (TBPGraph::TNodeI NI = OrigGraph->BegLNI(); NI < OrigGraph->EndLNI(); NI++) {
00274     const int NId = NI.GetId();
00275     for (int e = 0; e < NI.GetOutDeg(); e++) {
00276       EdgeSet.AddKey(TIntPr(NId, NI.GetOutNId(e))); } // edges left-->right
00277     Graph.AddNode(NI.GetId(), true); } // left nodes
00278   for (TBPGraph::TNodeI NI = OrigGraph->BegRNI(); NI < OrigGraph->EndRNI(); NI++) {
00279     Graph.AddNode(NI.GetId(), false); } // right nodes
00280   IAssert(EdgeSet.Len() == Edges);
00281   // edge switching
00282   uint skip=0;
00283   for (uint swps = 0; swps < 2*uint(Edges)*uint(NSwitch); swps++) {
00284     const int keyId1 = EdgeSet.GetRndKeyId(Rnd);
00285     const int keyId2 = EdgeSet.GetRndKeyId(Rnd);
00286     if (keyId1 == keyId2) { skip++; continue; }
00287     const TIntPr& E1 = EdgeSet[keyId1];
00288     const TIntPr& E2 = EdgeSet[keyId2];
00289     TIntPr NewE1(E1.Val1, E2.Val2), NewE2(E2.Val1, E1.Val2);
00290     if (NewE1!=NewE2 && NewE1.Val1!=NewE1.Val2 && NewE2.Val1!=NewE2.Val2 && ! EdgeSet.IsKey(NewE1) && ! EdgeSet.IsKey(NewE2)) {
00291       EdgeSet.DelKeyId(keyId1);  EdgeSet.DelKeyId(keyId2);
00292       EdgeSet.AddKey(TIntPr(NewE1));
00293       EdgeSet.AddKey(TIntPr(NewE2));
00294     } else { skip++; }
00295     if (swps % Edges == 0) {
00296       printf("\r  %uk/%uk: %uk skip [%s]", swps/1000u, 2*uint(Edges)*uint(NSwitch)/1000u, skip/1000u, ExeTm.GetStr());
00297       if (ExeTm.GetSecs() > 2*3600) { printf(" *** Time limit!\n"); break; } // time limit 2 hours
00298     }
00299   }
00300   printf("\r  total %uk switchings attempted, %uk skiped  [%s]\n", 2*uint(Edges)*uint(NSwitch)/1000u, skip/1000u, ExeTm.GetStr());
00301   for (int e = 0; e < EdgeSet.Len(); e++) {
00302     Graph.AddEdge(EdgeSet[e].Val1, EdgeSet[e].Val2); }
00303   return GraphPt;
00304 }
00305 
00310 PUNGraph GenPrefAttach(const int& Nodes, const int& NodeOutDeg, TRnd& Rnd) {
00311   PUNGraph GraphPt = PUNGraph::New();
00312   TUNGraph& Graph = *GraphPt;
00313   Graph.Reserve(Nodes, NodeOutDeg*Nodes);
00314   TIntV NIdV(NodeOutDeg*Nodes, 0);
00315   // first edge
00316   Graph.AddNode(0);  Graph.AddNode(1);
00317   NIdV.Add(0);  NIdV.Add(1);
00318   Graph.AddEdge(0, 1);
00319   TIntSet NodeSet;
00320   for (int node = 2; node < Nodes; node++) {
00321     NodeSet.Clr(false);
00322     while (NodeSet.Len() < NodeOutDeg && NodeSet.Len() < node) {
00323       NodeSet.AddKey(NIdV[TInt::Rnd.GetUniDevInt(NIdV.Len())]);
00324     }
00325     const int N = Graph.AddNode();
00326     for (int i = 0; i < NodeSet.Len(); i++) {
00327       Graph.AddEdge(N, NodeSet[i]);
00328       NIdV.Add(N);
00329       NIdV.Add(NodeSet[i]);
00330     }
00331   }
00332   return GraphPt;
00333 }
00334 
00335 namespace TSnapDetail {
00337 void GetSphereDev(const int& Dim, TRnd& Rnd, TFltV& ValV) {
00338   if (ValV.Len() != Dim) { ValV.Gen(Dim); }
00339   double Length = 0.0;
00340   for (int i = 0; i < Dim; i++) {
00341     ValV[i] = Rnd.GetNrmDev();
00342     Length += TMath::Sqr(ValV[i]); }
00343   Length = 1.0 / sqrt(Length);
00344   for (int i = 0; i < Dim; i++) {
00345     ValV[i] *= Length;
00346   }
00347 }
00348 } // namespace TSnapDetail
00349 
00355 PUNGraph GenGeoPrefAttach(const int& Nodes, const int& OutDeg, const double& Beta, TRnd& Rnd) {
00356   PUNGraph G = TUNGraph::New(Nodes, Nodes*OutDeg);
00357   TFltTrV PointV(Nodes, 0);
00358   TFltV ValV;
00359   // points on a sphere of radius 1/(2*pi)
00360   const double Rad = 0.5 * TMath::Pi;
00361   for (int i = 0; i < Nodes; i++) {
00362     TSnapDetail::GetSphereDev(3, Rnd, ValV);
00363     PointV.Add(TFltTr(Rad*ValV[0], Rad*ValV[1], Rad*ValV[2]));
00364   }
00365   const double R2 = TMath::Sqr(log((double) Nodes) / (pow((double) Nodes, 0.5-Beta)));
00366   TIntV DegV, NIdV;
00367   int SumDeg;
00368   for (int t = 0; t < Nodes; t++) {
00369     const int pid = t;
00370     const TFltTr& P1 = PointV[pid];
00371     // add node
00372     if (! G->IsNode(pid)) { G->AddNode(pid); }
00373     // find neighborhood
00374     DegV.Clr(false);  NIdV.Clr(false);  SumDeg=0;
00375     for (int p = 0; p < t; p++) {
00376       const TFltTr& P2 = PointV[p];
00377       if (TMath::Sqr(P1.Val1-P2.Val1)+TMath::Sqr(P1.Val2-P2.Val2)+TMath::Sqr(P1.Val3-P2.Val3) < R2) {
00378         NIdV.Add(p);
00379         DegV.Add(G->GetNI(p).GetDeg()+1);
00380         SumDeg += DegV.Last();
00381       }
00382     }
00383     // add edges
00384     for (int m = 0; m < OutDeg; m++) {
00385       const int rnd = Rnd.GetUniDevInt(SumDeg);
00386       int sum = 0, dst = -1;
00387       for (int s = 0; s < DegV.Len(); s++) {
00388         sum += DegV[s];
00389         if (rnd < sum) { dst=s;  break; }
00390       }
00391       if (dst != -1) {
00392         G->AddEdge(pid, NIdV[dst]);
00393         SumDeg -= DegV[dst];
00394         NIdV.Del(dst);  DegV.Del(dst);
00395       }
00396     }
00397   }
00398   return G;
00399 }
00400 
00406 PUNGraph GenSmallWorld(const int& Nodes, const int& NodeOutDeg, const double& RewireProb, TRnd& Rnd) {
00407   THashSet<TIntPr> EdgeSet(Nodes*NodeOutDeg);
00408   
00409   IAssertR(Nodes > NodeOutDeg, TStr::Fmt("Insufficient nodes for out degree, %d!", NodeOutDeg));
00410   for (int node = 0; node < Nodes; node++) {
00411     const int src = node;
00412     for (int edge = 1; edge <= NodeOutDeg; edge++) {
00413       int dst = (node+edge) % Nodes;      // edge to next neighbor
00414       if (Rnd.GetUniDev() < RewireProb) { // random edge
00415         dst = Rnd.GetUniDevInt(Nodes);
00416         while (dst == src || EdgeSet.IsKey(TIntPr(src, dst))) {
00417           dst = Rnd.GetUniDevInt(Nodes); }
00418       }
00419       EdgeSet.AddKey(TIntPr(src, dst));
00420     }
00421   }
00422   PUNGraph GraphPt = TUNGraph::New();
00423   TUNGraph& Graph = *GraphPt;
00424   Graph.Reserve(Nodes, EdgeSet.Len());
00425   int node;
00426   for (node = 0; node < Nodes; node++) {
00427     IAssert(Graph.AddNode(node) == node);
00428   }
00429   for (int edge = 0; edge < EdgeSet.Len(); edge++) {
00430     Graph.AddEdge(EdgeSet[edge].Val1, EdgeSet[edge].Val2);
00431   }
00432   Graph.Defrag();
00433   return GraphPt;
00434 }
00435 
00436 PNGraph GenForestFire(const int& Nodes, const double& FwdProb, const double& BckProb) {
00437   return TForestFire::GenGraph(Nodes, FwdProb, BckProb);
00438 }
00439 
00447 PNGraph GenCopyModel(const int& Nodes, const double& Beta, TRnd& Rnd) {
00448   PNGraph GraphPt = TNGraph::New();
00449   TNGraph& Graph = *GraphPt;
00450   Graph.Reserve(Nodes, Nodes);
00451   const int startNId = Graph.AddNode();
00452   Graph.AddEdge(startNId, startNId);
00453   for (int n = 1; n < Nodes; n++) {
00454     const int rnd = Graph.GetRndNId();
00455     const int NId = Graph.AddNode();
00456     if (Rnd.GetUniDev() < Beta) {
00457       Graph.AddEdge(NId, rnd); }
00458     else {
00459       const TNGraph::TNodeI NI = Graph.GetNI(rnd);
00460       const int rnd2 = Rnd.GetUniDevInt(NI.GetOutDeg());
00461       Graph.AddEdge(NId, NI.GetOutNId(rnd2));
00462     }
00463   }
00464   return GraphPt;
00465 }
00466 
00472 PNGraph GenRMat(const int& Nodes, const int& Edges, const double& A, const double& B, const double& C, TRnd& Rnd) {
00473   PNGraph GraphPt = TNGraph::New();
00474   TNGraph& Graph = *GraphPt;
00475   Graph.Reserve(Nodes, Edges);
00476   IAssert(A+B+C < 1.0);
00477   int rngX, rngY, offX, offY;
00478   int Depth=0, Collisions=0, Cnt=0, PctDone=0;
00479   const int EdgeGap = Edges / 100 + 1;
00480   // sum of parameters (probabilities)
00481   TVec<double> sumA(128, 0), sumAB(128, 0), sumAC(128, 0), sumABC(128, 0);  // up to 2^128 vertices ~ 3.4e38
00482   for (int i = 0; i < 128; i++) {
00483     const double a = A * (Rnd.GetUniDev() + 0.5);
00484     const double b = B * (Rnd.GetUniDev() + 0.5);
00485     const double c = C * (Rnd.GetUniDev() + 0.5);
00486     const double d = (1.0 - (A+B+C)) * (Rnd.GetUniDev() + 0.5);
00487     const double abcd = a+b+c+d;
00488     sumA.Add(a / abcd);
00489     sumAB.Add((a+b) / abcd);
00490     sumAC.Add((a+c) / abcd);
00491     sumABC.Add((a+b+c) / abcd);
00492   }
00493   // nodes
00494   for (int node = 0; node < Nodes; node++) {
00495     IAssert(Graph.AddNode(-1) == node);
00496   }
00497   // edges
00498   for (int edge = 0; edge < Edges; ) {
00499     rngX = Nodes;  rngY = Nodes;  offX = 0;  offY = 0;
00500     Depth = 0;
00501     // recurse the matrix
00502     while (rngX > 1 || rngY > 1) {
00503       const double RndProb = Rnd.GetUniDev();
00504       if (rngX>1 && rngY>1) {
00505         if (RndProb < sumA[Depth]) { rngX/=2; rngY/=2; }
00506         else if (RndProb < sumAB[Depth]) { offX+=rngX/2;  rngX-=rngX/2;  rngY/=2; }
00507         else if (RndProb < sumABC[Depth]) { offY+=rngY/2;  rngX/=2;  rngY-=rngY/2; }
00508         else { offX+=rngX/2;  offY+=rngY/2;  rngX-=rngX/2;  rngY-=rngY/2; }
00509       } else
00510       if (rngX>1) { // row vector
00511         if (RndProb < sumAC[Depth]) { rngX/=2; rngY/=2; }
00512         else { offX+=rngX/2;  rngX-=rngX/2;  rngY/=2; }
00513       } else
00514       if (rngY>1) { // column vector
00515         if (RndProb < sumAB[Depth]) { rngX/=2; rngY/=2; }
00516         else { offY+=rngY/2;  rngX/=2;  rngY-=rngY/2; }
00517       } else { Fail; }
00518       Depth++;
00519     }
00520     // add edge
00521     const int NId1 = offX;
00522     const int NId2 = offY;
00523     if (NId1 != NId2 && ! Graph.IsEdge(NId1, NId2)) {
00524       Graph.AddEdge(NId1, NId2);
00525       if (++Cnt > EdgeGap) {
00526         Cnt=0;  printf("\r  %d%% edges", ++PctDone); }
00527       edge++;
00528     } else {
00529       Collisions++; }
00530   }
00531   printf("\r  RMat: nodes:%d, edges:%d, Iterations:%d, Collisions:%d (%.1f%%).\n", Nodes, Edges,
00532     Edges+Collisions, Collisions, 100*Collisions/double(Edges+Collisions));
00533   Graph.Defrag();
00534   return GraphPt;
00535 }
00536 
00541 PNGraph GenRMatEpinions() {
00542   return GenRMat(75888, 508837, 0.550, 0.228, 0.212);
00543 }
00544 
00545 } // namespace TSnap