I believe I'm making progress on the exact problem. I have come up with a flawed algorithm that gives answers that are "close to" correct. The algorithm is for a general m-by-n board (it is easier to test that way). Also I wrote a brute-forcer for smaller boards, to test against. The flawed algorithm is simple (in a manner of speaking), concise, and very fast, but the answers given are always a bit high (except for 2-by-2 board), for example, to your 8-by-8 question the answer given is 0.00030602613...(more decimal places). I'm working on a way to fix it.

Monte Carlo is at about 1.2 * 10^11 trials (!) and the answer seems to be in the range [0.000305835, 0.000305875] again as a rough estimate.

I'll hold off on explaining how it works for the time being; I'm just focused on trying to fix it.

Here is the complete code up to this point:

Code:import java.util.Random; public class ChessBlackSquare22 { static Random g=new Random(); static double g1; public static void main(String[] args) { //monteCarlo(); solve(0.05,3,4); naive(0.05,3,4); } static void solve(double p,int m,int n) { int i,j; long t=time(); double z[][]=new double[m][n],q,r,a,b,c,x=0,y=1; for(i=0;i<m;i++) z[i][0]=p; for(i=1;i<n;i++) z[0][i]=p; for(i=1;i<m;i++) for(j=1;j<n;j++) { a=z[i-1][j-1]; b=z[i-1][j]; c=z[i][j-1]; x+=y*(q=a*b*c*p); z[i][j]=(p-q)/(r=1-q); z[i-1][j]=(b-q)/r; z[i][j-1]=(c-q)/r; y*=r; } System.out.println(x); System.out.println("Elapsed: "+(time()-t)/1000.0+" s"); } static void naive(double p,int m,int n) { long t=time(); if(m*n>25) return; g1=0; naive2(m,n,p,0,0,1,new boolean[m][n]); System.out.println(g1); System.out.println("Elapsed: "+(time()-t)/1000.0+" s"); } static void naive2(int m,int n,double p,int a,int b,double q,boolean[][] z) { int i,j,a2,b2; boolean[][] y; if(a==m) { outer: for(i=0;i<m-1;i++) for(j=0;j<n-1;j++) if(z[i][j]&&z[i][j+1]&&z[i+1][j]&&z[i+1][j+1]) { g1+=q; break outer; } return; } if(b<n-1) { a2=a; b2=b+1; } else { a2=a+1; b2=0; } naive2(m,n,p,a2,b2,q*(1-p),deepcopy(z)); y=deepcopy(z); y[a][b]=true; naive2(m,n,p,a2,b2,q*p,y); } static boolean[][] deepcopy(boolean[][] a) { int i; boolean[][] b=new boolean[a.length][]; for(i=0;i<a.length;i++) b[i]=a[i].clone(); return b; } static void monteCarlo() { int i,j; long s=0,t,t1=time(); boolean b,c[][]; for(t=0;;t++) { if(t>0&&(t%500000)==0) System.out.println("Trials: "+t+ ", Probability: "+(s/(double)t)+", Elapsed:"+((time()-t1)/1000.0)); b=false; c=new boolean[8][8]; for(i=0;i<8;i++) for(j=0;j<8;j++) if(g.nextInt(20)==0) c[i][j]=true; for(i=0;i<7;i++) for(j=0;j<7;j++) if(c[i][j]&&c[i][j+1]&&c[i+1][j]&&c[i+1][j+1]) { b=true; break; } if(b) s++; } } static long time() { return System.currentTimeMillis(); } }