Data hosted with ♥ by Pastebin.com - Download Raw - See Original
  1. #Modified from: http://xianblog.wordpress.com/2012/03/17/simulated-annealing-for-sudokus-2
  2.  
  3. #Calculates the error from a proposed solution
  4. target=function(s){
  5.   tar=sum(apply(s,1,duplicated)+apply(s,2,duplicated))
  6.   for (r in 1:9){
  7.     bloa=(1:3)+3*(r-1)%%3
  8.     blob=(1:3)+3*trunc((r-1)/3)
  9.     tar=tar+sum(duplicated(as.vector(s[bloa,blob])))
  10.   }
  11.   return(tar)
  12. }
  13.  
  14.  
  15. cost=function(i,j,s){
  16.   #constraint violations in cell (i,j)
  17.   cos=sum(s[,j]==s[i,j])+sum(s[i,]==s[i,j])
  18.   boxa=3*trunc((i-1)/3)+1;
  19.   boxb=3*trunc((j-1)/3)+1;
  20.   cos+sum(s[boxa:(boxa+2),boxb:(boxb+2)]==s[i,j])
  21. }
  22.  
  23. #Creates a starter suggestion substituting ZEROS by sampled numbers!
  24. entry=function(to_solve){
  25.   s=to_solve
  26.   pop=NULL
  27.   for (i in 1:9)
  28.     pop=c(pop,rep(i,9-sum(to_solve==i)))
  29.   s[s==0]=sample(pop)
  30.   return(s)
  31. }
  32.  
  33. # generates a step following simulated annealing strategy!
  34. move=function(tau,s,con){
  35.   pen=(1:81)
  36.   for (i in pen[con==0]){
  37.     a=((i-1)%%9)+1
  38.     b=trunc((i-1)/9)+1
  39.     pen[i]=cost(a,b,s)
  40.   }
  41.   #It seems that Evolutionary Algorithms or Tabu Search would work better
  42.   #as it is statying a lot of time in the same solution!
  43.   #Maybe simply storing the best from the past and trying go back after stuck!
  44.   wi=sample((1:81)[con==0],2,prob=exp(pen[(1:81)[con==0]]))
  45.   prop=s
  46.   prop[wi[1]]=s[wi[2]]
  47.   prop[wi[2]]=s[wi[1]]
  48.  
  49.   if (target(prop) < target(s))
  50.     s=prop
  51.   else if (runif(1)<exp((target(s)-target(prop)))/tau)
  52.     s=prop
  53.   return(s)
  54. }
  55.  
  56.  
  57. run_annealing=function(num_it, num_try, tau, to_solve, temp_degradation=0.99, draw=FALSE, point_color="black") {
  58.   s=entry(to_solve)
  59.   for(t in 1:num_it) {
  60.     for(v in 1:num_try) s=move(tau,s,to_solve)
  61.     tau=tau*temp_degradation
  62.     if (draw) points(t,target(s),col=point_color, cex=.3,pch=19)
  63.     if(target(s)==0) break()
  64.   }
  65. }
  66.  
  67. #Too easy one
  68. s=matrix(0,ncol=9,nrow=9)
  69. s[1,c(6,8)]=c(6,4)
  70. s[2,c(1:3,8)]=c(2,7,9,5)
  71. s[3,c(2,4,9)]=c(5,8,2)
  72. s[4,3:4]=c(2,6)
  73. s[6,c(3,5,7:9)]=c(1,9,6,7,3)
  74. s[7,c(1,3:4,7)]=c(8,5,2,4)
  75. s[8,c(1,8:9)]=c(3,8,5)
  76. s[9,c(1,7,9)]=c(6,9,1)
  77.  
  78. plot(c(0,0),col="white",xlim=c(1,100),ylim=c(0,target(s)),xlab="iterations",ylab="penalty")
  79. run_annealing(100,100,100,s,0.99,TRUE,"blue")
  80. run_annealing(100,100,50,s,0.99,TRUE,"red")
  81. run_annealing(100,100,100,s,0.999,TRUE,"green")
  82. run_annealing(100,100,50,s,0.999,TRUE,"black")