Мне нужно было решить некоторые двоичные задачи линейного программирования, матрицы ограничений которых были очень велики (до такой степени, что они даже не могли быть сгенерированы); поэтому я подумал: почему бы не сделать разреженные матрицы ограничений? Однако ни один из известных мне LP-решателей (lpSolve, Rsymphony), по-видимому, не мог обрабатывать матрицы разреженных ограничений.
Итак, я нашел и установил Rglpk, который обрабатывает разреженные матрицы через slam.
Затем я смог создать свои матрицы, используя xtabs с опцией sparse , преобразовать их в формат slam и отправить их в решатель. Я был впечатлен улучшением производительности при использовании Rglpk в отличие от lpSolve (который также был медленнее, чем Rsymphony при некоторых других проблемах). И в Rglpk я мог бы также использовать опцию presolve , которую lpSolve тоже имел, но фактически игнорировал.
Однако теперь у меня есть два сомнения.
Когда я решил большую проблему с Rglpk, я сначала попытался с presolve=TRUE, и это закончилось менее чем за 20 минут.
Но когда я попробовал то же самое с presolve=FALSE, не только он еще не закончился через несколько часов, но установив verbose=TRUE, я увидел, что целевая функция достигла еще лучших значений по сравнению с «оптимальным» значением, найденным presolve=TRUE пробег.
Q1: не следует ли presolve=TRUE просто ускорить процесс, фактически не изменяя проблему таким образом, чтобы получить менее оптимальное решение?
при поиске дополнительной информации по Rglpk я обнаружил этот отчет , в котором сравнивалось несколько решателей LP.
В докладе, похоже, делается вывод о том, что CLP является лучшим решателем с открытым исходным кодом. Я искал R реализации, но я смог найти только clpAPI ссылку , которая, похоже, не имеет команды высокого уровня, как Rglpk, Rsymphony и lpSolve do. Похоже, что он может обрабатывать разреженные матрицы: команда loadMatrixCLP, кажется, принимает индексы строк и столбцов + значения для ненулевых элементов; но ясно, что это потребует некоторой работы. И тогда я все равно не знаю, как обращаться с различными низкоуровневыми решателями, опциями и т. Д. И я не могу найти никаких ссылок на методы предварительного решения.
Q2: кто-нибудь знает о реализации высокого уровня CLP в R?
Rsymphony должно быть довольно хорошим (COIN-OR и т. Д.), И в его справочном руководстве дано краткое упоминание о разреженных матрицах, но тогда я не могу найти четкого объяснения того, как разреженные матрицы могут быть переданы в решатель.
Я проконсультировался с этой веб-страницей , но не смог найти ничего, что поразило цель.
Есть идеи?
Спасибо!
РЕДАКТИРОВАТЬ: пример проблемы, решаемой с помощью Rsymphony, lpSolve и Rglpk
После дальнейших поисков по совету Саши выяснилось, что Rsymphony может обрабатывать разреженные матрицы.
Поэтому я попробовал это и обнаружил, что это сработало; однако lpSolve нашел дополнительное решение, которое Rsymphony не смогло.
В приведенном ниже примере начальный data.frame t0 представляет 7 объектов (идентифицируемых свойством 'ID'), каждый из которых связан со свойством 'SC' (уникальным для каждого объекта) и свойством массива 'FP'.
Поскольку я не знаю, как правильно обрабатывать свойства массива в R, я распаковал исходные данные отдельно для свойств 'SC' и 'FP' и пометил соответствующие строки, используя свойство 'prop_to_handle'.
Целью оптимизации является выбор 1 объекта с SC == "A", 2 объектов с SC == "B" и 1 объекта с SC == "C", чтобы число уникальных FP в решение максимально.
Я обсуждал эту тему (считая уникальные возможности) в других статьях StackExchange; до сих пор никто не придумал альтернативных решений тому, который я здесь реализую, который состоит в добавлении 1 двоичной переменной для каждой функции FP, а затем соответствующих ограничений, обеспечивающих отсутствие противоречий.
В этом случае имеется 7 идентификаторов и 12 различных функций FP, поэтому необходимо назначить 19 двоичных переменных.
#Initial data
t0 <- structure(list(ID = c(1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 4,
4, 4, 4, 5, 5, 5, 5, 6, 6, 6, 6, 7, 7, 7, 7), prop_to_handle = c("FP",
"FP", "FP", "SC", "FP", "FP", "FP", "SC", "FP", "FP", "FP", "SC",
"FP", "FP", "FP", "SC", "FP", "FP", "FP", "SC", "FP", "FP", "FP",
"SC", "FP", "FP", "FP", "SC"), SC = c(NA, NA, NA, "A", NA, NA,
NA, "A", NA, NA, NA, "B", NA, NA, NA, "B", NA, NA, NA, "B", NA,
NA, NA, "C", NA, NA, NA, "C"), N = c(NA, NA, NA, 1, NA, NA, NA,
1, NA, NA, NA, 2, NA, NA, NA, 2, NA, NA, NA, 2, NA, NA, NA, 1,
NA, NA, NA, 1), FP = c(100, 200, 300, NA, 100, 400, 600, NA,
500, 200, 300, NA, 500, 400, 600, NA, 500, 900, 700, NA, 250,
150, 300, NA, 250, 175, 350, NA)), .Names = c("ID", "prop_to_handle",
"SC", "N", "FP"), row.names = c(NA, -28L), class = "data.frame")
#Installing and loading packages
if (length(find.package(package="Rsymphony",quiet=TRUE))==0) install.packages("Rsymphony")
if (length(find.package(package="slam",quiet=TRUE))==0) install.packages("slam")
if (length(find.package(package="lpSolve",quiet=TRUE))==0) install.packages("lpSolve")
if (length(find.package(package="Rglpk",quiet=TRUE))==0) install.packages("Rglpk")
if (length(find.package(package="slam",quiet=TRUE))==0) install.packages("slam")
require(Rsymphony)
require(slam)
require(lpSolve)
require(Rglpk)
#making sure SC is not a factor
t0["SC"] <- as.character(t0$SC)
#separating the SC and FP handling
t0_FP <- t0[t0$prop_to_handle=="FP",]
t0_SC <- t0[t0$prop_to_handle=="SC",]
#sparse constraint matrix for SC
SC.ID.mat <- xtabs(~SC+ID,t0_SC,sparse=T)
#sparse constraint matrix for FP
FP.ID.mat <- xtabs(~FP+ID,t0_FP,sparse=T)
NIDs_vs_FP <- table(t0_FP[["FP"]])
N.FP <- length(NIDs_vs_FP)
N.ID <- length(dimnames(SC.ID.mat)[[2]])
N.SC <- length(dimnames(SC.ID.mat)[[1]])
N.ID.SC <- as.vector(xtabs(N~SC,t0_SC[!duplicated(t0_SC$SC),]))
NIDs_vs_SC <- table(t0_SC[["SC"]])
N.ID.SC[(unname(NIDs_vs_SC - N.ID.SC) < 0)] <- NIDs_vs_SC[(unname(NIDs_vs_SC - N.ID.SC) < 0)]
#converting the sparse matrices from xtabs to the format accepted by Rglpk
FP.ID.mat <- as.simple_triplet_matrix(FP.ID.mat)
SC.ID.mat <- as.simple_triplet_matrix(SC.ID.mat)
#putting together the overall constraint matrix
cm <- rbind(cbind(FP.ID.mat,-simple_triplet_diag_matrix(NIDs_vs_FP)),cbind(FP.ID.mat,-simple_triplet_diag_matrix(NIDs_vs_FP)),cbind(SC.ID.mat,simple_triplet_zero_matrix(N.SC,N.FP)))
#making a non-sparse overall constraint matrix for lpSolve
cm2 <- as.matrix(cm)
#directions and rhs of the constraints
cdir <- c(rep("<=",N.FP),rep(">=",N.FP),rep("==",N.SC))
crhs <- unname(c(rep(0,N.FP),1-NIDs_vs_FP,N.ID.SC))
#objective vector, assigning -1 to each FP feature
cobj <- c(rep(0,N.ID),rep(-1,N.FP))
#solution using Rsymphony
sol.Rsymphony <- Rsymphony_solve_LP (cobj, cm, cdir, crhs, types = "B", max = FALSE, verbosity=0, first_feasible = FALSE)$solution
#solution using lpSolve
sol.lpSolve <- lp("min", cobj, cm2, cdir, crhs, all.bin=TRUE, num.bin.solns = 2, use.rw = FALSE, transpose.constraints = TRUE)$solution
#solution using Rglpk
sol.Rglpk <- Rglpk_solve_LP (cobj, cm, cdir, crhs, types = "B", max = FALSE, control=list(verbose=TRUE,presolve=FALSE))$solution
#comparison of the solutions
sol.Rsymphony
#[1] 1 0 0 1 1 0 1 1 0 1 1 1 1 1 1 1 1 1 1
sol.lpSolve
#[1] 1 0 0 1 1 0 1 1 0 1 1 1 1 1 1 1 1 1 1 0 1 1 0 1 0 1 1 0 1 1 1 1 1 1 1 1 1 1 1
sol.Rglpk
#[1] 0 1 1 0 1 0 1 1 0 1 1 1 1 1 1 1 1 1 1
И Rsymphony, и Rglpk находят одно решение, и в каждом случае оно отличается; lpSolve находит оба решения.
К сожалению, я не могу использовать lpSolve для моей реальной проблемы, потому что она слишком велика, поэтому не только фактическое решение займет слишком много времени, но я даже не смог создать исходную матрицу. Я мог бы попытаться поиграть с параметром «density.const» ...
Любые комментарии или предложения приветствуются.