蒙特卡罗方法求解有约束的非线性规划问题的matlab程序.doc

蒙特卡罗方法求解有约束的非线性规划问题的matlab程序首先我们说明一下：观察函数f(x)，它恰好可以应用雅克比矩阵来计算得到：clear;symsx1x2x3x4x5;y=(1+(x1*x2+x3+x4^2)^(1/2)/(x5+x3^2))^(1/2);f=simple(jacobian(y)*[x1;x2;x3;x4;x5]);然后将得到的f表达式前面加上@(x1,x2,x3,x4,x5)：f=@(x1,x2,x3,x4,x5)-1/4*x3*(2*x1*x2*x3+x5+3*x3^2+2*x4^2*x3)/(x5+x3^2)/((x5+x3^2+(x1*x2+x3+x4^2)^(1/2))*(x5+x3^2))^(1/2)/(x1*x2+x3+x4^2)^(1/2);作为下面程序中的f即可，已经验证过结果仍然差不多。下面程序中的f仍为原程序所给的f。matlab蒙特卡罗方法程序（相关原理参见附件文献）：clear;f=@(x1,x2,x3,x4,x5)-(x3*(3*x3^2+2*x3*x4^2+2*x1*x2*x3+x5))./(4*((x5+(x4^2+x3+x1*x2)^(1/2)+x3^2)/(x3^2+x5))^(1/2)*(x3^2+x5)^2*(x4^2+x3+x1*x2)^(1/2));MIN=inf;LIMIT=10000;whileLIMIT>0 x(3)=4.5*rand+0.5;%将【0,1】区间上的随机数转化到【0.5,5】上的随机数，下面其余数类同x(4)=2*rand+1;x(5)=3*rand+1;x(1)=4.5*rand+0.5;x(2)=4.5*rand+0.5;whilex(1)+x(2)^2>=1x(4)=2*rand+1;x(5)=3*rand+1;x(1)=4.5*rand+0.5;x(2)=4.5*rand+0.5;ifx(1)+x(2)^2>=1x(4)=2*rand+1;x(5)=3*rand+1;x(1)=4.5*rand+0.5;x(2)=4.5*rand+0.5;whilex(1)+x(2)^210 x(1)=4.5*rand+0.5;x(2)=4.5*rand+0.5;end;x1=x(1);x2=x(2);x3=x(3);x4=x(4);x5=x(5);temp=f(x1,x2,x3,x4,x5);iftemp10时，再次令设：x(1)=4.5*rand+0.5;x(2)=4.5*rand+0.5;重新进行while循环，只有当条件上面条件不满足才完成内层的while循环，接着执行后面的语句。所以内层while循环完全等效于如下if判断：clear;f=@(x1,x2,x3,x4,x5)-(x3*(3*x3^2+2*x3*x4^2+2*x1*x2*x3+x5))./(4*((x5+(x4^2+x3+x1*x2)^(1/2)+x3^2)/(x3^2+x5))^(1/2)*(x3^2+x5)^2*(x4^2+x3+x1*x2)^(1/2));MIN=inf;LIMIT=10000;k=0;whileLIMIT>0 x(3)=4.5*rand+0.5;x(4)=2*rand+1;x(5)=3*rand+1;x(1)=4.5*rand+0.5;x(2)=4.5*rand+0.5;ifx(1)+x(2)^210continue;end;x1=x(1);x2=x(2);x3=x(3);x4=x(4);x5=x(5);temp=f(x1,x2,x3,x4,x5);iftemp