生产计划多阶段动态规划LINGO代码

网上的文章附带的代码有误，我花时间修改了下。

model:
max=11328*(1-k2)-0.55*x1^2-0.55*x2^2-0.6*x3^2-0.6*x4^2-0.65*x5^2-0.65*x6^2
!Storage cost in january;
-0.1*28
!Storage cost in february;
-0.1*@if(x1#gt#25,x1-9,(x1+7)^2/64)  
!Storage cost in march;
-0.1*@if(x1#gt#25#and#x1+x2#gt#66,x1+x2-45.5,
@if(x1#gt#25,(x1+x2-25)^2/82,
@if(x1#le#25#and#2*x1+x2#gt#91,2*x1+x2-70.5,(2*x1+x2-50)^2/82)))
!Storage cost in april;
-0.1*@if(x1#le#25#and#2*x1+x2#gt#91#and# x1+x2+x3#gt#133, x1+x2+x3-99.5,
@if(x1#le#25#and#2*x1+x2#gt#91, (x1+x2+x3-66)^2/134,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224, (3*x1+2*x2+x3-190.5), 
@if(x1#le#25,(3*x1+2*x2+x3-157)^2/134,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133, x1+x2+x3-99.5,
@if(x1#gt#25#and#x1+x2#gt#66, (x1+x2+x3-66)^2/134,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199,2*x1+2*x2+x3-165.5,(2*x1+2*x2+x3-132)^2/134)))))))
!Storage cost in may;
-0.1*@if(x1#le#25#and#2*x1+x2#gt#91#and#x1+x2+x3#gt#133#and#x1+x2+x3+x4#gt#158,x1+x2+x3+x4-145.5,
@if(x1#le#25#and#2*x1+x2#gt#91#and#x1+x2+x3#gt#133, (x1+x2+x3+x4-133)^2/50,
@if(x1#le#25#and#2*x1+x2#gt#91#and#2*x1+2*x2+2*x3+x4#gt#291,2*x1+2*x2+2*x3+x4-278.5,
@if(x1#le#25#and#2*x1+x2#gt#91, (2*x1+2*x2+2*x3+x4-266)^2/50,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224#and#x1+x2+x3+x4#gt#158, x1+x2+x3+x4-145.5,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224, (x1+x2+x3+x4-133)^2/50,
@if(x1#le#25#and#4*x1+3*x2+2*x3+x4#gt#382, 4*x1+3*x2+2*x3+x4-369.5,
@if(x1#le#25, (4*x1+3*x2+2*x3+x4-357)^2/50,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133#and#x1+x2+x3+x4#gt#158,x1+x2+x3+x4-145.5,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133, (x1+x2+x3+x4-133)^2/50,
@if(x1#gt#25#and#x1+x2#gt#66#and#2*x1+2*x2+2*x3+x4#gt#291,2*x1+2*x2+2*x3+x4-278.5, 
@if(x1#gt#25#and#x1+x2#gt#66,(2*x1+2*x2+2*x3+x4-266)^2/50,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199#and#x1+x2+x3+x4#gt#158,x1+x2+x3+x4-145.5,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199, (x1+x2+x3+x4-133)^2/50,
@if(x1#gt#25#and#3*x1+3*x2+2*x3+x4#gt#357,3*x1+3*x2+2*x3+x4-344.5,(3*x1+3*x2+2*x3+x4-332)^2/50)))))))))))))))
!Storage cost in june;
-0.1*@if(x1#le#25#and#2*x1+x2#gt#91#and#x1+x2+x3#gt#133#and#x1+x2+x3+x4#gt#158#and#x1+x2+x3+x4+x5#gt#187, x1+x2+x3+x4+x5-172.5,
@if(x1#le#25#and#2*x1+x2#gt#91#and#x1+x2+x3#gt#133#and#x1+x2+x3+x4#gt#158,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#le#25#and#2*x1+x2#gt#91#and#x1+x2+x3#gt#133#and#2*x1+2*x2+2*x3+2*x4+x5#gt#345, 2*x1+2*x2+2*x3+2*x4+x5-330.5,
@if(x1#le#25#and#2*x1+x2#gt#91#and#x1+x2+x3#gt#133,(2*x1+2*x2+2*x3+2*x4+x5-316)^2/58,
@if(x1#le#25#and#2*x1+x2#gt#91#and#2*x1+2*x2+2*x3+x4#gt#291#and#x1+x2+x3+x4+x5#gt#187, x1+x2+x3+x4+x5-172.5,
@if(x1#le#25#and#2*x1+x2#gt#91#and#2*x1+2*x2+2*x3+x4#gt#291,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#le#25#and#2*x1+x2#gt#91#and#3*x1+3*x2+3*x3+2*x4+x5#gt#478,3*x1+3*x2+3*x3+2*x4+x5-463.5,
@if(x1#le#25#and#2*x1+x2#gt#91, (3*x1+3*x2+3*x3+2*x4+x5-449)^2/58,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224#and#x1+x2+x3+x4#gt#158#and#x1+x2+x3+x4+x5#gt#187, x1+x2+x3+x4+x5-172.5,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224#and#x1+x2+x3+x4#gt#158,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224#and#2*x1+2*x2+2*x3+2*x4+x5#gt#345,2*x1+2*x2+2*x3+2*x4+x5-330.5,
@if(x1#le#25#and#3*x1+2*x2+x3#gt#224, (2*x1+2*x2+2*x3+2*x4+x5-316)^2/58,
@if(x1#le#25#and#4*x1+3*x2+2*x3+x4#gt#382#and#x1+x2+x3+x4+x5#gt#187,x1+x2+x3+x4+x5-172.5,
@if(x1#le#25#and#4*x1+3*x2+2*x3+x4#gt#382, (x1+x2+x3+x4+x5-158)^2/58,
@if(x1#le#25#and#5*x1+4*x2+3*x3+2*x4+x5#gt#569, 5*x1+4*x2+3*x3+2*x4+x5-554.5,
@if(x1#le#25, (5*x1+4*x2+3*x3+2*x4+x5-540)^2/58,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133#and#x1+x2+x3+x4#gt#158#and#x1+x2+x3+x4+x5#gt#187, x1+x2+x3+x4+x5-172.5,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133#and#x1+x2+x3+x4#gt#158,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133#and#2*x1+2*x2+2*x3+2*x4+x5#gt#345,2*x1+2*x2+2*x3+2*x4+x5-330.5,
@if(x1#gt#25#and#x1+x2#gt#66#and#x1+x2+x3#gt#133,(2*x1+2*x2+2*x3+2*x4+x5-316)^2/58,
@if(x1#gt#25#and#x1+x2#gt#66#and#2*x1+2*x2+2*x3+x4#gt#291#and#x1+x2+x3+x4+x5#gt#187,x1+x2+x3+x4+x5-172.5,
@if(x1#gt#25#and#x1+x2#gt#66#and#2*x1+2*x2+2*x3+x4#gt#291,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#gt#25#and#x1+x2#gt#66#and#3*x1+3*x2+3*x3+2*x4+x5#gt#478,3*x1+3*x2+3*x3+2*x4+x5-463.5,
@if(x1#gt#25#and#x1+x2#gt#66, (3*x1+3*x2+3*x3+2*x4+x5-449)^2/58,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199#and#x1+x2+x3+x4#gt#158#and#x1+x2+x3+x4+x5#gt#187, x1+x2+x3+x4+x5-172.5,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199#and#x1+x2+x3+x4#gt#158,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199#and#2*x1+2*x2+2*x3+2*x4+x5#gt#345,2*x1+2*x2+2*x3+2*x4+x5-330.5,
@if(x1#gt#25#and#2*x1+2*x2+x3#gt#199,(2*x1+2*x2+2*x3+2*x4+x5-316)^2/58,
@if(x1#gt#25#and#3*x1+3*x2+2*x3+x4#gt#357#and#x1+x2+x3+x4+x5#gt#187,x1+x2+x3+x4+x5-172.5,
@if(x1#gt#25#and#3*x1+3*x2+2*x3+x4#gt#357,(x1+x2+x3+x4+x5-158)^2/58,
@if(x1#gt#25#and#4*x1+4*x2+3*x3+2*x4+x5#gt#544,4*x1+4*x2+3*x3+2*x4+x5-529.5, (4*x1+4*x2+3*x3+2*x4+x5-515)^2/58)))))))))))))))))))))))))))))));
x1+x2+x3+x4+x5+x6=187;
x1>0;  
x2>0;  
x3>0;
x4>0;
x5>0;
x6>0;
@gin(x1);
@gin(x2);
@gin(x3);
@gin(x4);
@gin(x5);
@gin(x6);
k2=0.1;
end

网上文章的结果也可能有问题，我的代码要勾选住use global solver选项进行运行，运行结果和一篇用matlab做的结果一样，为33 34 31 31 29 29