Using R and the library lpSoveAPI to find solutions to FlowFree boards. Note that the side has to be specified. There are 3 problems (*.csv) files, included for testing the testing. Create your own csv for other problems.

flowFree.R

library(lpSolveAPI)
source("~/RStats/FlowFree/flowfree_lib.R")

problemfile = "problem5_33.csv" ; side<- 5
problemfile = "problem.csv"
problemfile = "problem6_4.csv"; side <- 6
problemfile = "problem6_20.csv"; side <- 6
problemfile = "problem8_17.csv"; side <- 8

#load one specific Puzzle
terminal.cells <- read.csv(problemfile, header=T, stringsAsFactors=FALSE)
terminal.cells$tcell <- (terminal.cells$Y -1)* side + terminal.cells$X # cell serial number
num.colors <- length(unique(terminal.cells$color))
str(terminal.cells)
colorpalette <- unique(terminal.cells$palette)

num.colors
init(side, num.colors)
n.row; n.col
#df, const.type.vec have been initialized
df <- populate.Amatrix()
num.cells; num.edges; num.colors; length(rhs.vector); length(const.type.vec); num.nt
rhs.vector <- create_rhs_vector(rhs.vector, terminal.cells)
const.type.vec <- createConstraintTypeVector(const.type.vec)
length(rhs.vector); length(const.type.vec)
dim(df)

# actual problem definition
lpff <- make.lp(nrow=n.row, ncol=n.col)
defineIP()
lpff
solve(lpff)
sol <- get.variables(lpff)
#sol #take a quick look at the solution
sum(unlist(sol[1:num.edges]))
######################################################
### Done with the IP, now on to plotting
colorpalette
plotSol(terminal.cells, colorpalette)
terminal.cells

Flowfree_lib.R

rm(list=ls())
library(lpSolveAPI)

#utility functions
getColnum <- function(cell, color, edge) {
#  return(4*(cell-1)+edge)
  celloffset <- (cell-1)*4*num.colors
  coloroffset <- (color-1)*4
  return(celloffset+coloroffset+edge)
}

#Important: Even the non-terminal cells are included here, with all 0's.
#This is to avoid ugly if nt type checks.
#Thus the total number of Ycols = Num.Cell x Num.Colors
getYColnum <- function(cell, k, terminal.cells) {
  skip <- num.edges
  celloffset <- (cell-1)*num.colors
  coloroffset <- k
  return(skip+celloffset+coloroffset)
}



getColnumFromXYE <- function(xpos,ypos, color, edge) {
  cellnum <- side*(ypos-1) + xpos
  celloffset <- (cellnum-1)*4*num.colors
  coloroffset <- (color-1)*4
  return(celloffset+coloroffset+edge)  
}

getCellandEdgeGiveColnum <- function(col) {
  cell <- ((col-1) %/% (num.colors*4)) + 1 #serial number of the cell
  k <- ((col-1) %% (num.colors*4) %/% 4) + 1 # color value of the column
  edge <- col %% 4 #get the edge value, 1 through 4
  if (edge==0) {edge=4}
  return(list(cell, k, edge))
}



getSolutionGivenCell <- function(cell, sol) {
  skip <- (cell-1)*4*num.colors  
  range <- (skip+1):(skip+4*num.colors)
  cell.sol <- sol[range]
  print(sum(cell.sol))
  return(cell.sol)
}


IsCellNonTerminal <- function(c, terminal.cells) {
  for (i in 1:nrow(terminal.cells)) {
    if(terminal.cells$tcell[i] == c) {
      return(0)
    }
  }
  return(1) #default
}


### End of Utility Functions

createXcolumnDF <- function() {
  cells <- 1:num.cells
  colors <- 1:num.colors
  x <- expand.grid(1:4, colors, cells)
  x.col <- data.frame(1:num.edges, x[c(3,2,1)]) #re-arranging
  names(x.col) <- c("col","cell","color","edge")
  x.col #return 
}



