工程结构优化设计—编程作业.doc

《工程结构优化设计—编程作业.doc》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、工程结构优化设计—编程作业1.用黄金分割法求方程在区间[-1,1]上的解。在MATLAB中没有专门的函数实现黄金分割法求解线性方程，可通过编写glodf.m函数实现黄金分割法求解，其代码如下:functionx=glodf(f,a,b,eps)if(nargin==3)eps=1.0e-4;endf1=subs(sym(f),findsym(sym(f)),a);f2=subs(sym(f),findsym(sym(f)),b);if(f1==0)x=a;endif(f2==0)x=b;endi

2、f(f1*f2>0)disp(’两端点函数值乘积大于0!’);return;elset1=a+(b-a)*0.382;t2=a+(b-a)*0.618;f_1=subs(sym(f),findsym(sym(f)),t1);f_2=subs(sym(f),findsym(sym(f)),t2);tol=abs(t1-t2);while(tol>eps)if(f_1*f_2<0)a=t1;b=t2;elsefa=subs(sym(f),findsym(sym(f)),a);if(f_1*fa>0)

3、a=t2;elseb=t1;endendt1=a+(b-a)*0.382;t2=a+(b-a)*0.618;f_1=subs(sym(f),findsym(sym(f)),t1);f_2=subs(sym(f),findsym(sym(f)),t2);tol=abs(t2-t1);endx=(t1+t2)/2;end其实现的MATLAB代码如下:>>clearall;>>tic>>x=glodf('x^3-3*x+1',-1,1,1e-6)x=0.34731.用抛物线法求解的极值点，取值范围[0

4、.5，1]，允许误差为1e-6。在MATLAB中没有专门的函数实现抛物线法求解线性方程，可通过编写paowuxian.m函数实现抛物线法求解，其代码如下:functionx=paowuxian(f,x0,x1,x2,e)ifnargin<4,e=1e-4;endx=x2;y=x1;z=x0;whileabs(x-y)>eh1=y-z;h2=x-y;c1=(feval(f,y)-feval(f,z))/h1;c2=(feval(f,x)-feval(f,y))/h2;d=(c1-c2)/(h2+h

5、1);w=c2+h2*d;xi=x-(2*feval(f,x))/(w+(w/abs(w))*sqrt(w^2-4*feval(f,x)*d));z=y;y=x;x=xi;end其实现的MATLAB代码如下:>>fun=inline('x^3+2*x-0.2');>>paowuxian(fun,0,0.5,1,1e-6),formatshortans=0.09952.用齿形法对下图的结构进行杆件截面优化设计：3.节点处作用两种工况：（1）P1=2000kN，P2=0kN；（2）P1=0kN，P2

6、=2000kN。各杆采用同一材料制成，弹性模量E为常数，容许应力分别为：[σ+]=2╳107kPa，[σ-]=1.5╳107kPa。此结构为对称结构，工况也对称，所以A1=A3，计算时只取两个变量A1和A2，一种工况P1=2000kN，P2=0kN。目标函数取为结构总重量：W(22A1A2)l在MATLAB中没有专门的函数实现齿形法求解线性方程，可通过编写Zigzag.m函数实现齿形法求解，其代码如下:function[newA1,newA2,newA3,W,exitflag]=Zigzag

7、(A1,A2,maxstep)ifnargin<3maxstep=96;endA=[A1;A2];U=[(1e-4)*(sqrt(2)*A1+A2)/(sqrt(2)*A1^2+2*A1*A2);(1e-4)*(sqrt(2)*A1)/(sqrt(2)*A1^2+2*A1*A2)];w=sqrt(2)*2*A1+A2;w1=w;w2=w;newA1=A1;newA2=A2;newA3=A1;symspl;W=vpa(w*p*l);A11=A1;A13=A1;A22=A2;A24=A2;exitf

8、lag=-1;fork=1:maxstepifw2>w1newA1=A11;newA2=A22;newA3=A11;W=vpa(w1*p*l,4);exitflag=1;return;endu1=max(U);A=A.*u1;A11=A(1);A22=A(2);U=U./u1;U=roundn(U,-4);w=sqrt(2)*2*A11+A22;w1=w;ifw1>w2newA1=A13;newA2=A24;newA3=A13;W=vpa(w2*p*l,4);exitflag=0;return;