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

拓扑优化的99行程序学写报告余艳辉 SY1005525拓扑优化的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或者Emax的趋势也更明显。三、 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）,返回位移列阵UC. 循环遍历所有单元，对离散的单元逐一进行单元节点编号；定义出目标函数c和目标函数的变化率dc。D. 上一步计算得到的目标函数变化率dc作为输入参数之一，调用网格过滤子程序check，得到输出dc(new)；E. 上一步计算得到的dc(new)作为输入参数之一，调用OC优化准则子程序OC，得到输出xnew；F. 输出结果；G. 对结果进行可视化；H. 当xnew和xold之间达到要求的精度后，停止循环，结束迭代；否则转到B继续循环。设计域的离散化及程序参数的初始化主循环调用有限元分析求解子程序FE U作为输入定义目标函数，计算目标函数变化率dc 调用网格过滤子程序checkdc作为输入 dc(new)作为输入调用OC优化准则子程序OC否是否满足精度要求 输出x(new)是结束循环程序流程示意图四、 算例及结果分析这一部分主要通过若干简单的算例来研究各个输入参数、载荷施加方式以及边界约束条件对拓扑优化结果的影响，进而对拓扑优化问题加深理解和认识。受力模型如上图所示，带双孔的悬臂梁，在右下端点处受垂直向下的单位载荷。 长宽比为6:2，两孔圆心的位置分别为(L/4,W/2)，(3L/4,W/2) 处。A. nelx和nely对结果的影响：nelx*nely=60*20 nelx*nely=72*24nelx*nely=90*30nelx*nely=120*40单元划分数对结果的影响由以上结果可以大致看出，随着单元数的增加，用于分割材料的最小单元尺寸减小，拓扑优化后的图形显示上锯齿现象有所缓解，局部细节更加清晰；拓扑结构也有一定程度的改变，但总的说来还是有一致的趋势。从理论上说，拓扑优化过程中，在进行有限元离散时，进行离散的单元数越多，越接近于材料真实的无穷自由度情况，得到的有限元解（单元位移）越接近于真实水平，设计变量的增加可以对结构边界进行更精细尺度的描述，使拓扑结果中出现更多的细小分支结构，因此进行描述材料属性的参数也就越真实，在其他输入参数不变的情况下，优化结果越优。但另一方面，过于密集的网格划分容易造成拓扑优化结果中过多的分支结构，现过多的孔洞，使得结构的几何复杂性增加，结构制造成本提高，降低了结构局部刚度和强度。而拓扑结果形式上的变化我认为是在单元数较少时，描述整体材料属性时的精度不足（一个单元所描述的区域较大，不够真实的反应该区域的力学特性变化）引起的误差，不得不在某些位置将刚度改变，从而改变了整体结构的刚度，使的整体的传力路线呈现一定程度上的差异，但是传力路线的趋势仍大致相同。以上结果也是拓扑优化的网格依赖性的一种体现。B. 材料体积保留分数volfrac的影响： volfrac=0.4 volfrac=0.5volfrac=0.6 volfrac=0.7 volfrac=0.8由以上的优化结果可以看出，随着体积保留分数volfrac的增大，优化结果中保留的材料也逐渐增多。这很容易理解，体积保留分数的增大，意味着保留的材料也增多，在基本的传力路线的基础上增加了一些额外的材料。另一方面，过小的volfrac取值，可能造成改变结构传力的路线，或者迭代多次而收敛速度却较慢。volfrac对优化结果的影响从程序中也可以体现：volfrac在程序中主要作为OC优化准则中利用二分算法搜索拉格朗日乘子的判定准则。C. 惩罚因子P对结果的影响：P=1.5p=2.5P=3.0由优化结果可以看到：随着惩罚因子取值的增大，呈现出的拓扑结果呈现日趋清晰的趋势，这是因为在SIMP材料插值模型中，参数p较小时，对中间密度的惩罚程度较小，或者说惩罚效果不够明显，因此呈现的结果会有模糊的灰色区域；随着参数p取值增大，对中间密度的惩罚程度也越大，单元等效弹性模量逼近0或者Emax的趋势也更明显，因此，优化结果也愈加清晰。另一方面，取值越大，迭代次数也越多，计算量增大，甚至无解（例如，在该算例中，取p=4.5时，主循环的循环中止条件change值=0.200长时间无变化，目标函数的变化也很小，也呈现震荡趋势，可能是解收敛太慢的原因）。总的说来，各个取值情况下的优化结构大致相似。因此，在实际的拓扑优化中，要注意选择适宜的惩罚因子，达到平衡迭代时间和较优的优化效果。 D. 过滤半径rmin对拓扑结果的影响 r=0.8r=1.0r=1.5r=2.0r=3.0r=4.5文章中提到过：网格过滤半径rmin的选取如果小于1的话，被过滤的目标函数变化率将会和原始的变化率相等，使得过滤无效。在结果上表现出来的是棋盘格现象，这在工程上不可能实现的，对工程实践没有意义。从上面的优化结果分析，随着rmin取值的增大，过滤效果随着其增大而增强，在局部出现模糊的灰色区域。在实际的操作中，要综合考虑结果的精度要求以及计算时间等各方面因素，综合考量，选取适宜的网格过滤半径，这样才能够得到既满足设计优化要求，计算量又不会太大的优化结果。E. 边界约束条件和加载方式对拓扑结果的影响 1) 悬臂梁上端受均布载荷2) 双工况情况： 3) 左右两侧悬臂，上端受均布载荷的情况： 从上面的结果可以看出，随着边界约束条件和加载方式的变化，拓扑优化结果发生明显变化，这是因为外部约束方式和加载方式的变化引起结构内部传力路线的变化。F. 工程实际应用 电动机轴承支承座的优化处理问题问题简要描述： 电动机轴承支座部分受力问题，孔的圆心在矩形上水平边的中点处，承载方式简化为对其有点接触的集中力，垂直向下，约束方式为简支。试求在该带孔的矩形板上的拓扑优化结构。利用程序将力学模型修改完毕后，运行程序得到优化结果如上右图所示，得到的大致轮廓与实际中的电动机轴承支座形状基本一致。五、 结论及心得 研究分析了程序输入参数nelx,nely,volfrac,p,rmin对优化结果的影响，通过改变加载矩阵F和边界约束条件，实现了对不同加载方式和不同边界约束下的简单结构静力学问题的优化。结果看到，输入参数的变化对优化结果影响显著，但是不同的输入参数变化对结果影响也不同。在实际的应用中应综合考虑，选取适宜的优化参数。 通过对 A 99 line topology optimization code written in Matlab 一文的学习，对拓扑优化的理论原理以及实际的计算机方法实现都有了一定的认识，了解了拓扑优化的一般步骤，拓扑优化问题中的理论基础等。通过对99行程序的学习，结合实际的上机操作，得到了一些简单的结构问题的拓扑优化求解结果，并且分析对比了该程序主要输入参数的变化对拓扑优化结果的影响；另一方面，由于文章仅是用于教学目的，为了处理问题的简便，在