北航拓扑优化程序学习报告

拓扑优化的

99行程序学习报告4月19日2011

《结构优化设计》课程学习报告

任课教师：李书

一、前言：

在最近的结构优化设计课程上学习了O.Sigmund的《A 99 line topology optimization code written in Matlab》一文，对拓扑优化的理论原理与实际的计算机程序实现都有了一定的理解，文章主要是通过拓扑优化的原理来实现对简单结构的静力学问题的优化求解，而编写的代码仅有99行，包括36行的主程序，12行的OC优化准则代码，16行的网格过滤代码和35行的有限元分析代码。

自1988 年丹麦学者Bendsoe与美国学者Kikuchi提出基于均匀化方法的结构拓扑优化设计基本理论以来,均匀化方法应用到具有周期性结构的材料分析中,近几年该方法已经成为分析夹杂、纤维增强复合材料、混凝土材料等效模量,以及材料的细观结构拓扑优化常用的手段之一。其基本思想是在组成拓扑结构的材料中引入微结构，优化过程中以微结构的几何尺寸作为设计变量，以微结构的消长实现其增删，并产生介于由中间尺寸微结构组成的复合材料，从而实现了结构拓扑优化模型与尺寸优化模型的统一。文章就是通过均匀化的基础，结合拓扑结构优化的工程实际，以计算机模拟的方法将拓扑优化的一般过程呈现出来，有助于初涉拓扑优化的读者对拓扑优化有个基础的认识。

二、拓扑优化问题描述

为了简化问题的描述，文中假设设计域是简单的矩形形式，且在进行有限元离散的时候采用正方形单元对其进行离散。这样不仅便于进行单元离散和单元编号，也利于对结构进行几何外形的描述。

一般说来，基于指数逼近法的拓扑优化最小化的问题可作如下描述：

文中采用的对结构材料属性的描述是所谓的“指数逼近法”或者称为SIMP 逼近法，即（Solid Isotropic Material with Penalization带惩罚因子的各项同性材料模型法），该方法是拓扑优化中常用的变密度材料插值模型中最具代表性的一种。

基于SIMP格式的材料插值模型为

其中ρ表示材料的相对密度(0≤ ρ ≤1)，p代表惩罚因子，共同描述材料的力学属性。SIMP材料插值模型中，随着惩罚因子p取值的增大，对中间密度的惩罚程度越大，单元等效弹性模量逼近0或者E max的趋势也更明显。

三、Matlab代码实现

99行程序代码主要包括：主程序，OC优化准则代码，网格过滤代码，有限元分析代码等5个部分，而主函数的调用方式为top(nelx,nely,volfrac,penal,rmin)，共有5个输入参数，分别是：nelx，x方向（水平方向）单元划分数；nely，y方向（垂直方向）上的单元划分数；volfrac，结构体积保留分数；penal,是惩罚因子，即指数逼近模型中的指数因子，用于对材料属性按照指数模型进行逼近描述；rmin是网格过滤半径。

通过在调用主函数时改变上述5个输入参数，以及修改程序中的外载荷矩阵F和边界约束条件fixeddofs矩阵等代码段语句可以建立不同的输入模型，求解不同约束条件，不同外载荷下的拓扑优化解。

程序流程：

A.设计域的离散化

B.主循环，通过调用有限元分析子程序（Finete Element

subroutine）,返回位移列阵U

C.循环遍历所有单元，对离散的单元逐一进行单元节点编号；定

义出目标函数c和目标函数的变化率dc。

D.上一步计算得到的目标函数变化率dc作为输入参数之一，调

用网格过滤子程序check，得到输出dc(new)；

E.上一步计算得到的dc(new)作为输入参数之一，调用OC优化

准则子程序OC，得到输出xnew；

F.输出结果；

G.对结果进行可视化；

H.当xnew和xold之间达到要求的精度后，停止循环，结束迭代；

否则转到B继续循环。

程序流程示意图

四、算例及结果分析

这一部分主要通过若干简单的算例来研究各个输入参数、载荷施加方式以及边界约束条件对拓扑优化结果的影响，进而对拓扑优化问题加深理解和认识。

受力模型如上图所示，带双孔的悬臂梁，在右下端点处受垂直向下的单位载荷。长宽比为6:2，两孔圆心的位置分别为(L/4,W/2)，(3L/4,W/2) 处。

A.nelx和nely对结果的影响：

nelx*nely=60*20

nelx*nely=72*24

nelx*nely=90*30

nelx*nely=120*40

单元划分数对结果的影响

由以上结果可以大致看出，随着单元数的增加，用于分割材料的最小单元尺寸减小，拓扑优化后的图形显示上锯齿现象有所缓解，局部细节更加清晰；拓扑结构也有一定程度的改变，但总的说来还是有一致的趋势。从理论上说，拓扑优化过程中，在进行有限元离散时，进行离散的单元数越多，越接近于材料真实的无穷自由度情况，得到的有限元解（单元位移）越接近于真实水平，设计变量的增加可以对结构边界进行更精细尺度的描述，使拓扑结果中出现更多的细小分支结构，因此进行描述材料属性的参数也就越真实，在其他输入参数不变的情况下，优化结果越优。但另一方面，过于密集的网格划分容易造成拓扑优化结果中过多的分支结构，现过多的孔洞，使得结构的几何复杂性增加，结构制造成本提高，降低了结构局部刚度和强度。

而拓扑结果形式上的变化我认为是在单元数较少时，描述整体材料属性时的精度不足（一个单元所描述的区域较大，不够真实的反应该区域的力学特性变化）引起的误差，不得不在某些位置将刚度改变，从而改变了整体结构的刚度，使的整体的传力路线呈现一定程度上的差异，但是传力路线的趋势仍大致相同。

以上结果也是拓扑优化的网格依赖性的一种体现。

B.材料体积保留分数volfrac的影响：