现代设计理论及方法优化设计实验报告

西安交通大学实验报告

课程名称：现代设计理论与方法实验名称：优化设计上机实验

学院：实验日期：

班级：姓名：学号：

一、实验要求

1. 采用MA TLAB等编程语言，编写优化程序，计算优化结果；

2. 完成大作业书面报告，对每个题目进行分析建模，包括：

①设计变量的选择；

②优化目标函数的确定；

③约束条件的确定。

二、优化分析

1. 镗刀杆（销轴）结构参数优化

①设计变量的选择

题目要求“试在满足强度、刚度条件下，设计一个用料最省的方案”，即在满足性能要求的前提下，使设计方案的质量（体积）最小。

最直接的思路为，控制长度L和直径d最小。而根据条件分析，亦可通过改变截面形状（改变轮廓形状、使截面空心等）、改变不同L处截面形状等复杂的空间质量分布模式等，来到达最优的目的。

为便于分析，此处选择设计变量为刀杆直径d、长度L（实际可直接取最小值）为设计变量。

②优化目标函数的确定

刀杆用料最省，即体积最小：

V=1

4

πd2L→min

设

x1=d,x2=L 则目标函数为

min f X=V=1

πd2l=

1

πx12x2=0.785x12x2

③约束条件的确定

根据材料力学知识，应有：

σmax<σ

τmax<τ

f max

σmax=32

πd3

M Z2+0.75M n2≤[σ]

τmax=16M n

πd3

≤[τ]

f max=PL3

≤[f]

L≥80

2. 梯形截面管道参数优化

①设计变量的选择

题目要求设计管道的参数，即取管道截面的高度h、底边长度c、底边与侧边所夹锐角θ为设计变量。

②优化目标函数的确定

已知管道内液体的流速与管道截面的周长S的倒数成正比例关系，当液体流速最大时，则管道截面周长最小：

S=2h

+2c+

2h

→min

③约束条件的确定

管道的截面面积一定，则有：

A=h c+h

=64516

管道截面为梯形，则有：

0<θ<

π

h>0

c>0

3. 厂址选择问题

①设计变量的选择

设位于i地的工厂从j地购买的原材料质量为B ij，位于i地的工厂向j地销售的产品质量为S ij。即设计变量为B ij ,S ij (i,j=A,B,C)

②优化目标函数的确定

记i、j两地距离为L ij(i,j=A,B,C)，记i地生产费用为P i(i=A,B,C)题目要求总费用最小：

C=5000B ij L ij+5000S ij L ij+S i P i→min

③约束条件的确定

从各地购买的原材料总和小于等于其产量

在各地销售的产品总和小于等于其销售量

(S ij?1 3

j

B ij)<0

S ij

j

?5<0

4.桁架结构优化

①设计变量的选择

题目要求“对该桁架优化设计使其质量最轻”，题目给出的限制条件较少。同样可认为优化算法中，可以改变每根杆的截面形状、改变每根杆上不同位置处的截面形状，或者仅改变每根杆的直径，使得体积最小，即质量最轻。

所以，取每根杆的直径d i(i=1,2,3,…,10)为设计变量。

②优化目标函数的确定

桁架质量最轻，即总体积最小。忽略连接处体积变化，则有：

V=

1

πL i d i2→min 10

i=1

③约束条件的确定

各根杆在同一节点处位移相等

相应杆在节点1、6处位移为零

