Lingo与线性规划

Lingo与线性规划线性规划的标准形式是11nnMinzcxcx1111111..0,1,2,,nnmmnnmiaxaxbstaxaxbxin(1)其中11nnzcxcx称为目标函数，自变量ix称为决策变量,不等式组(1)称为约束条件.满足不等式组(1)的所有1(,,)nxx的集合称为可行域，在可行域里面使得z取最小值的**1(,,)nxx称为最优解，最优解对应的函数值称为最优值。求解优化模型的主要软件有Lingo、Matlab、Excel等。其中Lingo是一款专业求解优化模型的软件，有其他软件不可替代的方便功能。本文将简要介绍其在线性规划领域的应用。一、基本规定1、目标函数输入格式max=函数解析式；或者min=函数解析式；2、约束条件输入格式利用：、、=、=等符号。但是与=没有区别。Lingo软件默认所以自变量都大于等于0.3、运算加(+),减(-),乘(*),除(/),乘方(x^a)，要注意乘号(*)不能省略。4、变量名不区分大小写字母，不超过32个字符，必须以字母开头。5、标点符号每个语句以分号“；”结束，感叹号“！”开始的是说明语句（说明语句也需要以分号“;”结束）。但是，model，sets，data以“：”结尾。endsets，enddata，end尾部不加任何符号。6、命令不考虑先后次序7、MODEL语句一般程序必须先输入MODEL：表示开始输入模型，以“END”结束。对简单的模型，这两个语句也可以省略。8、改变变量的取值范围@bin(变量名)；限制该变量为0或1.@bnd(a,变量名,b)；限制该变量介于a,b之间.@free(变量名)；允许该变量为负数.@gin(变量名)；限制该变量为整数.例1求目标函数1223zxx的最小值，约束条件为1211212..35010022600,0stxxxxxxx输入Lingo程序：min=2*x1+3*x2;x1+x2=350;x1=100;2*x1+x2=600;有两种运行方式：1、点击工具条上的按钮即可。2、点击菜单：LINGO→Solve运行结果如下：下面对其各个部分进行说明：Globaloptimalsolutionfound：表示已找到全局最优解。Objectivevalue：表示最优值的大小。可见本题函数最小值minz800。Infeasibilities：矛盾约束的数目。Totalsolveriterations:迭代次数。Variable：变量。本题有两个变量。Value：变量对应的最优解，即12250,100xx。ReducedCost：变量ix在最优解的基础上增加一个单位，目标函数值的改变量。例如，一个变量的ReducedCost值为8，那么当该变量增加一个单位，在最大化（最小化）问题中目标函数值将减少（增大）8个单位。SlackorSurplus：表示接近等于的程度，即约束离相等还差多少。在约束条件是=中，表示松弛程度，在约束条件是=中，不是过剩程度。如果约束条件是=，则SlackorSurplus为0，该约束是个紧约束(或有效约束)。如果一个约束是矛盾的，即模型无可行解则Slackorsurplus的值是负数。知道SlackorSurplus的值，可以帮助我们发现优化模型中错误的约束条件。在上例中第2和第4行松弛变量均为0,说明对于最优解来讲,两个约束(第2和4行)均取等号,即都是紧约束，第3行为150，即最优解使得第3行过剩150.DualPrice：对偶价格的值，它表示约束条件中的常数，每增加一个单位，目标函数值改变的数量（在最大化问题中目标函数值是增加，在最小化问题中目标函数值是减少）。比如，在上一个Min模型中第四行的1，表示2*x1+x2=600增加一个单位到2*x1+x2=601，可以使目标值增加-1(因为第一行是目标函数的DualPrice是-1)，即Objectivevalue=799;如果增加-1个单位到599会使目标值增加到801。例2求目标函数2221234212zxxx的最小值，约束条件为123123123..3291,,stxxxxxxxxxR输入Lingo程序：min=4*x1^2-x2^2+2*x3^2+12;3*x1+2*x2+x3=9;x1+x2+x3=-1;@free(x1);@free(x2);@free(x3);运行结果：即当1231,8,10xxx时，min152z。二、灵敏度分析灵敏度分析是指：找出模型变量系数的一个变化范围，使得最优基(即最优解)保持不变。一般只对线性规划模型做灵敏度分析。1、灵敏度分析操作步骤第一步：菜单lingo--options--generalsolver--dualcomputations:prices&ranges--ok.第二步：菜单lingo--range2、灵敏度报告中常见的词汇Currentcoefficient：当前目标函数系数Allowableincrease：允许增加量Allowabledecrease：允许减少量CurrentRHS：当前右边常数项INFINITY：表示正无穷。例1求解下列模型：max=72*x1+64*x2;x1+x2=50;12*x1+8*x2=480;3*x1=100;并做灵敏度分析。求解报告：灵敏度分析报告：灵敏度分析报告的解读：x1的系数变化范围是（72-8，72+24）=（64，96）；x2的系数变化范围是（64-16，64+8）=（48，72）。注意：x1系数的允许范围需要x2系数64不变，反之亦然。由于目标函数的费用系数变化并不影响约束条件，因此此时最优基不变可以保证最优解也不变，但最优值变化。右边常数项中，第2行原来为50，当它在[50-6.67，50+10]=[43.33，60]范围变化时，最优基保持不变。第3行可以类似解释。对第4行，原来为100，当它在[100-40，100+∞]=[60，+∞]范围变化时，最优基保持不变。不过由于此时约束发生变化，最优基即使不变，最优解、最优值也会发生变化。三、数据输入对于大型的优化问题，即自变量比较多的时候，还像上两节那样输入目标函数和约束条件就比较麻烦了。一般输入数据的方法有两种：一、建立向量、矩阵输入；二、调用外部数据。这里仅介绍第一种方法。1、建立向量命令格式：集合名称/集合维数/：向量名称例如：sets：set1/1..9/:x;set2/1..5/:a,b;endsets表示建立了两类集合。第一类集合set1，维数为9，x和y是向量名。向量x=（x(1),…,x(9)）,其中x(i)是x的元素。第二类集合set2，维数为5，a和b都是向量名。向量a=（a(1),…,a(5)）,其中a(i)是a的元素。向量b=（b(1),…,b(5)）,其中b(i)是b的元素。2、建立矩阵命令格式：集合名称（集合1，集合2）/：矩阵名称例如：sets：set1/1..3/:x;set2/1..4/:a;link(set1,set2):A;endsets表示建立了一个矩阵类link，其矩阵的阶数为34，A是具体的矩阵名。有两个命令是比较常见的：求和语句:@sum(集合名(i):含集合名(i)的语句);循环语句:@for(集合名(i):循环的语句);例3：求目标函数121115zxx的最小值，约束条件为12121212..2030360302520003035300,0stxxxxxxxx输入Lingo程序：model:sets:set1/1..2/:c,x;set2/1..3/:b;link(set2,set1):A;endsetsmax=@sum(set1(i):c(i)*x(i));@for(set2(i):@sum(link(i,j):A(i,j)*x(j))=b(i));data:c=1115;A=203030253025;b=3602000300;enddataend运行结果报告：例4、某地区有三个蔬菜生产基地，估计每年可供应本地区的蔬菜量表为：生产基地ABC蔬菜生产量(吨)783有四个地市需要该类蔬菜，需求表为：地区甲乙丙丁蔬菜生产量(吨)6633如果从各蔬菜生产基地到各地市的每吨蔬菜的运价表(单位：万元/吨)为:地市生产基地甲乙丙丁A5879B49107C8429为了降低运输费，需要合理调拨资源.(1)根据以上资料表制订一个使总的运费为最少的蔬菜调拨方案．(2)如果有机会增加生产基地的产量1吨，问应当优先增加那个基地的产量？(3)如果将A到乙市的运价减少为5万元/吨，问这会影响最优的调拨方案吗？设iA：第i个蔬菜生产基地，1,2,3i，分别对应生产基地A，B，C；iB：第i个蔬菜需求地，1,2,3,4i，分别对应蔬菜需求地市甲、乙、丙、丁；Q：总运输费用；ijx：表示的是从第i个生产基地向第j个地市运输的蔬菜数量；ijc：表示的是从第i个生产基地向第j个地市运输蔬菜的运价；ib：第i个蔬菜生产基地的蔬菜产量；iq：第j个地市的蔬菜需求量；那么有优化模型：3411minijijijQcx..st11121314121222324231323334311213111222322132333314243440,1,2,3;1,2,3,4ijxxxxbxxxxbxxxxbxxxqxxxqxxxqxxxqxij输入Lingo程序求解模型：model:sets:set1/1..3/:b;set2/1..4/:q;link(set1,set2):c,x;endsetsmin=@sum(link(i,j):c(i,j)*x(i,j));@for(set1(i):@sum(link(i,j):x(i,j))=b(i));@for(set2(j):@sum(link(i,j):x(i,j))=q(j));data:c=5,8,7,94,9,10,78,4,2,9;b=7,8,3;q=6,6,3,3;enddataend运行结果如下：Globaloptimalsolutionfound.Objectivevalue:100.0000Infeasibilities:0.000000Totalsolveriterations:6VariableValueReducedCostB(1)7.0000000.000000B(2)8.0000000.000000B(3)3.0000000.000000Q(1)6.0000000.000000Q(2)6.0000000.000000Q(3)3.0000000.000000Q(4)3.0000000.000000C(1,1)5.0000000.000000C(1,2)8.0000000.000000C(1,3)7.0000000.000000C(1,4)9.0000000.000000C(2,1)4.0000000.000000C(2,2)9.0000000.000000C(2,3)10.000000.000000C(2,4)7.0000000.000000C(3,1)8.0000000.000000C(3,2)4.0000000.000000C(3,3)2.0000000.000000C(3,4)9.0000000.000000X(1,1)1.0000000.000000X(1,2)6.0000000.000000X(1,3)0.0000001.000000X(1,4)0.0000001.000000X(2,1)5.0000000.000000X(2,2)0.0000002.000000X(2,3)0.0000005.000000X(2,4)3.0000000.000000X(3,1)0.0000007.000000X(3,2)0.0000000.000000