Estimate the probability that the triangle formed by three uniformly random points on the sphere contains its circumcentre; see https://math.stackexchange.com/questions/4842652.

Question4840933.java

public class Question4840933 {
    final static long ntrials = 2000000;
    
    public static void main(String [] args) {
        double [] [] x = new double [3] [3];
        double [] c = new double [3];
        
        long count = 0;
        
        for (long n = 0;n < ntrials;n++) {
            // the first point can be fixed arbitrarily
            x [0] [0] = 0;
            x [0] [1] = 0;
            x [0] [2] = 1;

            // the azimuth of the second point can be fixed arbitrarily
            double cos1 = 1 - 2 * Math.random();
            double sin1 = Math.sqrt (1 - cos1 * cos1);
            x [1] [0] = 0;
            x [1] [1] = sin1;
            x [1] [2] = cos1;

            // the third point has to be in general position
            double cos2 = 1 - 2 * Math.random();
            double sin2 = Math.sqrt (1 - cos2 * cos2);
            double phi = 2 * Math.PI * Math.random();
            double cos = Math.cos(phi);
            double sin = Math.sin(phi);
            x [2] [0] = sin2 * cos;
            x [2] [1] = sin2 * sin;
            x [2] [2] = cos2;
            
            // form the differences from the first point
            x [1] [2]--;
            x [2] [2]--;
            
            // form the cross product of the differences
            for (int i = 0;i < 3;i++)
                c [i] = x [1] [(i + 1) % 3] * x [2] [(i + 2) % 3] - x [2] [(i + 1) % 3] * x [1] [(i + 2) % 3];
            normalize (c);

            // restore the original points
            x [1] [2]++;
            x [2] [2]++;

            // project the points along the cross product
            for (int j = 0;j < 3;j++) {
                double dot = dot (x [j],c);
                for (int i = 0;i < 3;i++)
                    x [j] [i] -= dot * c [i];
                normalize (x [j]);
            }
            
            // check whether one of the angles is reflex
            double angle = 0;
            for (int j = 0;j < 3;j++)
                angle += Math.acos(dot (x [j],x [(j + 1) % 3]));
            if (angle < 2 * Math.PI - 1e-11)
                count++;
        }
        
        System.out.println(count / (double) ntrials);
    }
    
    static void normalize (double [] x) {
        double norm = 0;
        for (int i = 0;i < 3;i++)
            norm += x [i] * x [i];
        norm = 1 / Math.sqrt(norm);
        for (int i = 0;i < 3;i++)
            x [i] *= norm;
    }
    
    static double dot (double [] x,double [] y) {
        double dot = 0;
        for (int i = 0;i < 3;i++)
            dot += x [i] * y [i];
        return dot;
    }
}