R语言求解线性规划和非线性规划

第七章线性规划与非线性规划

例1m a x z=10x

1+5x2

s.t.5x1+2x2<=8

3x1+4x2=9

x1+x2>=1

x1,x2>=0

首先可化为标准形式：min - z = -10x1 -5x2

s.t. 5x1+2x1<=8

-x1-x2<=-1

3x1+4x2=9

x1,x2>=0

library(Rglpk)

obj<-c(-10,-5)

mat<-matrix(c(5,2,-1,-1,3,4),3,2,T)

dir<-c("<=","<=","==")

rhs<-c(8,-1,9)

Rglpk_solve_LP(obj,mat,dir,rhs)

#直接求解

library(Rglpk)

obj<-c(10,5)

mat<-matrix(c(5,2,1,1,3,4),3,2,T)

dir<-c("<=",">=","==")

rhs<-c(8,1,9)

Rglpk_solve_LP(obj,mat,dir,rhs,max=T)

非线性规划求解（Rdonlp2）

例2 有如下的条件约束最优化问题：

22min(sin cos )

1001001001002133

2sin cos 3z x y y x x y x y x y xy x y =+-<<⎧⎪-<<⎪⎪≤+⎨≤-≤⎪⎪=⎪≤⎩

library (Rdonlp2) p = c(10,10) #迭代初始值

#对求解问题进行描述

22min(sin cos )z x y y x =+

fn = function (x){

x[1]^2*sin(x[2])+x[2]^2*cos(x[1])

}

#对x,y 值域描述

100100100100

x y -<<-<< ## par.l 和par.u 分别为约束的左边和右边

par.l = c(-100,-100); par.u = c(100,100) ## 目标值域

#对线性约束进行描述

2133

x y x y ≤+≤-≤ A = matrix(c(1,1,3,-1),2,byrow=TRUE ) ##线性约束系数

lin.l = c(2,1); lin.u = c(+Inf ,3) ## 分别为约束的左边和右边

#对非线性约束进行描述

2sin cos 3

xy x y =≤ nlcon1 = function (x){

x[1]*x[2] ##公式 x*y

}

nlcon2 = function (x){

sin(x[1])*cos(x[2]) ##公式 sin(x)*cos(y)

}

## 两个非线性约束的左右边

## x*y=2 等价于 2<=x*y<=2

nlin.l = c(2,-Inf ) ; nlin.u = c(2,0.6) #将参数输入donlp2函数中进行求解

## 输入参数第一行： x,y 值域及目标函数

## 输入参数第二行： 线性约束条件

## 输入参数第三，四行： 非线性约束条件

ret = donlp2(p, fn, par.u=par.u, par.l=par.l,

A, lin.l=lin.l,lin.u=lin.u,

nlin=list(nlcon1,nlcon2),

nlin.u=nlin.u, nlin.l=nlin.l)

## 输出结果

ret$par

ret$par

例3 解下列二次规划

[][]111222111()26122x x f x x x x x -⎡⎤⎡⎤⎡⎤=-⎢⎥⎢⎥⎢⎥-⎣⎦⎣⎦⎣⎦

二次规划的优化问题，这是一种特殊形式的非线性约束优化问题。

二次规划在许多领域都有运用，比如投资组合优化、求解支持向量机(SVM)分类问题等。

想要用quadprog 包求解二次规划，我们需要同时转化我们的目标函数和约束条件为矩阵形式。

quadprog 包默认是求解最小化问题，目标函数二次，约束一次。所以，我们的约束条件默认的形式也就是AX>=bvec 。

通常我们需要把一些原来是求极大值的问题或者<=约束通过乘以负号来转化 library(quadprog)

Dmat <-matrix(c(1,-1,1,2),2,2,T)

Dmat

dvec <- c(2,6)

A<- matrix(-c(1,1,-1,2,2,1),3,2,T)

A

bvec <- c(-2,-2,-3)

Amat <- t(A)

sol <- solve.QP(Dmat, dvec, Amat, bvec)

sol

参数Dmat 表示海赛矩阵

参数dvet 表示一阶向量，和Dmat 的维数要相对应。

参数Amat 表示约束矩阵，默认的约束都是>=。

参数bvet 表示右边值，由向量，和Amat 的维数要相对应。

参数 meq 表示从哪一行开始Amat 矩阵中的约束是需要被当作等式约束的。

例4 假设以决策变量x1、x2、x3分别表示甲、乙、丙、丁4种肥料的用量，得线性规划模型

1234124134141234min 0.040.150.10.125..0.030.30.15320.050.20.1240.140.0742,,,0z x x x x s t x x x x x x x x x x x x =+++⎧⎪++≥⎪⎪++=⎨⎪+≤⎪⎪≥⎩

library(Rglpk)

obj<-c(0.04,0.15,0.1,0.125)

mat<-matrix(c(0.03,0.3,0,0.15,0.05,0,0.2,0.1,0.14,0,0,0.07),3,4,T)

mat

dir<-c(">=","==","<=")

rhs<-c(32,24,42)

Rglpk_solve_LP(obj,mat,dir,rhs)