init <- function(side, num.colors) {  
  num.cells <<- side*side
  num.edges <<- num.cells* 4 * num.colors
  num.nt <<- num.cells - (2*num.colors) #number of non-terminal cells
  num.horiz <<- side * (side-1) * num.colors
  num.vert <<- side * (side-1)  * num.colors
  num.limitedge.constraints <<- num.cells * 4
  num.boundary.constraints <<- side * 4 * num.colors #first row, top row, first col, last col
  num.noncolor.constraints <<- num.colors * 2 #one for each terminal of the segment of that color
  num.pick1 <<- num.nt
  num.samecolor <<- (num.nt * num.colors)
  n.row <<- num.cells + num.limitedge.constraints + num.horiz + num.vert + num.boundary.constraints + num.noncolor.constraints + num.pick1 + num.samecolor
  n.col <<- num.edges + (num.cells * num.colors) # Xcke + Yck variables
  
  start.pickone <<- num.cells + num.limitedge.constraints + num.vert + num.horiz + num.boundary.constraints + num.noncolor.constraints + 1
  start.samecolor <<- start.pickone + num.pick1  
  
  #Ax <= B
  const.type.vec <<- rep("=", n.row)
  rhs.vector <<- rep(0, n.row)
  x.col <<- createXcolumnDF()  
}


defineIP <- function() {
  set.objfn(lpff, rep(1, n.col))
  set.constr.type(lpff, const.type.vec, 1:n.row) #horiz, vert, corners
  set.rhs(lpff, b=rhs.vector, constraints=1:n.row) # assign rhs values   
  #Set all the columns at once
  for (col in 1:n.col) {
    set.column(lpff, col, df[ ,col])  
    set.type(lpff, col, "binary")      
  }
  #assemble it all
  dimnames(lpff) <- setRowAndColNames()
  write.lp(lpff, "flowfreeIP.lp", "lp")#write it out
  
}


# Order of Constraints
#Cover (1 for each cell)
#Limit Each edge to at most 1 (1 for each edge)
#Horizontal Connectivity
#Vert Connectivity
#Boundary Removal (boundary edges are set to 0)
#Set Non-colors in TERMINAL NODES to 0 (One constraint for each terminal node)
populate.Amatrix <- function() {  
  df <- data.frame(matrix(0, n.row, n.col))  #create a dataframe shell  
  #Every Cell has 2 nz edges out of the 4 possible ones
  for (row in 1:num.cells) {
        skip = (row-1)*num.colors*4
        range <- 1:(num.colors*4)
        df[row, skip+range] <- 1
  }
  print("Done with cover Constraints")
  ##############
  df <- limit_edge_flow_to_atmost_one(df)
  print("Done with LimitEdge Constraints")
  
  
  # Horizontal connectivity constratints
  start.constraint.num <- num.cells + num.limitedge.constraints
  crow <- start.constraint.num
  for (row in 1:side) {
    for (link in 1:(side-1)){
      for ( k in 1:num.colors) {      
        # x13 = x21; x23 = x31; ...; x43=x51
        cell <- (row-1)*side + link 
        cnum <-  x.col[ x.col$cell==cell & x.col$color==k & x.col$edge==3  , 1] #lookup       
        cnum2 <-  x.col[ x.col$cell==(cell+1) & x.col$color==k & x.col$edge==1  , 1]  #lookup               
        #print(paste(cnum, cnum2))
        crow <- crow + 1
        df[crow, cnum] <- 1
        df[crow, cnum2] <- -1
      }
    }
  }
  print("Done with Horiz df")
  
  # Vertical connectivity constratints
  start.constraint.num <- num.cells + num.limitedge.constraints + num.horiz  
  crow <- start.constraint.num
  for (row in 1:(side-1)) {
    for (link in 1:(side)) {
      for ( k in 1:num.colors) {      
        bottom.cell <- (row-1)*side + link 
        top.cell <- (row)*side + link 
        #cnum <- getColnum(bottom.cell,k,2)
        #cnum2 <- getColnum(top.cell,k,4)
        cnum <-  x.col[ x.col$cell==bottom.cell & x.col$color==k & x.col$edge==2  , 1] #lookup column        
        cnum2 <-  x.col[ x.col$cell==top.cell & x.col$color==k & x.col$edge==4  , 1]  #lookup      
        
        crow <- crow + 1
        #print(paste(cnum, cnum2, crow))
        df[crow, cnum] <- 1
        df[crow, cnum2] <- (-1)
      }
    }
  }
  print("Done with Vert df")
  
  df <- boundary.edge.removal(df) 
  df <- set_other_colors_in_terminals_to_zero(terminal.cells, df)   
  df <- pick_one_color_for_nonterminals_cells(terminal.cells,df)
  print(paste("Dimension", dim(df)))
  df <- single_color_per_cell_for_nonterminals_cells(terminal.cells,df)
  print(paste("Dimension", dim(df)))
  return(df) # updated a-matrix
}

