工程结构优化设计编程作业.docx
- 文档编号:29098701
- 上传时间:2023-07-20
- 格式:DOCX
- 页数:16
- 大小:25.72KB
工程结构优化设计编程作业.docx
《工程结构优化设计编程作业.docx》由会员分享,可在线阅读,更多相关《工程结构优化设计编程作业.docx(16页珍藏版)》请在冰豆网上搜索。
工程结构优化设计编程作业
工程结构优化设计—编程作业
1.用黄金分割法求方程
在区间[-1,1]上的解。
在MATLAB中没有专门的函数实现黄金分割法求解线性方程,可通过编写glodf.m函数实现黄金分割法求解,其代码如下:
functionx=glodf(f,a,b,eps)
if(nargin==3)
eps=1.0e-4;
end
f1=subs(sym(f),findsym(sym(f)),a);
f2=subs(sym(f),findsym(sym(f)),b);
if(f1==0)
x=a;
end
if(f2==0)
x=b;
end
if(f1*f2>0)
disp(’两端点函数值乘积大于0!
’);
return;
else
t1=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;
else
fa=subs(sym(f),findsym(sym(f)),a);
if(f_1*fa>0)
a=t2;
else
b=t1;
end
end
t1=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);
end
x=(t1+t2)/2;
end
其实现的MATLAB代码如下:
>>clearall;
>>tic
>>x=glodf('x^3-3*x+1',-1,1,1e-6)
x=
0.3473
2.用抛物线法求解
的极值点,取值范围[0.5,1],允许误差为1e-6。
在MATLAB中没有专门的函数实现抛物线法求解线性方程,可通过编写paowuxian.m函数实现抛物线法求解,其代码如下:
functionx=paowuxian(f,x0,x1,x2,e)
ifnargin<4,e=1e-4;
end
x=x2;y=x1;z=x0;
whileabs(x-y)>e
h1=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+h1);
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),formatshort
ans=
0.0995
3.用齿形法对下图的结构进行杆件截面优化设计:
节点处作用两种工况:
(1)P1=2000kN,P2=0kN;
(2)P1=0kN,P2=2000kN。
各杆采用同一材料制成,弹性模量E为常数,容许应力分别为:
[σ+]=2╳107kPa,[σ-]=1.5╳107kPa。
此结构为对称结构,工况也对称,所以A1=A3,计算时只取两个变量A1和A2,一种工况P1=2000kN,P2=0kN。
目标函数取为结构总重量:
在MATLAB中没有专门的函数实现齿形法求解线性方程,可通过编写Zigzag.m函数实现齿形法求解,其代码如下:
function[newA1,newA2,newA3,W,exitflag]=Zigzag(A1,A2,maxstep)
ifnargin<3
maxstep=96;
end
A=[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;
exitflag=-1;
fork=1:
maxstep
ifw2>w1
newA1=A11;
newA2=A22;
newA3=A11;
W=vpa(w1*p*l,4);
exitflag=1;
return;
end
u1=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>w2
newA1=A13;
newA2=A24;
newA3=A13;
W=vpa(w2*p*l,4);
exitflag=0;
return;
end
A=A.*U;
A1=A
(1);
A2=A
(2);
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)];
u1=max(U);
A=A.*u1;
A13=A
(1);
A24=A
(2);
U=U./u1;
U=roundn(U,-4);
w=sqrt
(2)*2*A13+A24;
w2=w;
A=A.*U;
A1=A
(1);
A2=A
(2);
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)];
end
end
其实现的MATLAB代码如下:
>>clear;
A1=1e-4;
A2=1e-4;
[Anew1,Anew2,Anew3,W,exitflag]=Zigzag(A1,A2)
Anew1=
7.7346e-05
Anew2=
4.5309e-05
Anew3=
7.7346e-05
W=
0.0002641*l*p
exitflag=
0
4.用修改齿形法对结构进行杆件截面优化设计。
在MATLAB中没有专门的函数实现修改齿形法求解线性方程,可通过编写IZigzag.m函数实现修改齿形法求解,其代码如下:
function[newA1,newA2,newA3,W,exitflag]=IZigzag(A1,A2,v,maxstep)
ifnargin<5
maxstep=96;
end
A=[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;
exitflag=-1;
fork=1:
maxstep
ifw2>w1
newA1=A11;
newA2=A22;
newA3=A11;
W=vpa(w1*p*l,4);
exitflag=1;
return;
end
u1=max(U);
A=A.*u1;
A11=A
(1);
A22=A
(2);
U=U./u1;
x1=U
(1);
x2=U
(2);
y1=1-v*(1-x1);
y2=1-v*(1-x2);
U=[y1;y2];
U=roundn(U,-4);
w=sqrt
(2)*2*A11+A22;
w1=w;
ifw1>w2
newA1=A13;
newA2=A24;
newA3=A13;
W=vpa(w2*p*l,4);
exitflag=0;
return;
end
A=A.*U;
A1=A
(1);
A2=A
(2);
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)];
u1=max(U);
A=A.*u1;
A13=A
(1);
A24=A
(2);
U=U./u1;
x1=U
(1);
x2=U
(2);
y1=1-v*(1-x1);
y2=1-v*(1-x2);
U=[y1;y2];
U=roundn(U,-4);
w=sqrt
(2)*2*A13+A24;
w2=w;
A=A.*U;
A1=A
(1);
A2=A
(2);
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)];
end
end
其实现的MATLAB代码如下:
>>clear;
A1=1e-4;
A2=1e-4;
v=0.5;
[Anew1,Anew2,Anew3,W,exitflag]=IZigzag(A1,A2,v)
Anew1=
7.8168e-05
Anew2=
4.2839e-05
Anew3=
7.8168e-05
W=
0.0002639*l*p
exitflag=
0
5.用满应力法对结构进行杆件截面优化设计
在MATLAB中没有专门的函数实现满应力法求解线性方程,可通过编写Full_Stressed.m函数实现满应力法求解,其代码如下:
function[newA1,newA2,newA3,exitflag]=Full_Stressed(A1,A2,wucha,maxstep)
ifnargin<4
maxstep=96;
end
ifnargin<3
wucha=1e-6;
end
A=[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)];
newA1=A1;
newA2=A2;
newA3=A1;
exitflag=-1;
fork=1:
maxstep
m=U
(1);
h=1-m;
n=U
(2);
j=1-n;
ifabs(h) newA1=A1; newA2=A2; newA3=A1; exitflag=0; return; end 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)]; U=roundn(U,-4); A=A.*U; A1=A (1); A2=A (2); A=[A1;A2]; end newA1=A1; newA2=A2; newA3=A1; end 其实现的MATLAB代码如下: >>clearall;closeall; A1=1e-4; A2=1e-4; [newA1,newA2,newA3,exitflag]=Full_Stressed(A1,A2) newA1= 9.8988e-05 newA2= 1.4373e-06 newA3= 9.8988e-05 exitflag= -1 6.用加松弛因子满应力法对结构进行杆件截面优化设计 在MATLAB中没有专门的函数实现加松弛因子满应力法求解线性方程,可通过编写SRF_Full_Stressed.m函数实现满应力法求解,其代码如下: function[newA1,newA2,newA3,exitflag]=SRF_Full_Stressed(A1,A2,B,wucha,maxstep) ifnargin<5 maxstep=96; end ifnargin<3 wucha=1e-6; end A=[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)]; newA1=A1; newA2=A2; newA3=A1; exitflag=-1; fork=1: maxstep m=U (1); h=1-m; n=U (2); j=1-n; ifabs(h) newA1=A1; newA2=A2; newA3=A1; exitflag=0; return; end 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)]; U=roundn(U,-4); U=U.^B; A=A.*U; A1=A (1); A2=A (2); A=[A1;A2]; end newA1=A1; newA2=A2; newA3=A1; end 其实现的MATLAB代码如下: >>clearall;closeall; A1=1e-4; A2=1e-4; B=1.1; wucha=0.05; [newA1,newA2,newA3,exitflag]=SRF_Full_Stressed(A1,A2,B,wucha) newA1= 9.5527e-05 newA2= 6.3262e-06 newA3= 9.5527e-05 exitflag= 0
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 工程 结构 优化 设计 编程 作业
![提示](https://static.bdocx.com/images/bang_tan.gif)