R/write_MPSILP.R

In DoE.MIParray: Creation of Arrays by Mixed Integer Programming

write_MPSILP <- function (file, qco=NULL, qcoinfo=NULL, cvec=NULL, bvec=NULL, Amat=NULL, bsense=NULL, boundmat=NULL, name = "EnforceResolution", commentline="* no comment"){
  if (!is.character(file)) stop("file must be a character string that provides the storage location for the MPS file")
  
  if (!is.null(qco)){
    ## cvec, bvec, Amat, boundmat and commentline are ignored
    ## qco is either a qco object or a list produced by create_ILP_list
    stopifnot("qco" %in% class(qco))
    ## a single problem
    cvec <- qco$c
    ## needed for incorporating forced runs
    ## "G" occurs for distinct=FALSE
    bsense <- rep("E", ncol(qco$bc))
    bsense[qco$bc[2,]>qco$bc[1,]] <- "G"
    bvec <- qco$bc[1,]
    Amat <- qco$A
    boundmat <- qco$bx
  }
  if (!is.null(qcoinfo)){
    ## qcoinfo allows to provide the tailored comment line
    stopifnot(is.list(qcoinfo))
    
    commentline <- "* "
    if (!is.null(qcoinfo$nlev)) commentline <- paste0(commentline, "levels: ", paste(qcoinfo$nlev, collapse=", "),"; ")
    if (!is.null(qcoinfo$nruns)) commentline <- paste0(commentline, qcoinfo$nruns, " runs; ")
    if (!is.null(qcoinfo$nruns)) commentline <- paste0(commentline, "requested resolution: ", qcoinfo$reso)
  }

  if (any(is.null(cvec), is.null(bvec), is.null(Amat))){
    stop("Essential specifications are missing")
  }

  ## adapted from package linprog
  nCon <- length(bvec)
  nVar <- length(cvec)
  
  stopifnot(is.numeric(cvec), is.numeric(bvec))
  stopifnot(is.numeric(as.matrix(Amat)))
  if (!is.null(boundmat)){
    stopifnot(is.matrix(boundmat))
    if (!nrow(boundmat)==2) stop("boundmat must have two rows")
    stopifnot(ncol(boundmat)==nVar)
  }
  else boundmat <- matrix(c(0,1), 2, nVar, byrow=FALSE)
  
  if (!is.null(bsense)){
    stopifnot(is.character(bsense))
    stopifnot(length(bsense)==nCon)
    stopifnot(all(bsense %in% c("E","G","L","N")))
  }
    
  
  if (is.null(names(bvec))) {
    blab <- paste0("R_", 1:nCon)
  }
  else {
    blab <- names(bvec)
    blab <- substr(gsub(" ", "", blab), 1, 8)
    j <- 2
    if (length(unique(blab)) < nCon){
      for (i in 1:nCon)
        while (i > 1 & blab[i] %in% blab[1:(i - 1)]) {
          blab[i] <- paste(substr(blab[i], 1, 7 - nchar(as.character(j))), 
                           "_", as.character(j), sep = "")
          j <- j + 1
        }
    }
  }
  if (is.null(names(cvec))) {
    clab <- rep("", nVar)
    clab <- paste0("C_", 1:nVar)
  }
  else {
    clab <- gsub(" ", "", names(cvec))
    clab <- substr(clab, 1,8)
    j <- 2
    if (length(unique(clab)) < nVar){
      for (i in 1:nVar)
        while (i > 1 & clab[i] %in% clab[1:(i - 1)]) {
          clab[i] <- paste(substr(clab[i], 1, 7 - nchar(as.character(j))), 
                           "_", as.character(j), sep = "")
          j <- j + 1
        }
    }
  }
  nc <- nchar(clab) ## for column adjustment
  cv <- as.character(signif(cvec, 10))
  ncv <- nchar(cv)
  nb <- nchar(blab) ## for column adjustment
  bv <- as.character(signif(bvec, 10))
  nbv <- nchar(bv)
  Amatv <- matrix(as.character(signif(Amat, 10)), nCon, nVar)
  ncA <- nchar(Amatv)
  message("start writing ...")
  cat(paste("NAME          ", name, sep = ""), file=file, fill=TRUE)
  cat(commentline, file=file, append=TRUE, fill=TRUE)
  cat("OBJSENSE", file=file, append=TRUE, fill=TRUE)
  cat("    MIN", file=file, append=TRUE, fill=TRUE)
  message("writing ROWS ...")
  cat("ROWS", file=file, append = TRUE, fill=TRUE)
  cat(" N  obj", file=file, append = TRUE, fill=TRUE)
  for (i in 1:nCon) {
    cat(paste0(" ", bsense[i], "  ", blab[i]), file=file, 
          append = TRUE, fill=TRUE)
  }
  message("writing COLUMNS ...")
  cat("COLUMNS", file=file, append = TRUE, fill=TRUE)
  line <- "    COUNTS    'MARKER'                 'INTORG'"
  cat(line, file=file, append = TRUE, fill=TRUE)
  for (i in 1:nVar) {
    line <- paste("    ", clab[i], sep = "")
    line <- paste(line, paste(rep(" ", 10 - nc[i]), 
                              collapse = ""), "obj", sep = "")
    ## -3 for nchar("obj")
    line <- paste(line, paste(rep(" ", 22 - 3 - ncv[i]), collapse = ""), 
                  cv[i], sep = "")
    cat(line, file=file, append = TRUE, fill=TRUE)
    
    for (j in 1:nCon) {
      if (Amat[j, i] != 0) {
        line <- paste0("    ", clab[i], sep = "")
        line <- paste0(line, paste(rep(" ", 10 - nc[i]), collapse = ""), blab[j], 
                       sep = "")
        line <- paste0(line, paste(rep(" ", 22 - 
                                         nb[j] - ncA[j,i]), collapse = ""), Amatv[j, i], sep = "")
        cat(line, file=file, append = TRUE, fill=TRUE)
      }
    }
  }
  line <- "    COUNTS    'MARKER'                 'INTEND'"
  cat(line, file=file, append = TRUE, fill=TRUE)
  ## entire COLUMNS group enclosed in markers
  
  message("writing RHS ...")
  cat("RHS", file=file, append = TRUE, fill=TRUE)
  for (i in 1:nCon) {
    line <- paste("    RHS       ", blab[i], sep = "")
    line <- paste(line, paste(rep(" ", 22 - nb[i] - 
                                    nbv[i]), collapse = ""), 
                  bv[i], sep = "")
    cat(line, file=file, append = TRUE, fill=TRUE)
  }
  ### defaults for upper integer bounds (UI) are solver-specific (e.g. 1 for gurobi)
  ### Therefore, bounds are always written.
  ### Infinity bounds are currently replaced by nruns. Is that wise?
  message("writing BOUNDS ...")
  cat("BOUNDS", file=file, append = TRUE, fill=TRUE)
  for (i in 1:nVar){ 
    cat(paste0(" ", ifelse(boundmat[2,i]<Inf, "BV", "PL"), " BND       ", 
                 clab[i]), file=file, append=TRUE, fill=TRUE)
  }
  cat("ENDATA", file=file, append = TRUE, fill=TRUE)
}