##################################


#We don't want the same edge to have two colors.
# Limit that by taking each edge and forcing it to be atmost 1. <= constraint
# Total number of constraints = number of edges = number of cells*4
limit_edge_flow_to_atmost_one<- function(df) {  
  crow <- num.cells
  for (c in 1:num.cells){
      #cnum <-  x.col[ x.col$cell==c & x.col$color==1 & x.col$edge==1, 1]
      cnum<- ((c-1)*num.colors*4)+1
      end <- cnum + (num.colors-1)*4
      range <- seq(cnum, end, by=4) #every 4th column is for the same edge, different color
      df[crow+1, range]   <- 1  # limit cell C, edge 1
      df[crow+2, range+1] <- 1  # limit cell C, edge 2
      df[crow+3, range+2] <- 1  # limit cell C, edge 3
      df[crow+4, range+3] <- 1  # limit cell C, edge 4
      crow <- crow+4            #jump ahead to next 4 rows
    }
    return(df)
}




#########
# boundary Edge removals - 4 * side number of them
#find the correct column numbers and make the rhs  = 0 
boundary.edge.removal <- function(df) {
  crow <- num.cells + num.limitedge.constraints + num.vert + num.horiz
  for(index in 1:side) {
    for ( k in 1:num.colors) {      
      aa <- getColnumFromXYE(index, 1, k,4) #1st row, 4th edges to be removed
      crow <- crow+1; df[crow,aa] <- 1
      bb <- getColnumFromXYE(index, side, k,2) #top row, 2nd edges to be removed
      crow <- crow+1; df[crow,bb] <- 1
      cc <- getColnumFromXYE(1, index, k,1) #first column, 1st edges to be removed
      crow <- crow+1; df[crow,cc] <- 1
      dd <- getColnumFromXYE(side, index, k, 3) #last column, 3rd edges to be removed
      crow <- crow+1; df[crow,dd] <- 1
      #print(paste(aa,bb,cc,dd))
    }  
  }
  return(df)
}

# In Terminal nodes, only one edge can be 1.
# That edge color has to be the same as the terminal color.
#everything else is set to 0
set_other_colors_in_terminals_to_zero <- function(terminal.cells, df) {  
  crow <- num.cells + num.limitedge.constraints + num.vert + num.horiz + num.boundary.constraints
  for (i in 1:nrow(terminal.cells)) {
    cellnum <- terminal.cells[i,1] + (terminal.cells[i,2]-1)*side
    crow <- crow+1
    for(k in 1:num.colors) {
      allowed.color <- terminal.cells[i,3] # get hold of its color, so that it can be spared. Everything else is set to be 0
      if(k != allowed.color){
        for (e in 1:4) {
          cnum <-  x.col[ x.col$cell==cellnum & x.col$color==k & x.col$edge==e  , 1]        
          #cnum <- getColnum(cellnum, k, e)
          df[crow, cnum] <- 1  
        }
      }
    }
  }
  return(df)
}

  



#This constraints ensures that non-terminal cells do not take on multiple colors.
#We introduce a new type of coloumn Y_cell_color. and set only one of them to be 1
pick_one_color_for_nonterminals_cells <- function(terminal.cells, df) {  
  crow <- start.pickone - 1
  for (c in 1:num.cells) {
    if(IsCellNonTerminal(c, terminal.cells)) {
      crow <- crow + 1 
      for(k in 1:num.colors) { 
        cnum <- getYColnum(c, k, terminal.cells)
        df[crow, cnum] <- 1  
        rhs.vector[crow] <<- 1
      }
    }
  }
  return(df)
}


#for any given cell, there can only be one color that is non-zero
# the sum of those colored-edges should total 2.
#We equate it to 2 * Y_ck and pick one of the colors
single_color_per_cell_for_nonterminals_cells <- function(terminal.cells, df) {  
  crow <- start.samecolor - 1
  for (c in 1:num.cells) {
    if(IsCellNonTerminal(c, terminal.cells)) {
      #print(paste("terminal cell", c))
      for(k in 1:num.colors) { 
        crow <- crow + 1 
        for (e in 1:4) {
          #cnum <- getColnum(c, k, e)
          cnum <-  x.col[ x.col$cell==c & x.col$color==k & x.col$edge==e  , 1]        
          
          df[crow, cnum] <- 1  
        }        
        cnum <- getYColnum(c, k, terminal.cells)
        df[crow, cnum] <- (-2)  
      }
    }
  }
  return(df)
}

