mardi 13 août 2019

Solving 2008 Google Code Jam Fly Swatter with Monte Carlo Simulation

Looking back over the 2008 Qualifying round Question C, I came up with the idea that this could be solved using a Monte Carlo Simulation. So far my attempts have not succeeded. However my results consistently over-estimate the probablility of hitting the fly. So I suspect it is a systematic issue: Any ideas would be great...

Running over the basic Solution:

  1. Run monte carlo batches of 100000 tests at a time, 5 or more times where the stopper is when the probability of hitting a fly does not change by EPSILON or more (e.g. .0000001) between runs (monte carlo answers are not too sensitive to epsilon within reason, but the run time goes up... I have run many different values of epsilon)

  2. Each case (read the 2008 Qualifying round Question C) is run by generating a random number for radius, and another for "theta" the angle... Determine the following: A: did the Rim of the Racket hit the fly, B: is the fly bigger than the gap, or C: did part of the fly catch part of a string.

please: note the comments contain // [test cases * 100000]: error/delta between runs: average of the random numbers generated for the monte carlo simulation (0-1.0)

Case #1: 1.000000     // 5: 0.0 : drand_avg=0.4999613
Case #2: 0.917519     // 36: 7.857142847456089E-8 : drand_avg=0.49994657
Case #3: 0.000000     // 5: 5.0E-8 : drand_avg=0.4999613373
Case #4: 0.002615     // 26: 1.8461538461595622E-8 : drand_avg=0.499946629
Case #5: 0.575228     // 104: 2.949962663034711E-8 : drand_avg=0.4999600

Comparing this to to the Expected Results from Google:

Case #1: 1.000000
Case #2: 0.910015
Case #3: 0.000000
Case #4: 0.002371
Case #5: 0.573972

Finally here is my Code:

    class PCase extends Case {

        double flyRadius;
        double racketRadius;
        double racketThickness;
        double stringRadius;
        double stringGap;
        int caseNumbert;

        PCase(int caseNumber, 
                double f,
                double R, 
                double t, 
                double r, 
                double g) {
            super(caseNumber);

            //f, R, t, r and g
            flyRadius = f;
            racketRadius = R;
            racketThickness = t;
            stringRadius= r;
            stringGap=g;
            caseNumbert=caseNumber;

        }
        String result = "";
        Double EPSILON = 0.0000001;
        DecimalFormat df2 = new DecimalFormat( "#0.000000" );

        @Override
        protected String caseResult() {
            return " "+result;
        }

        @Override
        protected Case processCase() {

            int hits=0;
            int total = 0;

            double percent = 1;
            double delta = 1;
            int c = 0;
            do{
                ++c;
                for(int i=0;i<1000000;++i){
                    ++total;
                    if(monteCarloTest()){
                        ++hits;
                    }
                }
                double t = (double)hits/(double)total;
                delta = Math.abs(t-percent);
                percent =t;
            }while(delta>EPSILON || c<5);
            result = String.valueOf(df2.format(percent));
            System.out.println("Case #"+caseNumbert+":"+this.caseResult()+"\n\t\t// "+ c +": "+delta 
                    + " : drand_avg="
                    +(totalRand/(double)total) 
                    + " : drand2_avg="
                    +(totalRand2/(double)total) 
            );
            return this;

        }
        Random r = new Random(7777777);
        Random r2 = new Random(123456);
        double totalRand =0;
        double totalRand2 =0;

        private boolean monteCarloTest(){
            //boolean hitRacket = false;

            double drand  = r.nextDouble();
            double drand2 = r2.nextDouble();
            totalRand+=drand;
            totalRand2+=drand2;
            double pathRadius = drand*racketRadius;
            double pathAngle  = drand2*(Math.PI/2.0);
            if(pathRadius>racketRadius - racketThickness-flyRadius){
                return true;
            }else{
                double pitch = 2.0d*stringRadius+stringGap;// - 2*stringRadius;
                if(2*flyRadius>stringGap){
                    return true;
                }
                // Convert to cartesian coordinates and try...
                double x = Math.cos(pathAngle) * pathRadius;
                double y = Math.sin(pathAngle) * pathRadius;

                double minX = x - flyRadius;
                double maxX = x + flyRadius;
                double minY = y - flyRadius;
                double maxY = y + flyRadius;

                int leftX   = (int)(minX/pitch);
                int rightX  = (int)(maxX/pitch);
                int bottomY = (int)(minY/pitch);
                int topY    = (int)(maxY/pitch);
                if(leftX!=rightX  || bottomY!=topY){
                    return true;
                }

                double leftStringEdge 
                        = pitch * leftX + stringRadius;
                if(leftStringEdge>minX){
                    return true;
                }


                double rightStringEdge 
                        = pitch * (leftX+1) - stringRadius;
                assert (rightStringEdge-leftStringEdge -stringGap)<EPSILON;
                if(rightStringEdge<maxX){
                    return true;
                }


                double bottomStringEdge =
                        pitch * bottomY + stringRadius;
                if(bottomStringEdge>minY){
                    return true;

                }


                double topStringEdge =
                        pitch * (bottomY+1) - stringRadius;
                // Confirm the Top String's Bottom Edge is the same via both calcs...
                assert (topStringEdge - bottomStringEdge-stringGap)<EPSILON;
                if(topStringEdge < maxY){
                    return true;
                }


            }
            return false;
        }
    }




Aucun commentaire:

Enregistrer un commentaire