**A** rather different Le Monde mathematical puzzle:

A two-person game is played on an nxn grid filled with zero. Each player pick an entry, which is increased by one as well as all adjacent entries. The game stops when all entries are equal. For n=3,4,5, what are the possible grids with identical values all over?

**I**f I define an R neighbourhood function

nighbo=function(i,n){ neigh=i if (i%%n!=1) neigh=c(i-1,neigh) if (i%%n>0) neigh=c(i+1,neigh) if (i+n<=n*n) neigh=c(i+n,neigh) if (i-n>0) neigh=c(i-n,neigh) return(neigh)}

and try a brute force filling of the grid

while ((min(grid)==0)||(length(unique(grid))>1)){ ent=sample(1:(n*n),1,prob=1/as.vector(grid+1)^10) grid[nighbo(ent,n)]=grid[nighbo(ent,n)]+1}

the loop never stops. When thinking of the case n=3 [while running in the early morning], I wondered whether or not reaching an equal value on all entries was at all possible. Indeed, it is impossible to update one of the four corners without updating at least one of the neighbours, while the converse is false. Experimenting further with simulated annealing to optimise the probabilities of picking different entries in the table when n=4,5 seems to indicate this is also the case for larger values of n, in that all attempts lead to larger values of neighbours to the four corners :

outer=c(1,n,n*n,n*n-n+1) border=sort(unique(c(2:(n-1),(n*n-n+2):(n*n-1),1+n*(1:(n-2)),n+n*(1:(n-2))))) inner=(1:(n*n))[-c(outer,border)] # target=function(weit){ grid=matrix(0,n,n) for (t in 1:1e4){ cas=sample(1:3,1,prob=weit) if (cas==1) ent=sample(outer,1,prob=max(grid[outer])-grid[outer]+1) if (cas==2) ent=sample(border,1,prob=max(grid[border])-grid[border]+1) if (cas==3) ent=sample(inner,1,prob=max(grid[inner])-grid[inner]+1) ent=nighbo(ent,n) grid[ent]=grid[ent]+1} ave=c(mean(grid[outer]),mean(grid[border]),mean(grid[inner])) return(list(dive=max(diff(sort(ave))),ave=ave,ava=mean(grid))) } # weit=rep(1,3) cur=target(weit) T=100 while (cur$dive>0){ ind=sample(1:3,1,prob=1/cur$ave) peit=weit peit[ind]=weit[ind]-(cur$ave[ind]-cur$ava)*runif(1)/(cur$ava) while(peit[ind]<0) peit[ind]=weit[ind]-(cur$ave[ind]-cur$ava)*runif(1)/(cur$ava) prop=target(peit) if (log(runif(1))*1e4/T<prop$dive-cur$dive){ weit=peit;cur=prop} T=T*1.00001 print(cur$ave)}