### CONSTRAINT TYPES
createConstraintTypeVector <- function(const.type.vec) {
  #all equality constraints by default; set in init() function
  #overwrite the const.type vector, for the LimitEdge to 1 constraints. 
  # it is a less than or equal to
  const.type.vec[(num.cells+1):(num.cells+num.limitedge.constraints)] <- "<="  
  const.type.vec
}


# # # # # #
#  R H S # 
# # # # # 


#Problem specific adjustments to RHS
# Terminal cells have only one edge as opposed to all other edges that have 2 that are on
create_rhs_vector <- function(rhs.vector, terminal.cells) {
  rhs.vector[1:num.cells] <- 2 #the first n^2 rows are set to 2
  for (i in 1:nrow(terminal.cells)) {
    cellnum <- terminal.cells[i,1] + (terminal.cells[i,2]-1)*side
    rhs.vector[cellnum] <- 1 #Cover constraint for TerminalNodes is 1 (not 2)
  }
  rhs.vector[(num.cells+1):(num.cells*5)] <- 1 #limitedge to one color
  return(rhs.vector)
  #  rhs.vector[(num.cells+1):(num.cells+num.limitedge.constraints)] <- 1 #Limit each edge to be at most 1
  
}

  
getRowNames <- function() {  
  ######### Horizontal Connectivity
  row.cover.names<- paste("CoverCell", 1:num.cells, sep="_") #Cover constraints
  limitedge.row.names <- paste("LimitEdge_", 1:num.limitedge.constraints, sep="") # specify  row names
  
  start.constraint.num <- num.cells+num.limitedge.constraints+1 #Horiz start
  constraint.set <- seq(start.constraint.num, start.constraint.num+num.horiz-1)
  rowh.names <- paste("Horiz_", constraint.set, sep="") 
  
  # Vertical connectivity constratints
  start.constraint.num <- num.cells+num.horiz+1
  constraint.set <- seq(start.constraint.num, start.constraint.num + num.vert-1)
  rowv.names <- paste("Vert_", constraint.set, sep="") 
  
  boundary.row.names <- paste("Boundary_", 1:num.boundary.constraints, sep="") # specify  row names
  noncolor.row.names <- paste("Noncolor_", 1:num.noncolor.constraints, sep="") # specify  row names
  pick1.row.names <- paste("PickOneNT_", 1:num.nt, sep="") # One for each NT row
  samecolor.row.names <- paste("SameColorNT_", 1:(num.nt*num.colors), sep="") # One for each NT row and color combination
  
  row.names <- c(row.cover.names, limitedge.row.names, rowh.names, rowv.names, boundary.row.names, noncolor.row.names, pick1.row.names, samecolor.row.names)
  return(row.names)
}

#set rownames & colnames
setRowAndColNames<- function() {
  abc<- t(outer(1:num.cells, 1:num.colors, paste, sep="_"))
  abcd <- paste(abc, sep="")
  def <- t(outer(abcd, 1:4, paste, sep="-"))
  c.names<- paste("x_",def, sep="")
  ycolumn.names <- paste("y_", abcd, sep="")
  col.names <- c(c.names, ycolumn.names)
  
  row.names <- getRowNames() 
  return( list(row.names, col.names)  )
}

###################






### Plotting Related Functions

# First, get the x,y coords of the CENTER of each cell.
getCellCentersXY <- function() {
  x <- 1:side
  beat<- 2*x-1
  cx <- rep(beat, side)
  cy <- unlist(lapply(beat, rep, side))  
  return(list(cx,cy))  
}



segmentForEdge <- function(cell, edge, centers) {
  x <- centers[[1]][cell]
  y <- centers[[2]][cell]
  xoffset <- (switch(edge, -1,0,1,0)) #using switch instead of multiple if's
  yoffset <- (switch(edge, 0,1,0,-1)) #using switch instead of multiple if's
  xe <- x + xoffset
  ye <- y + yoffset
  return(list(x,y,xe,ye))
}

