From efc6996e97a236b14e1375a1953cdb1ddbdcf024 2017-03-13 17:28:41 From: Tom Bannink Date: 2017-03-13 17:28:41 Subject: [PATCH] Add adaptive mixing time --- diff --git a/cpp/switchchain.cpp b/cpp/switchchain.cpp index 559d9fa4814dbd2f820d3479cb823654384822fc..23d96eda05bceffb77cd0a7ed52c0f6a575ae204 100644 --- a/cpp/switchchain.cpp +++ b/cpp/switchchain.cpp @@ -36,21 +36,17 @@ class SwitchChain { } bool doMove() { - Edge e1, e2; int e1index, e2index; int timeout = 0; // Keep regenerating while conflicting edges do { e1index = edgeDistribution(mt); e2index = edgeDistribution(mt); - e1 = g.getEdge(e1index); - e2 = g.getEdge(e2index); - ++timeout; - if (timeout % 100 == 0) { + if (++timeout % 100 == 0) { std::cerr << "Warning: sampled " << timeout << " random edges but they keep conflicting.\n"; } - } while (edgeConflicts(e1, e2)); + } while (edgeConflicts(g.getEdge(e1index), g.getEdge(e2index))); // Consider one of the three possible permutations // 1) e1.u - e1.v and e2.u - e2.v (original) @@ -117,12 +113,10 @@ int main() { std::cout << "Running n = " << numVertices << ", tau = " << tau << ". \t" << std::flush; - constexpr int mixingTime = 40000; - constexpr int measureTime = 20000; + int mixingTime = (32.0f - 26.0f*(tau - 2.0f)) * numVertices; //40000; + constexpr int measurements = 50; constexpr int measureSkip = 200; // Take a sample every ... steps - constexpr int measurements = - (measureTime - 1) / measureSkip + 1; int movesDone = 0; int triangles[measurements]; @@ -131,15 +125,18 @@ int main() { if (chain.doMove()) ++movesDone; } - for (int i = 0; i < measureTime; ++i) { - if (chain.doMove()) - ++movesDone; - if (i % measureSkip == 0) - triangles[i / measureSkip] = chain.g.countTriangles(); + for (int i = 0; i < measurements; ++i) { + for (int j = 0; j < measureSkip; ++j) + if (chain.doMove()) + ++movesDone; + triangles[i] = chain.g.countTriangles(); } - std::cout << movesDone << '/' << mixingTime + measureTime - << " moves succeeded." << std::endl; + std::cout << movesDone << '/' << mixingTime + measurements * measureSkip + << " moves succeeded (" + << 100.0f * float(movesDone) / + float(mixingTime + measurements * measureSkip) + << "%)." << std::endl; if (outputComma) outfile << ','; diff --git a/showgraphs.m b/showgraphs.m index 727c8ecdb75c486f50a8c841da4e78a4d549d792..2d71c7db010b24428e83e09aab7709ae7920c6c3 100644 --- a/showgraphs.m +++ b/showgraphs.m @@ -64,7 +64,7 @@ Grid[Partition[gs,10],Frame->All] (*Data import and data merge*) -gsraw=Import[NotebookDirectory[]<>"data/graphdata_merged.m"]; +gsraw=Import[NotebookDirectory[]<>"data/graphdata2.m"]; gsraw=SortBy[gsraw,#[[1,1]]&]; (* Sort by n *) @@ -109,7 +109,7 @@ Sum[tmp[[i,j-i]]*tmp[[j,k-j]]*tmp[[i,k-i]],{i,1,n-2},{j,i+1,n-1},{k,j+1,n}] (* gdata[[ tau index, n index, run index , {ntau, #tris, ds} ]] *) -avgAndProp=ParallelMap[{getProperty[#[[3]]],Mean[#[[2,1;;-1]]]}&,gdata[[2,4,1;;100]]]; +avgAndProp=ParallelMap[{getProperty[#[[3]]],Mean[#[[2,1;;-1]]]}&,gdata[[2,2,1;;100]]]; Show[ListPlot[avgAndProp,AxesOrigin->{0,0},AxesLabel->{"degree-sequence-property","<#triangles>"},AspectRatio->1],Plot[x,{x,1,1000}]] @@ -120,17 +120,47 @@ Show[ListPlot[avgAndProp,AxesOrigin->{0,0},AxesLabel->{"degree-sequence-property numPlots=20; -selectedData=gsraw[[-numPlots;;-1]]; +selectedData=gdata[[5,-1]][[-numPlots;;-1]]; minCount=Min[Map[Min[#[[2]]]&,selectedData]]; maxCount=Max[Map[Max[#[[2]]]&,selectedData]]; maxTime=Max[Map[Length[#[[2]]]&,selectedData]]; skipPts=Max[1,Round[maxTime/100]]; (* Plotting every point is slow. Plot only once per `skipPts` timesteps *) coarseData=Map[#[[2,1;;-1;;skipPts]]&,selectedData]; labels=Map["{n,tau} = "<>ToString[#[[1]]]&,selectedData]; -ListPlot[coarseData,Joined->True,PlotRange->{minCount,maxCount},DataRange->{0,maxTime},PlotLegends->labels] +ListPlot[coarseData,Joined->True,PlotRange->{minCount,maxCount},DataRange->{0,200*maxTime},PlotLegends->labels] (* Map[ListPlot[#,Joined->True,PlotRange\[Rule]{minCount,maxCount},DataRange\[Rule]{0,maxTime}]&,coarseData] *) +(* ::Subsection:: *) +(*Compute 'mixing time'*) + + +(* Compute average of last part and check when the value drops below that for the first time *) +getMixingTime[values_]:=Module[{avg}, + (* average over the last 20 percent *) + avg=Mean[values[[-Round[Length[values]/5];;-1]]]; + FirstPosition[values,_?(#<=avg&)][[1]] +] +(* gdata[[ tau index, n index, run index , {ntau, #tris, ds} ]] *) +mixingTimes=Map[{#[[1,1]],(1/#[[1,1]])200 * getMixingTime[#[[2]]]}&,gdata,{3}]; +mixingTimesBars=Map[ + {{#[[1,1]],Mean[#[[All,2]]]},ErrorBar[StandardDeviation[#[[All,2]]]/Sqrt[Length[#]]]}& +,mixingTimes,{2}]; +ErrorListPlot[mixingTimesBars,Joined->True,PlotMarkers->Automatic,AxesLabel->{"n","~~mixing time divided by n"},PlotLegends->taulabels] + + +(* For n fixed, look at function of tau *) +mixingTimes=Map[{#[[1,2]],(1/#[[1,1]])200 * getMixingTime[#[[2]]]}&,gdata,{3}]; +mixingTimesBars=Map[ + {{#[[1,1]],Mean[#[[All,2]]]},ErrorBar[StandardDeviation[#[[All,2]]]]}& +,mixingTimes[[All,-1]],{1}]; + + +Show[ +ErrorListPlot[mixingTimesBars,Joined->True,PlotMarkers->Automatic,AxesLabel->{"tau","~~mixing time divided by n, for n = 1000"},PlotRange->{{2,3},{0,30}}] +,Plot[(32-26(tau-2)),{tau,2,3}]] + + (* ::Subsection:: *) (*Plot average #triangles vs n*) @@ -143,11 +173,11 @@ nlabels=Map["n = "<>ToString[#]&,averagesGrouped[[1,All,1,1,1]]]; taulabels=Map["tau = "<>ToString[#]&,averagesGrouped[[All,1,1,1,2]]]; averagesErrorBars=Map[ {{#[[1,1,1]],Mean[#[[All,2]]]}, -ErrorBar[StandardDeviation[#[[All,2]]]/Sqrt[Length[#]]] +ErrorBar[StandardDeviation[#[[All,2]]]] }&,averagesGrouped,{2}]; -ErrorListPlot[averagesErrorBars,Joined->True,PlotMarkers->Automatic,AxesLabel->{"n","\[LeftAngleBracket]triangles\[RightAngleBracket]"},PlotLegends->taulabels] +ErrorListPlot[averagesErrorBars,Joined->True,PlotMarkers->Automatic,PlotRange->All,AxesLabel->{"n","\[LeftAngleBracket]triangles\[RightAngleBracket]"},PlotLegends->taulabels] ListLogLogPlot[averagesErrorBars[[All,All,1]],Joined->True,PlotMarkers->Automatic,AxesLabel->{"n","\[LeftAngleBracket]triangles\[RightAngleBracket]"},PlotLegends->taulabels] @@ -166,7 +196,7 @@ Show[ListLogLogPlot[averagesErrorBars[[All,All,1]],PlotMarkers->Automatic,AxesLa tauValues=averagesGrouped[[All,1,1,1,2]]; exponents=Transpose[{tauValues,fits[[All,2,1]]}]; -Show[ListPlot[exponents,Joined->True,PlotMarkers->Automatic,AxesLabel->{"tau","triangle-law-exponent"}],Plot[3/2(3-tau),{tau,2,3}]] +Show[ListPlot[exponents,Joined->True,PlotMarkers->Automatic,AxesLabel->{"tau","triangle-law-exponent"},PlotRange->{{2,3},{0,1.6}}],Plot[3/2(3-tau),{tau,2,3}]] (* ::Subsection:: *)