R Notebook
1.参考资料
1.1 Modern optimization with R
Cortez, P. (2021). Modern optimization with R.Springer.
cover_modern_optimization_with_R.jpg
1.2 Modeling and Solving Linear Programming with R
1.3 Linear Programming with R: Exploring the “lpSolve” R package
Roberto Salazar,Nov 17, 2019
Linear Programming with R: Exploring the “lpSolve” R package
1.4 lpSolveAPI Package Users Guide by Kjell Konis
2. prepraration
2.1 import libraries
# Import lpSolve package
library(lpSolve)
library(XLConnect)
2.2 the problem: 班次安排问题
## 2.3 连接Excel文件
mybook<-loadWorkbook("D:/kedu/teaching_datasets/Excel_model/sm_solver.xls")
rsheets<-getSheets(mybook)
rsheets
[1] "table-chair" "beef" "transport" "staff"
3. linear programming
3.1 Set coefficients of the objective function
设定目标函数,因为是求和,所以矩阵是1,1,1,1,1,1 matrix在这里非常重要 研究一下如何生成重复的相同元素
#f.obj <- c(4, 2)
f.obj <-as.matrix(rep(1,each=7))
3.2 Set matrix corresponding to coefficients of constraints by rows
Do not consider the non-negative constraint; it is automatically assumed
就是员工工作的可用矩阵 注意:这里面通常包含标题,默认header=TRUE,那么startRow=7
#f.con <- matrix(c(5, 15,20, 5), nrow = 2, byrow = TRUE)
f.con <-as.matrix(readWorksheet(mybook, sheet = "staff", startRow = 8, endRow = 14,
startCol = 17, endCol = 23,header=FALSE))
3.3 Set unequality signs
设置不等式符号 可以每行设置不同的符号
f.dir <- as.matrix(rep(">=",each=7))
3.4 Set right hand side coefficients
约束条件的范围
#f.rhs <- c(50,40)
f.rhs <-as.matrix(readWorksheet(mybook, sheet = "staff", startRow = 8, endRow =14,
startCol = 27, endCol = 27,header=FALSE))
3.5 设定变量取整数 f.intvec,例如f.intvec <- c(1,2)表示x1,x2取整数
如果自变量全都是整数,那么 all.int = TRUE就可以了
#f.intvec <- c(1,2)
3.6 Final value (z)
计算结果
#report_lp<-lp("min", f.obj, f.con, f.dir, f.rhs,int.vec = f.intvec)
# 部分变量取整数:int.vec = f.intvec
report_lp<-lp("min", f.obj, f.con, f.dir, f.rhs, all.int = TRUE)
report_lp
Success: the objective function is 9
report_lp$objval
[1] 9
# output the final value
writeWorksheet(mybook,report_lp$objval,sheet = "staff", startRow =16,startCol = 25,header = FALSE)
saveWorkbook(mybook)
3.7 Variables final values
report_solution<-report_lp$solution
report_solution
[1] 0 1 0 0 1 3 4
## 保存变量取值
## 结果是一个矩阵(1列),所以为了在excel变为1行,需要转置,t()
writeWorksheet(mybook,t(report_solution),sheet = "staff", startRow =16,startCol = 17,header = FALSE)
saveWorkbook(mybook)
3.8 Sensitivities
敏感度分析
report_lp<-lp("min", f.obj, f.con, f.dir, f.rhs, all.int = TRUE,compute.sens=TRUE)
report_lp$sens.coef.from
[1] 1e+00 -1e+30 0e+00 0e+00 0e+00 0e+00 0e+00
report_lp$sens.coef.to
[1] 1.000000e+30 1.000000e+30 1.333333e+00 1.000000e+00 1.333333e+00
[6] 1.333333e+00 1.333333e+00
3.9 Dual Values (first dual of the constraints and then dual of the variables)
Duals of the constraints and variables are mixed
report_lp$duals
[1] 0 1 0 0 0 0 0 0 1 0 0 0 0 0
3.10 Duals lower and upper limits
report_lp$duals.from
[1] -1e+30 -1e+30 7e+00 6e+00 6e+00 5e+00 4e+00 -1e+30 -1e+30 -1e+30
[11] -1e+30 -1e+30 -1e+30 -1e+30
report_lp$duals.to
[1] 1.000000e+30 1.000000e+30 8.333333e+00 8.666667e+00 8.666667e+00
[6] 9.000000e+00 9.333333e+00 1.250000e+00 2.500000e-01 1.000000e+30
[11] 1.000000e+30 1.000000e+30 1.000000e+30 1.000000e+30
LS0tDQp0aXRsZTogIlIgTm90ZWJvb2siDQpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sNCi0tLQ0KDQojIDEu5Y+C6ICD6LWE5paZDQojIyMgMS4xIE1vZGVybiBvcHRpbWl6YXRpb24gd2l0aCBSDQpDb3J0ZXosIFAuICgyMDIxKS4gW01vZGVybiBvcHRpbWl6YXRpb24gd2l0aCBSLl0oaHR0cHM6Ly9saW5rLnNwcmluZ2VyLmNvbS9ib29rLzEwLjEwMDcvOTc4LTMtMzE5LTA4MjYzLTkpU3ByaW5nZXIuDQoNCltzcHJpbmdlciBsaW5rXShodHRwczovL2xpbmsuc3ByaW5nZXIuY29tL2Jvb2svMTAuMTAwNy85NzgtMy0zMTktMDgyNjMtOSkNCg0KW2FtYXpvbiBsaW5rXShodHRwczovL3d3dy5hbWF6b24uY29tL2RwLzMwMzA3MjgxODgpDQoNCltsaWJyYXJ5IGdlbmVzaXNdKGh0dHA6Ly9saWJnZW4ucnMvYm9vay9pbmRleC5waHA/bWQ1PTBGNDlCMkFCMjMxODE5M0NFOTE4QjM2RTBCMUUxNjQyKQ0KDQohW10oRDovd2VpeXVuX3N5bmMvcGhvdG9fYW5kX2ltYWdlL2thaXd1X2ltYWdlcy9wZGF0YV9pbWFnZXMvY292ZXJfbW9kZXJuX29wdGltaXphdGlvbl93aXRoX1IuanBnKQ0KY292ZXJfbW9kZXJuX29wdGltaXphdGlvbl93aXRoX1IuanBnDQoNCg0KIyMjIDEuMiBNb2RlbGluZyBhbmQgU29sdmluZyBMaW5lYXIgUHJvZ3JhbW1pbmcgd2l0aCBSDQpbaHR0cHM6Ly91cGNvbW1vbnMudXBjLmVkdS9iaXRzdHJlYW0vaGFuZGxlLzIxMTcvNzgzMzUvTW9kZWxpbmcrYW5kK1NvbHZpbmcrTGluZWFyK1Byb2dyYW1taW5nK3dpdGgrUi5wZGZdKGh0dHBzOi8vdXBjb21tb25zLnVwYy5lZHUvYml0c3RyZWFtL2hhbmRsZS8yMTE3Lzc4MzM1L01vZGVsaW5nK2FuZCtTb2x2aW5nK0xpbmVhcitQcm9ncmFtbWluZyt3aXRoK1IucGRmKQ0KDQohW10oRDovd2VpeXVuX3N5bmMvcGhvdG9fYW5kX2ltYWdlL2thaXd1X2ltYWdlcy9wZGF0YV9pbWFnZXMvY292ZXJfTW9kZWxpbmdfYW5kX1NvbHZpbmdfTGluZWFyX1Byb2dyYW1taW5nX3dpdGhfUi5qcGcpDQoNCg0KIyMjIDEuMyBMaW5lYXIgUHJvZ3JhbW1pbmcgd2l0aCBSOiBFeHBsb3JpbmcgdGhlIOKAnGxwU29sdmXigJ0gUiBwYWNrYWdlDQpSb2JlcnRvIFNhbGF6YXIsTm92IDE3LCAyMDE5DQoNCltMaW5lYXIgUHJvZ3JhbW1pbmcgd2l0aCBSOiBFeHBsb3JpbmcgdGhlIOKAnGxwU29sdmXigJ0gUiBwYWNrYWdlXShodHRwczovL3Rvd2FyZHNkYXRhc2NpZW5jZS5jb20vbGluZWFyLXByb2dyYW1taW5nLWluLXItNDQ0ZTljMTk5MjgwKQ0KDQoNCiMjIyAxLjQgbHBTb2x2ZUFQSSBQYWNrYWdlIFVzZXJzIEd1aWRlIGJ5IEtqZWxsIEtvbmlzDQpbbHBTb2x2ZUFQSSBQYWNrYWdlIFVzZXJzIEd1aWRlXShodHRwOi8vY2l2aWwuY29sb3JhZG8uZWR1L35iYWxhamlyL0NWRU41MzkzL1Itc2Vzc2lvbnMvc2VzczEvbHBTb2x2ZUFQSS12aWduZXR0ZXMucGRmKQ0KDQoNCg0KIyAyLiBwcmVwcmFyYXRpb24NCg0KIyMgMi4xIGltcG9ydCBsaWJyYXJpZXMNCg0KYGBge3J9DQojIDIuMSBJbXBvcnQgbHBTb2x2ZSBwYWNrYWdlDQpsaWJyYXJ5KGxwU29sdmUpDQpsaWJyYXJ5KFhMQ29ubmVjdCkNCmBgYA0KDQoNCg0KDQoNCiMjIDIuMiB0aGUgcHJvYmxlbTog54+t5qyh5a6J5o6S6Zeu6aKYDQoNCiFbaW1hZ2VuYW1lXShEOi93ZWl5dW5fc3luYy9waG90b19hbmRfaW1hZ2Uva2Fpd3VfaW1hZ2VzL3BkYXRhX2ltYWdlcy9yX3N0YWZmLmpwZykNCiMjICAyLjMg6L+e5o6lRXhjZWzmlofku7YNCmBgYHtyfQ0KbXlib29rPC1sb2FkV29ya2Jvb2soIkQ6L2tlZHUvdGVhY2hpbmdfZGF0YXNldHMvRXhjZWxfbW9kZWwvc21fc29sdmVyLnhscyIpDQpyc2hlZXRzPC1nZXRTaGVldHMobXlib29rKQ0KcnNoZWV0cw0KYGBgDQoNCiMgMy4gbGluZWFyIHByb2dyYW1taW5nDQoNCg0KIyMjIDMuMSBTZXQgY29lZmZpY2llbnRzIG9mIHRoZSBvYmplY3RpdmUgZnVuY3Rpb24NCuiuvuWumuebruagh+WHveaVsO+8jOWboOS4uuaYr+axguWSjO+8jOaJgOS7peefqemYteaYrzHvvIwx77yMMe+8jDHvvIwx77yMMQ0KbWF0cml45Zyo6L+Z6YeM6Z2e5bi46YeN6KaBDQrnoJTnqbbkuIDkuIvlpoLkvZXnlJ/miJDph43lpI3nmoTnm7jlkIzlhYPntKANCg0KDQpgYGB7cn0NCiNmLm9iaiA8LSBjKDQsIDIpDQpmLm9iaiA8LWFzLm1hdHJpeChyZXAoMSxlYWNoPTcpKQ0KYGBgDQoNCiMjIyAzLjIgU2V0IG1hdHJpeCBjb3JyZXNwb25kaW5nIHRvIGNvZWZmaWNpZW50cyBvZiBjb25zdHJhaW50cyBieSByb3dzDQoNCiMjIyBEbyBub3QgY29uc2lkZXIgdGhlIG5vbi1uZWdhdGl2ZSBjb25zdHJhaW50OyBpdCBpcyBhdXRvbWF0aWNhbGx5IGFzc3VtZWQNCuWwseaYr+WRmOW3peW3peS9nOeahOWPr+eUqOefqemYtQ0K5rOo5oSP77ya6L+Z6YeM6Z2i6YCa5bi45YyF5ZCr5qCH6aKY77yM6buY6K6kaGVhZGVyPVRSVUXvvIzpgqPkuYhzdGFydFJvdz03DQpgYGB7cn0NCiNmLmNvbiA8LSBtYXRyaXgoYyg1LCAxNSwyMCwgNSksIG5yb3cgPSAyLCBieXJvdyA9IFRSVUUpDQoNCmYuY29uIDwtYXMubWF0cml4KHJlYWRXb3Jrc2hlZXQobXlib29rLCBzaGVldCA9ICJzdGFmZiIsIHN0YXJ0Um93ID0gOCwgZW5kUm93ID0gMTQsDQpzdGFydENvbCA9IDE3LCBlbmRDb2wgPSAyMyxoZWFkZXI9RkFMU0UpKQ0KDQpgYGANCg0KDQojIyMgMy4zIFNldCB1bmVxdWFsaXR5IHNpZ25zDQrorr7nva7kuI3nrYnlvI/nrKblj7cNCuWPr+S7peavj+ihjOiuvue9ruS4jeWQjOeahOespuWPtw0KYGBge3J9DQpmLmRpciA8LSBhcy5tYXRyaXgocmVwKCI+PSIsZWFjaD03KSkNCg0KDQpgYGANCg0KDQoNCiMjIyAzLjQgU2V0IHJpZ2h0IGhhbmQgc2lkZSBjb2VmZmljaWVudHMNCue6puadn+adoeS7tueahOiMg+WbtA0KYGBge3J9DQojZi5yaHMgPC0gYyg1MCw0MCkNCmYucmhzIDwtYXMubWF0cml4KHJlYWRXb3Jrc2hlZXQobXlib29rLCBzaGVldCA9ICJzdGFmZiIsIHN0YXJ0Um93ID0gOCwgZW5kUm93ID0xNCwNCnN0YXJ0Q29sID0gMjcsIGVuZENvbCA9IDI3LGhlYWRlcj1GQUxTRSkpDQpgYGANCg0KIyMjIDMuNSDorr7lrprlj5jph4/lj5bmlbTmlbAgZi5pbnR2ZWPvvIzkvovlpoJmLmludHZlYyA8LSBjKDEsMinooajnpLp4Me+8jHgy5Y+W5pW05pWwDQrlpoLmnpzoh6rlj5jph4/lhajpg73mmK/mlbTmlbDvvIzpgqPkuYggYWxsLmludCA9IFRSVUXlsLHlj6/ku6XkuoYNCmBgYHtyfQ0KI2YuaW50dmVjIDwtIGMoMSwyKSAgIA0KYGBgDQoNCg0KIyMjIDMuNiBGaW5hbCB2YWx1ZSAoeikNCuiuoeeul+e7k+aenA0KDQpgYGB7cn0NCiNyZXBvcnRfbHA8LWxwKCJtaW4iLCBmLm9iaiwgZi5jb24sIGYuZGlyLCBmLnJocyxpbnQudmVjID0gZi5pbnR2ZWMpDQojIOmDqOWIhuWPmOmHj+WPluaVtOaVsO+8mmludC52ZWMgPSBmLmludHZlYw0KDQpyZXBvcnRfbHA8LWxwKCJtaW4iLCBmLm9iaiwgZi5jb24sIGYuZGlyLCBmLnJocywgYWxsLmludCA9IFRSVUUpDQpyZXBvcnRfbHANCnJlcG9ydF9scCRvYmp2YWwNCg0KIyBvdXRwdXQgdGhlIGZpbmFsIHZhbHVlDQp3cml0ZVdvcmtzaGVldChteWJvb2sscmVwb3J0X2xwJG9ianZhbCxzaGVldCA9ICJzdGFmZiIsIHN0YXJ0Um93ID0xNixzdGFydENvbCA9IDI1LGhlYWRlciA9IEZBTFNFKQ0Kc2F2ZVdvcmtib29rKG15Ym9vaykNCmBgYA0KDQoNCiMjIyAzLjcgVmFyaWFibGVzIGZpbmFsIHZhbHVlcw0KDQpgYGB7cn0NCnJlcG9ydF9zb2x1dGlvbjwtcmVwb3J0X2xwJHNvbHV0aW9uDQpyZXBvcnRfc29sdXRpb24NCg0KIyMg5L+d5a2Y5Y+Y6YeP5Y+W5YC8DQojIyDnu5PmnpzmmK/kuIDkuKrnn6npmLXvvIgx5YiX77yJ77yM5omA5Lul5Li65LqG5ZyoZXhjZWzlj5jkuLox6KGM77yM6ZyA6KaB6L2s572u77yMdCgpDQp3cml0ZVdvcmtzaGVldChteWJvb2ssdChyZXBvcnRfc29sdXRpb24pLHNoZWV0ID0gInN0YWZmIiwgc3RhcnRSb3cgPTE2LHN0YXJ0Q29sID0gMTcsaGVhZGVyID0gRkFMU0UpDQpzYXZlV29ya2Jvb2sobXlib29rKQ0KYGBgDQoNCg0KIyMjIDMuOCBTZW5zaXRpdml0aWVzDQrmlY/mhJ/luqbliIbmnpANCg0KYGBge3J9DQpyZXBvcnRfbHA8LWxwKCJtaW4iLCBmLm9iaiwgZi5jb24sIGYuZGlyLCBmLnJocywgYWxsLmludCA9IFRSVUUsY29tcHV0ZS5zZW5zPVRSVUUpDQpgYGANCg0KDQoNCmBgYHtyfQ0KcmVwb3J0X2xwJHNlbnMuY29lZi5mcm9tDQpyZXBvcnRfbHAkc2Vucy5jb2VmLnRvDQpgYGANCg0KDQoNCiMjIyAzLjkgRHVhbCBWYWx1ZXMgKGZpcnN0IGR1YWwgb2YgdGhlIGNvbnN0cmFpbnRzIGFuZCB0aGVuIGR1YWwgb2YgdGhlIHZhcmlhYmxlcykNCiMjIyBEdWFscyBvZiB0aGUgY29uc3RyYWludHMgYW5kIHZhcmlhYmxlcyBhcmUgbWl4ZWQNCmBgYHtyfQ0KcmVwb3J0X2xwJGR1YWxzDQpgYGANCg0KDQojIyMgMy4xMCBEdWFscyBsb3dlciBhbmQgdXBwZXIgbGltaXRzDQpgYGB7cn0NCnJlcG9ydF9scCRkdWFscy5mcm9tDQpyZXBvcnRfbHAkZHVhbHMudG8NCmBgYA0KDQoNCg0KDQoNCg==