CounterexampleSearch.prejava

#include "macros.h"

import com.donhatchsw.compat.Format;

public class CounterexampleSearch
{
    private static String repeat(String s, int n) { return new String(new char[n]).replace("\0", s); }

    public static void main(String args[])
    {
        if (args.length != 2)
        {
            System.err.println("Usage: CounterexampleSearch <lookahead> <nSamples>");
            System.exit(1);
        }
        int lookahead = Integer.parseInt(args[0]); // e.g. 5
        int nSamples = Integer.parseInt(args[1]); // e.g. 1000

        int verboseLevel = 1;
        long seed = 0;
        java.util.Random rng = new java.util.Random(seed);
        int nInitialVerts = 3; // must be 3 or 4, not sure if 4 even works
        int maxVerts = 100; // we'll never get that high

        ConvexNoiseMaker convexNoiseMaker = new ConvexNoiseMaker(nInitialVerts, maxVerts);
        convexNoiseMaker.pushRandomVertex(rng, verboseLevel); // always add that deterministic vertex in the middle

        while (true)
        {
            OUT("      top of loop");
            // add lookahead more verts, in each of nSamples ways
            double bestGoodFrac = Double.NEGATIVE_INFINITY; // for curiosity
            double worstGoodFrac = Double.POSITIVE_INFINITY;
            double worstLookahead[][] = new double[lookahead][3];

            // Save the mesh
            Mesh savedMesh = new Mesh(convexNoiseMaker.getInternalMesh());

            FORI (iSample, nSamples)
            {
                OUT("          iSample="+iSample+"/"+nSamples+":");
                int nVerts = convexNoiseMaker.getInternalMesh().verts.size();
                OUT("              pushing "+lookahead+" random vert"+(lookahead==1?"":"s")+": "+nVerts+" -> "+(nVerts+lookahead));
                FORI (iLookahead, lookahead)
                {
                    // note, verboseLevel=2 makes it print "." each time, currently
                    convexNoiseMaker.pushRandomVertex(rng, verboseLevel);
                }

                // Extract the mesh
                // and delaunayize it "the right way".
                // XXX wait, what? do I need this, since inside out is just a triangle?  think about this
                Mesh triangulatedMesh = MeshUtils.delaunayized(convexNoiseMaker.getInternalMesh(),
                                                               false,  // wrapAroundSphereFlagValue,
                                                               1.,     // wrapSphereCurvatureValue,
                                                               true,  // centerSphereFlagValue,
                                                               false,  // calcInsideOutDualVertsFlagValue (this is the important one)
                                                               false); // slightlyVerbose
                Mesh trivalentMesh = Mesh.makeDualMesh(triangulatedMesh,
                                                       false, // includeNonArity3
                                                       false, // includeInsideOut
                                                       /*zeroVerts=*/false,
                                                       false, // wrapAroundSphereFlag
                                                       false, // centerSphereFlag
                                                       0.);   // wrapSphereCurvature)
                OUT("                  triangulated mesh: "+triangulatedMesh.verts.size()+" verts, "+triangulatedMesh.edges.size()+" edges, "+trivalentMesh.verts.size()+" faces");
                OUT("                  trivalent mesh: "+trivalentMesh.verts.size()+" verts, "+triangulatedMesh.edges.size()+" edges, "+triangulatedMesh.verts.size()+" faces");
                CHECK_EQ(triangulatedMesh.verts.size(), nVerts+lookahead);
                CHECK_EQ(triangulatedMesh.edges.size(), trivalentMesh.edges.size());

                // Actually, it's a triangulated mesh, so can make some more assertions.
                // Can we get everything in terms of the original num verts?
                // For a triangulated mesh:
                //      V+F=E+2  by euler
                //      F*3=E*2  because triangulated
                // so:
                //      E = F*3/2
                //      V+F=F*3/2+2
                //      V=F/2+2
                //      F=2*(V-2)=2*V-4
                //      E = (2*V-4)*3/2 = (V-2)*3
                CHECK_EQ(trivalentMesh.verts.size()+1, 2*(triangulatedMesh.verts.size()-2));
                CHECK_EQ(trivalentMesh.edges.size(), 2*3*(triangulatedMesh.verts.size()-2));
                CHECK_EQ(trivalentMesh.verts.size()+1, 2*(nVerts+lookahead-2));
                CHECK_EQ(trivalentMesh.edges.size(), 2*3*(nVerts+lookahead-2));
                // Therefore the number of edges in trivalent's spanning tree
                // is 1 less than its number of vertices,
                // i.e. equal to its number of non-infinite vertices.
                int nSpanningTreeEdges = trivalentMesh.verts.size()+1;

                //PRINT(triangulatedMesh);
                //PRINT(trivalentMesh);
                BigInt nGood = new BigInt(-1);
                BigInt nBad = new BigInt(-1);
                BigInt nTotal = new BigInt(-1);
                int max = 1000*1000*1000;



                double goodFrac;
                double badFrac;
                if (true) // TODO: make this an option?
                {
                    // Estimate by sampling random spanning trees, instead of counting all.
                    //       if percentage isn't very near 0 or 100,
                    //       then don't need that many samples.
                    // TODO: also, could give it a percentage it needs to beat;
                    //       if it's nowhere near beating it, it doesn't need to get very accurate.
                    goodFrac = MeshUtils.estimateGoodFrac(trivalentMesh, triangulatedMesh,
                                                          10, // maxMin
                                                          1000*1000, // maxTotal
                                                          rng);
                    badFrac = 1. - goodFrac;
                    if (goodFrac == 0.)
                    {
                        OUT(repeat("!", 1000));
                        CHECK(false); // XXX what should I do here? notify the press?
                    }
                }
                else
                {
                    MeshUtils.countGoodBad(trivalentMesh, triangulatedMesh,
                                           max,
                                           nGood, nBad, nTotal,
                                           0); // verboseLevel
                    goodFrac = nGood.doubleValue()/nTotal.doubleValue();
                    badFrac = nBad.doubleValue()/nTotal.doubleValue();
                    {
                        double goodPercent = goodFrac * 100.;
                        double badPercent = badFrac * 100.;
                        String goodPercentString = Format.sprintf("%.5g", goodPercent);
                        String badPercentString = Format.sprintf("%.5g", badPercent);
                        OUT("                  (good:bad)/total = ("+nGood+":"+nBad+")/"+nTotal+" = "+goodPercentString+"%:"+badPercentString+"%");
                        if (nGood.eq(0))
                        {
                            OUT(repeat("!", 1000));
                            CHECK(false); // XXX what should I do here? notify the press?
                        }
                    }
                }
                OUT("              goodPercent = "+goodFrac*100.+"%");
                OUT("                          = "+Math.pow(goodFrac, 1./nSpanningTreeEdges)*100.+"% ^ ("+nSpanningTreeEdges+" edges)");
                if (goodFrac > bestGoodFrac)
                {
                    bestGoodFrac = goodFrac;
                }
                if (goodFrac < worstGoodFrac)
                {
                    worstGoodFrac = goodFrac;
                    FORI (iLookahead, lookahead)
                    {
                        int iVert = triangulatedMesh.verts.size() - lookahead + iLookahead;
                        Mesh.Vertex v = triangulatedMesh.getVert(iVert);
                        worstLookahead[iLookahead][0] = v.x();
                        worstLookahead[iLookahead][1] = v.y();
                        worstLookahead[iLookahead][2] = v.z();
                    }

                    MeshUtils.saveToOFF(convexNoiseMaker.getInternalMesh(),
                                        null, // net
                                        "DUMP.bestSoFarRadical.off",
                                        true); // saveHomo (as opposed to heights above paraboloid)
                }

                OUT("              popping "+lookahead+" vert"+(lookahead==1?"":"s")+".");
                FORI (iLookahead, lookahead)
                {
                    convexNoiseMaker.popVertex();
                }
                convexNoiseMaker.setInternalMesh(new Mesh(savedMesh));

                OUT("              bestGoodPercent = "+bestGoodFrac*100.+"%");
                OUT("                               = "+Math.pow(bestGoodFrac, 1./nSpanningTreeEdges)*100.+"% ^ ("+nSpanningTreeEdges+" edges)");
                OUT("              worstGoodPercent = "+worstGoodFrac*100.+"%");
                OUT("                               = "+Math.pow(worstGoodFrac, 1./nSpanningTreeEdges)*100.+"% ^ ("+nSpanningTreeEdges+" edges)");
            } // for iSample
            CHECK_LE(worstGoodFrac, 1.);
            convexNoiseMaker.pushSpecificVertex(
                worstLookahead[0][0],
                worstLookahead[0][1],
                worstLookahead[0][2],
                null);
            PRINT(convexNoiseMaker.getInternalMesh());
            MeshUtils.saveToOFF(convexNoiseMaker.getInternalMesh(),
                                null, // net
                                "DUMP.bestSoFarConservative.off",
                                true); // saveHomo (as opposed to heights above paraboloid)
        } // while true
    }
} // class CounterexampleSearch