getSegmentsgivenSol <- function(sol, centers){
  nz <- sum(unlist(sol))
  coords <- data.frame(matrix(0,nz,6)) #4 segment coords, cellnum and color
  r <- 0
  for (i in 1:length(sol)) {
    if(sol[i]==1){
      r <- r+1
      cell_edge <- getCellandEdgeGiveColnum(i) #returns a tuple (cell, k , edge)
      s <- segmentForEdge(cell_edge[[1]], cell_edge[[3]], centers) # we now have the 4 segment coords
      coords[r,5] <- cell_edge[[1]] #cellnum
      coords[r,6] <- cell_edge[[2]] #color
      for (j in 1:4) {
        coords[r,j] <- s[[j]] # add to the dataframe in the right places
      }
    }
  }    
  return(coords)
}



getGridlines <- function() {
  xs<- rep(0, side-1)
  ys <- seq(2,2*(side-1), by=2)
  xe<- rep(2*side, side-1)
  ye <- seq(2,2*(side-1), by=2)
  df1<- data.frame(xs,ys,xe,ye)  
  ys<- rep(0, side-1)
  ye<- rep(2*side, side-1)
  xs <- seq(2,2*(side-1), by=2)
  df2<- data.frame(xs, ys, xs, ye)
  names(df1) <- c("xs","ys", "xe", "ye")
  names(df2) <- c("xs","ys", "xe", "ye")
  df1 <- rbind(df1,df2)
  df1
}


plotSol <- function(terminal.cells,colorpalette) {  
  library(ggplot2)
  p <- NULL
  centers <- getCellCentersXY()
  xsol <- sol[1:num.edges]
  coords <- getSegmentsgivenSol(xsol, centers)
  p.df <- data.frame(centers[[1]], centers[[2]])
  p.df$cell <- 1:(side*side)
  term.df <- merge(p.df, terminal.cells[c(3,5)], by.x=c("cell"), by.y=c("tcell"))    
  grid.df <- getGridlines()
  names(p.df) <- c("cx", "cy", "cell")
  names(coords) <- c("x", "y", "xend", "yend", "cell", "color")
  names(term.df) <- c("cell", "cx", "cy", "color")
  p <- ggplot(coords) #+ geom_segment(aes(x, y, xend=xend, yend=yend, color=factor(color)), size=4)
  p <- p + geom_rect(aes(xmin=0, xmax=2*side, ymin=0, ymax=2*side), fill="gray10")  # box bg
  p <- p + geom_segment(data=grid.df, aes(x=xs, y=ys, xend=xe, yend=ye), color="gray50", size=2)
  p <- p + geom_segment(aes(x, y, xend=xend, yend=yend, color=factor(color)), size=4) #actual solution
  p <- p + geom_point(data=p.df, aes(x=cx, y=cy), color="grey50", size=4) #cell centers
  p <- p + geom_point(data=term.df, aes(x=cx, y=cy, color=factor(color)), size=8) + scale_colour_manual(values = colorpalette)
  p <- p + guides(color=FALSE) + theme(panel.background=element_rect(fill="white"), panel.grid.minor=element_blank())
  p <- p + scale_x_continuous(breaks=seq(0, 2*side, 2)) + scale_y_continuous(breaks=seq(0, 2*side, 2))
  return(p)
}

plotProb <- function(terminal.cells,colorpalette) {  
  p <- NULL
  centers <- getCellCentersXY()
  p.df <- data.frame(centers[[1]], centers[[2]])
  p.df$cell <- 1:(side*side)
  term.df <- merge(p.df, terminal.cells[c(3,5)], by.x=c("cell"), by.y=c("tcell"))    
  grid.df <- getGridlines()
  names(p.df) <- c("cx", "cy", "cell")
  names(term.df) <- c("cell", "cx", "cy", "color")
  p <- ggplot() #+ geom_segment(aes(x, y, xend=xend, yend=yend, color=factor(color)), size=4)
  p <- p + geom_rect(aes(xmin=0, xmax=2*side, ymin=0, ymax=2*side), fill="gray10")  # box bg
  p <- p + geom_segment(data=grid.df, aes(x=xs, y=ys, xend=xe, yend=ye), color="gray50", size=2)
  p <- p + geom_point(data=term.df, aes(x=cx, y=cy, color=factor(color)), size=12) + scale_colour_manual(values = colorpalette) # terminal nodes
  p <- p + guides(color=FALSE) + theme(panel.background=element_rect(fill="white"), panel.grid.minor=element_blank())
  p <- p + scale_x_continuous(breaks=seq(0, 2*side, 2)) + scale_y_continuous(breaks=seq(0, 2*side, 2))
  return(p)
}

problem5_2.csv

problem5_33.csv

problem6_20.csv