英语原文共 8 页，剩余内容已隐藏，支付完成后下载完整资料
用Matlab编写的99行拓扑优化代码
O. Sigmund
摘要：本文提出了一个简洁的可以使静载结构柔度最小化的拓扑优化程序的Matlab实现方法。Matlab输入行的总数为99行，其中包括了优化程序和有限元子程序。在这99行中有36行用于主程序，12行用于基于最优性准则的优化程序，16行用于网格独立滤波器，35行用于有限元代码。事实上，如果不包含与输出和有限元分析相关的注释行和行，求解一个适定的拓扑优化问题只需要49条Matlab的输入行。通过增加三条附加行，该程序还可以解决多个负载情况下的问题。本代码是为了教学目的而设计的。完整的Matlab代码在附录中给出，同时也可以从网站http://www.topopt.dtu.dk上下载。
关键词：拓扑优化，教学，最优性准则，万维网，Matlab代码
目录
1 背景介绍
本文给出的Matlab程序是为工程设计教学而设计的。拓扑优化领域的学生和新手可以从网站http://www.topopt.dtu.dk上下载代码。该代码可用于结构优化的课程，课上学生可以被安排进行扩展，如多个荷载工况、可替代的网格独立方案、无源区域等。另一种可能是利用这个程序来培养学生的直觉，从而进行优化设计。可以要求高等的学生猜测在给定边界条件和体积分数下的最佳拓扑，然后与程序运行得到的正确的最佳拓扑做对比。
在文献中，可以找到许多解决拓扑优化问题的方法。在原文Bendsoslash;e和Kikuchi (1988)中，基于研究解的存在性，所谓的微观结构或均匀化方法被使用。
基于均匀化的方法已经在许多文献中被采用，但其缺点是，如果问题没有被解决（对于不符合的问题），确定和评估最佳微观结构及其方向是很复杂的，而且，由于没有确定的长度尺度与微观结构相关，因此无法构建所得到的结构。然而，均匀化方法在拓扑优化中仍然很重要，因为在某种意义上它可以为结构的理论性能提供界限。
拓扑优化的另一种方法是所谓的“幂律方法”或SIMP方法（含有惩罚因子的固体各向同性材料）（Bendsoslash;e 1989；Zhou和Rozvany 1991；Mlejnek 1992）。这里，假设材料特性为用于离散设计域的每个单元内的常数，变量为单元相对密度。将材料特性建模设置为相对材料密度提高到一定的幂次，乘以固体材料的材料特性。这种方法已被批评，因为有人认为，没有幂律插值所描述的属性的物理材料存在。然而，Bendsoslash;e和Sigmund（1999）最近的一篇论文证明，只要满足关于幂的简单条件（例如，泊松比等于时pge;3），幂律方法在物理上是允许的。为了确保解决方案的存在，幂律方法必须与周长约束、梯度约束或滤波技术相结合（见Sigmund和Petersson 1998）。拓扑优化的幂律方法已被应用于多约束、多物理和多材料的问题。尽管上述方法的求解是基于数学规划技术和连续设计变量，但已有许多文献将拓扑优化问题作为一个整数问题来求解。Beckers（1999）成功地使用双重方法完全解决了大规模的柔度最小化问题，但其他基于遗传算法或半随机方法的其他方法需要数千个函数评估，即使是对于少量元素来说，这必须被视为是不切实际的。
除了上述方法都能解决定义明确的问题（例如柔度最小化）之外，一些启发式的或基于直觉的方法已经被证明可以降低柔度或其他目标函数。这些方法中有所谓的进化设计方法（参见Xie和Steven 1997；Baumgartner等人 1992）。除了非常容易理解和实现（至少对于柔度最小化的情况），进化方法的主要动机似乎是基于数学或连续变量的方法“涉及一些复杂的微积分运算和数学编程”（引自Li等人，1999），并且它们包含“一些复杂性的数学方法”（引自Zhao等人，1998），而进化方法“利用了强大的计算技术和自然界进化过程的直观概念”（引自Li等人，1999）。但是有两件事可以反驳这个。首先，进化方法本身变得很复杂，因为要考虑到比柔度最小化更复杂的目标。第二，如本文所示，柔度最小化的“基于数学的”方法也易于实现，并且在计算上同样有效。此外，基于数学规划的方法可以很容易地扩展到其他非柔度目标，例如非自伴问题和多物理量问题以及多约束问题（例如Sigmund 1999）。对这类问题的进化方法的扩展似乎更值得怀疑。
完整的Matlab代码在附录中给出。本文的剩余部分包括了优化问题的定义和讨论（第2章节），关于Matlab实现的评价（第3章节），关于扩展的讨论（第4章节）以及结论（第5章节）。
2 拓扑优化问题
为了简化Matlab程序，引入一些简化方法。首先，假设设计域为矩形，并用方形有限元进行离散。这样，单元和节点的编号很简单（从左上角开始逐列编号），结构的长宽比由水平方向（nelx）和垂直方向（nely）的单元比给出。一个基于幂律方法的拓扑优化问题，其目标是柔度最小化，可以写为
其中U和F分别为整体位移和力矢量，K为整体刚度矩阵，和分别为单元位移矢量和刚度矩阵，x为设计变量矢量，为最小相对密度矢量（为避免奇点，非零），N（=nelxtimes;nely）是用于离散设计域的元素数，p是惩罚幂（通常p=3），和分别是材料体积和设计域体积，f（volfrac）是规定的体积分数。
优化问题（1）可以用几种不同的方法来求解，如最优性准则（OC）方法、序列线性规划（SLP）方法或移动渐近线法（MMA by Svanberg 1987）等。为了简单方便，我们将在这里使用标准的OC方法。
根据Bendsoslash;e（1995），设计变量的启发式更新方案可以表述为
其中m（move）是正的移动极限，eta;（=1/2）是数值阻尼系数，从最优性条件中可以发现为
式中，lambda;是拉格朗日乘子，可通过双分段算法找到。
目标函数的灵敏度为
关于最优性标准方法的推导和实现的更多细节，请参考文献（例如Bendsoslash;e 1995）。
为了确保拓扑优化问题（1）的解的存在，必须引入对所得到的设计结果的某种限制（参见Sigmund和Peelsson 1998）。这里我们使用滤波技术（Sigmund 1994，1997）。必须强调的是，该滤波器还没有被证明能够确保解的存在，但是作者的许多应用已经证明了滤波器在实践中产生了网格独立设计。
网格独立过滤器通过修改单元灵敏度来工作，如下所示：
卷积算子（权重因子）写为
这里运算符（e，f）被定义为单元e的中心与单元f的中心之间的距离。卷积算子在滤波区域外为零，卷积算子随着单元f的距离线性衰减。在最优性准则更新（3）中使用修改后的灵敏度（5）代替原来的灵敏度（4）。
3 Matlab的实现
Matlab代码（见附录）是作为标准拓扑优化代码构建的。主程序是通过该行从Matlab提示符开始的
top(nelx,nely,volfrac,penal,rmin)
式中，nelx和nely分别是水平和垂直方向上的单元数量，volfrac是体积分数，penal是惩罚功率，rmin是滤波器尺寸（除以单元尺寸）。其他变量和边界条件在Matlab代码中定义，如果需要可以编辑。对于拓扑优化循环中的每次迭代，代码都会生成当前的密度分布图。图1显示了通过输入下面这行调用附录中给出的代码获得的密度分布
top(60,20,0.5,3.0,1.5)
图1 MBB梁的拓扑优化。
顶部：全设计域，
中间：具有对称边界条件的半设计域，
底部：生成拓扑优化梁（两半）
默认边界条件对应于“MBB梁”的一半（图1）。荷载在左上角垂直施加，沿左边缘有对称边界条件，结构在右下角水平支撑。下面的小节将讨论Matlab代码的一些重要细节。
3.1 主程序（第1-36行）
主程序（第1-36行）开始于在设计域中均匀分布材料（第4行）。在其他一些初始化之后，主循环从调用有限元子程序（第12行）开始，该子程序返回位移矢量U。由于实体材料的单元刚度矩阵对所有单元都是相同的，因此只调用一次单元刚度矩阵子程序（第14行）。随后，对所有单元循环（第16-24行）确定目标函数和灵敏度（4）。变量n1和n2表示全局节点号中的左上和右上单元节点号，用于从全局位移向量U中提取单元位移向量Ue。灵敏度分析之后调用网格独立滤波器（第26行）和最优性标准优化器（第28行）。当前的柔度和其他参数由第30-33行求出，并绘制出密度分布图（第35行）。如果设计变量的变化（第30行中确定的变化）小于1%，则主循环终止。否则重复上述步骤。
3.2 基于优化标准的优化器（第37-48行）
更新的设计变量由优化器找到（第37-48行）。知道材料体积（sum（sum（xnew））是拉格朗日乘子（lag）的单调递减函数，满足体积约束的拉格朗日乘子的值可以通过双分段算法（第40-48行）找到。通过猜测拉格朗日乘子的下限l1和上限l2（第39行），初始化双分段算法。限制拉格朗日乘子的区间重复减半，直到其大小小于收敛标准（第40行）。
<h2
剩余内容已隐藏，支付完成后下载完整资料</h2
资料编号：[240373]，资料为PDF文档或Word文档，PDF文档可免费转换为Word
以上是毕业论文外文翻译，课题毕业论文、任务书、文献综述、开题报告、程序设计、图纸设计等资料可联系客服协助查找。