σi max<σ (i=1,2,3, (10)

5.半圆弧拱结构截面优化

①设计变量的选择

题目要求“求截面应力小于20MPa的最优设计”。认为最优设计为用料最少，假定各处截面形状一致，即求截面面积最小。同时，假定截面形状固定为题目所示的“口”字形，则可以取设计变量为高度H、宽度B、壁厚T。

②优化目标函数的确定

用料最少，即截面面积最小：

S=BH?B?2T H?2T=2TB+2TH?4T2→min

③约束条件的确定

截面应力小于20MPa，即：

σmax<20 MPa

截面为“口”字形：

B>0,H>0,T>0

B>2T,H>2T

6. 正方形板拓扑优化

①设计变量的选择

题目要求“试用拓扑优化对其体积进行优化，使其质量最小，刚度最大”。则设计变量应为正方形板的质量在体积上的分布。考虑到其厚度相对宽度极小，可认为设计变量为其质量在面积上的分布X。

②优化目标函数的确定

设正方形板的质量为f M X，正方形板的刚度f S X，则有：

f M X→min

?f S X→min

所以优化目标函数可为两个函数的加权函数：

g f M,f S→min

③约束条件的确定

正方形板左边的变形为0

质量在边界范围内连续分布

σmax<[σ]

附：参考程序

1. 镗刀杆（销轴）结构参数优化

①优化目标函数：

function f=taget(x)

f=pi/4*x(1)^2*x(2);

②约束条件：

function[c,ceq]=mycon(x)

c(1)=((10000*x(2))^2+(0.75*100000)^2)^0.5/(pi*x(1)^3/32)-120;

c(2)=100000/(pi*x(1)^3/16)-80;

c(3)=10000*x(2)^3/(3*200000*pi*x(1)^4/32)-0.1;

ceq=[];

③主程序：

A=[0 -1];

b=[-80];

Aeq=[];

beq=[];

lb=[0;0];

ub=[inf;inf];

x0=[1;80];

options=optimset('LargeScale','off','display','iter');

[x,fval]=fmincon(@taget1,x0,A,b,[],[],lb,ub,@mycon1,options)

④计算求解：

2. 梯形截面管道参数优化

①主程序：

lb=[0;0;0];

A=[0 0 1];

b=[pi/2];

ub=[inf;inf;inf];

x0=[64516^0.5;64516^0.5;pi/2];

[x,fval]=fmincon(@taget2,x0,A,b,[],[],lb,ub,@mycon2)

②优化目标函数：

function f=taget2(x)

f=-1/(x(1)+2*x(2)/sin(x(3))+x(1)+2*x(2)/tan(x(3)));

③约束条件：

function[c,ceq]=mycon2(x)

ceq=(x(1)+x(1)+2*x(2)/tan(x(3)))*x(2)/2-64516;

c=[];

④计算求解：

3. 厂址选择问题

①优化目标函数：

function f=fun01(x)

f=75*x(1)+50*x(2)+75*x(3)+100*x(4)+50*x(5)+100*x(6)+150*x(7)+240*x(8)+210*x(9)+1 20*x(10)+160*x(11)+220*x(12);

②约束函数：

function[c,ceq]=con01(x)

c(1)=x(1)+x(2)+3*x(7)+3*x(8)-x(3)-x(5)-20; %A地的原料不能超过20万t

c(2)=x(3)+x(4)+3*x(9)+3*x(10)-x(1)-x(6)-16; %B地的原料不能超过16万t

c(3)=x(5)+x(6)+3*x(11)+3*x(12)-x(2)-x(4)-24; %C地的原料不能超过24万t

c(4)=x(9)+x(10)-5; %B地建厂规模不能超过5万t

ceq(1)=x(7)+x(9)+x(11)-7; %A地年消产品7万t

ceq(2)=x(8)+x(10)+x(12)-13; %B地年消产品13万t

③主程序：

x0=[0,0,0,0,0,0,0,0,0,0,0,0,] %初值

lb=[0,0,0,0,0,0,0,0,0,0,0,0,] %下限

ub=[Inf,Inf,Inf,Inf,Inf,Inf,Inf,Inf,Inf,Inf,Inf,Inf,] %上限

[x,fval]=fmincon(@fun01,x0,[],[],[],[],lb,ub,@con01) %调用函数fmincon，最优求解

④计算求解：

4.桁架结构优化

①空间转置矩阵

function T=TransformMatrix(ie)

global gElementgNode

gElement=[5 6;4 5;1 2;2 3;2 5;3 4;2 6;1 5;3 5;2 4];

gNode=[0 0;1 0;2 0;2 1;1 1;0 1];

xi = gNode( gElement ( ie, 1 ), 1 ) ;

yi = gNode( gElement ( ie, 1 ), 2 ) ;

xj = gNode( gElement ( ie, 2 ), 1 ) ;

yj = gNode( gElement ( ie, 2 ), 2 ) ;

L = sqrt( (xj -xi )^2 + (yj -yi )^2 ) ;

c = (xj -xi )/ L ;

s = (yj -yi )/ L ;

T=[ c -s 0 0;

s c 0 0;

0 0 c -s;

0 0 s c ];

return

②求解单刚阵

function k=StiffnessMatrix( ie )

syms A1A2A3A4A5A6A7A8A9A10;

gMaterial=[A1 A2 A3 A4 A5 A6 A7 A8 A9 A10];

global gNodegElement

k=zeros( 4,4 ) ;

E = 1 ;

A = gMaterial(ie);

xi = gNode(gElement (ie, 1 ), 1 ) ;

yi = gNode(gElement (ie, 1 ), 2 ) ;

xj = gNode(gElement (ie, 2 ), 1 ) ;

yj = gNode(gElement (ie, 2 ), 2 ) ;

L = ( (xj -xi )^2 + (yj -yi )^2 )^(1/ 2) ;

k=[ E*A/ L 0 -E*A/ L 0 ;0 0 0 0;-E*A/ L 0 E*A/ L 0;0 0 0 0] ; T = TransformMatrix(ie ) ;

k = T*k*transpose(T);

return

③合成总刚阵求解

clear;

gK=sym('gK',[12 12]);

for i =1:1:12

for j =1:1:12

gK(i,j)=0;

end

end

for ie=1:10

T=TransformMatrix(ie);

k=StiffnessMatrix(ie);

global gElement

for i =1:1:2

for j =1:1:2

for p=1:1:2

for q =1:1:2

m = (i -1)*2+p ;

n = (j -1)*2+q ;

M = (gElement (ie,i)-1)*2+p ;

N = (gElement (ie,j)-1)*2+q ;

k = StiffnessMatrix(ie);

gK(M,N) = gK(M,N)+ k(m,n);

end

end

end

end

end

gK

K=gK(3:10,3:10)

F=[0;-4000;0;-4000;0;0;0;0];

Q=K^-1

U=Q*F

syms A1A2A3A4A5A6A7A8A9A10;

ang=[180 180 0 0 90 90 135 45 135 45]*pi/180;

P(1,1)=[cos(ang(1)) sin(ang(1))]*[-U(7);-U(8)]*1/A1

P(2,1)=[cos(ang(2)) sin(ang(2))]*[U(7)-U(5);U(8)-U(6)]*1/A2;

P(3,1)=[cos(ang(3)) sin(ang(3))]*[U(1);U(2)]*1/A3;

P(4,1)=[cos(ang(4)) sin(ang(4))]*[U(3)-U(1);U(4)-U(2)]*1/A4;

P(5,1)=[cos(ang(5)) sin(ang(5))]*[U(7)-U(1);U(8)-U(2)]*1/A5;

P(6,1)=[cos(ang(6)) sin(ang(6))]*[U(5)-U(3);U(6)-U(4)]*1/A6;

P(7,1)=[cos(ang(7)) sin(ang(7))]*[-U(1);-U(2)]*1.414/A7;

P(8,1)=[cos(ang(8)) sin(ang(8))]*[U(7);U(8)]*1.414/A8;

P(9,1)=[cos(ang(9)) sin(ang(9))]*[U(7)-U(3);U(8)-U(4)]*1.414/A9;

P(10,1)=[cos(ang(10)) sin(ang(10))]*[U(5)-U(1);U(6)-U(2)]*1.414/A10;

P

I=[A1;A2;A3;A4;A5;A6;A7;A8;A9;A10]*200;

P=P-I

④优化目标函数：

function f=taget3(x)

f=x(1)+x(2)+x(3)+x(4)+x(5)+x(6)+1.414*x(7)+1.414*x(8)+1.414*x(9)+1.414*x(10);

⑤计算求解：

5. 半圆弧拱结构截面优化

①主程序：

finish

/clear

l=0.3 !截面长度

h=0.3 !截面宽度

t=0.04 !截面壁厚

/PREP7

et,1,188 !定义单元类型beam188

!定义材料参数

mp,ex,1,3e11 !弹性模量

mp,nuxy,1,0.17 !泊松比

mp,dens,1,2500 !密度

!定义截面类型（回字形）SECTYPE, 1, BEAM,hrec,,0 SECOFFSET, CENT

SECDATA,l,h,t,t,t,t

!建立模型

k,1,0,20,0

k,2,20,0,0

k,3,-20,0,0

larc,2,3,1,20

!划分网格

lsel,all

lesize,all,0.1

mat,1

type,1

secn,1

lmesh,1

alls

fini

!求解

/SOLU

ANTYPE,0

!加约束

d,1,all

d,2,all

!加节点集中载荷

nf=NODE(0,20,0)

f,nf,fy,-2e5

solve

fini

!后处理，观察应力云图/post1

/eshape,1

plesol,s,EQV,0,1

②计算求解：

6. 正方形板拓扑优化

①材料设定：

②划分网格：

③添加约束条件：

④添加载荷：

⑤计算求解：

（版权归上传者所有）