这是用户在 2025-3-8 1:27 为 https://app.immersivetranslate.com/pdf-pro/86b2a6a4-3cf8-481b-93b9-613b4555af41/ 保存的双语快照页面，由 沉浸式翻译 提供双语支持。了解如何保存？


一种针对实际设计应用的桁架结构新型进化拓扑优化方法


赖亚平 a , b , ∗ a , b , ∗ ^(a,b,**){ }^{\mathrm{a}, \mathrm{b}, *} , 蔡琦 c c ^(c){ }^{\mathrm{c}} , 李宇 d d ^(d)^{\mathrm{d}} , 陈家勇 b b ^(b){ }^{\mathrm{b}} , 谢怡敏 e , ∗ ∗ e , ∗ ∗ ^(e,^(****)){ }^{\mathrm{e},{ }^{* *}}
重庆大学
土木工程学院，中国重庆 400045
b b ^(b){ }^{\mathrm{b}} T.Y. Lin 国际工程咨询（中国）有限公司，重庆 401121，中国
c c  ^("c "){ }^{\text {c }} 河海大学机械与材料学院，中国南京 211189
同济大学建筑与城市规划学院，中国上海 200092
e e ^(e){ }^{\mathrm{e}} 创新结构与材料中心，工程学院，RMIT 大学，澳大利亚墨尔本 3001


文章信息

  关键词：

BESO
  多种材料
  结构优化
  大小优化
  地面结构

  摘要


本文提出了一种新的拓扑优化方法和一个全面的工作流程，用于桁架结构的实际设计。首先讨论了桁架拓扑优化的实际设计要求，然后介绍了一种创建适合复杂几何形状的任意形状基础结构的技术。为了解决当前桁架优化方法的局限性，我们提出了一种双材料桁架双向进化结构优化（DMT-BESO）方法。该方法利用两种在拉伸和压缩允许应力及弹性模量上差异显著的材料。DMT-BESO 方法将最小能量原理与全应力设计标准相结合，使用杆件截面积作为设计变量，以实现拓扑和尺寸的同时优化。通过考虑应力约束，该方法确保符合行业标准，提高了安全性和材料利用率。 此外，提出了一种结构复杂性控制策略，以生成近似最优的桁架设计并简化优化设计，同时保持效率，使其更适合实际应用。通过数值示例和复合材料拱桥的设计验证了 DMT-BESO 方法及其复杂性控制策略的有效性。

  1. 引言


拓扑优化（TO）是一种重要的基于性能的设计技术，旨在优化规定设计域内材料的分布，以最大化最终结构的性能。近年来，这项技术在连续和离散设计领域都引起了显著关注，导致方法论的快速进展。基于连续体的拓扑优化的最流行和成熟的方法，广泛应用于各个学科，可以大致分为四类：（1）基于微穿孔复合材料的均匀化方法 [1]；（2）材料密度分布方法，以固体各向同性材料惩罚（SIMP）方法 [2-4] 为代表；（3）进化方法，包括进化结构优化（ESO）方法 [ 5 , 6 ] [ 5 , 6 ] [5,6][5,6] 和双向 ESO（BESO）方法 [7]；以及（4）边界变异方法，如基于水平集的方法 [8,9]。

和相场方法[10]。最近的拓扑优化（TO）趋势集中在提高计算效率[11]和整合多物理场[12]、多尺度[13]和多目标[14]考虑，以应对现实工程挑战。考虑到材料属性和加载条件的不确定性，稳健的拓扑优化 [ 15 , 16 ] [ 15 , 16 ] [15,16][15,16] 对于实际应用至关重要，增强了优化结构的可靠性和性能。算法的持续改进使得拓扑优化技术成为设计师在概念设计过程中探索轻量化、创新和高性能结构的强大工具[17-22]。

毫无疑问，使用基于连续体的拓扑优化技术的设计师可以在建筑和土木工程领域获得具有独特形状和优越结构性能的概念解决方案。然而，尽管对拓扑优化算法进行了广泛而深入的学术研究，但它们在实际项目中的应用仍然有限。连续体拓扑优化通常会生成难以使用的不规则形状。

后续设计阶段的建设成本高，尤其是对于大规模应用。

相比之下，桁架结构因其优良的刚度与重量比、模块化、材料效率和承载能力而被广泛应用于工程领域，主要是由于其简单性和易于施工。桁架由通过铰链连接的直杆构成，使其非常适合优化。然而，基于连续体的拓扑优化方法不适用于铰链连接的桁架结构。因此，专门针对桁架结构的拓扑优化算法的研究对于实际设计应用至关重要。

在 19 世纪末和 20 世纪初，Maxwell [26] 和 Cilly [27] 提出了桁架结构拓扑优化的基本思想。后来，Michell [28] 提出了在应力约束下具有最小体积的 Michell 桁架模型。在应力约束优化模型中，他证明了在单一载荷情况下，具有最小体积的结构中的杆件处于完全应力状态。然而，尽管 Michell 桁架被证明是在应力约束下的最优设计，但由于其包含无限数量的杆件，难以在实际工程中应用。直到 1964 年，Dorn 等人 [29] 提出了地面结构（GS）方法，桁架拓扑优化的数值实现才变得可行，并获得了实际工程设计的重要意义。从那时起，许多学者围绕 GS 方法进行了广泛的研究。 尽管这种方法是桁架拓扑优化中最广泛使用的，但仍然存在一些弱点，包括：（1）最优拓扑高度依赖于初始几何形状；（2）最优解通常包含过多的杆件，并且大多数杆件过于细小，无法直接使用；（3）用于理论研究的设计域通常仅限于矩形等规则几何形状，在这些形状中，几何形状节点可以以方形网格布局排列，而复杂不规则形状的设计域在实际工程设计中应予以考虑。为了解决这些问题，Martinez 等人[30]、Mckeown[31]、Rule[32]、He 等人[33]、Gilbert 和 Tyas[34]、Wang 和 Arora[35]、Lu 和 Xie[36]提出了一系列有效的解决方案。然而，桁架拓扑优化在实际工程应用中还需要考虑复杂的结构几何形状、不规则的几何形状生成和多材料优化，这在以往的研究中讨论得较少。

鉴于此，本研究首先关注在任意形状设计域内创建桁架 GS 以及实现拓扑优化的工作流程。此外，基于 Cai 和 Feng [37] 提出的桁架结构 BESO 算法，我们的新算法结合了两种材料的拉伸和压缩特性的差异。本文的其余部分组织如下：第 2 节概述了桁架拓扑优化中的实际设计要求、创建任意 GS 的方法以及拓扑优化实施平台。第 3 节总结了桁架结构 BESO 方法的理论背景，并提出了双材料桁架 BESO 算法。第 4 节和第 5 节分别展示了 2D 和 3D 数值示例。最后，第 6 节得出结论。


2. GS 创建方法，TO 实施平台和工作流程


2.1. 桁架中的实际设计要求 T O T O TOT O


尽管当前拓扑优化（TO）技术具有强大的能力，并展示了显著的应用潜力，但在实际应用中常常会遇到挑战，包括由于结果过于复杂而导致的实用性差，以及与现有工作流程的不连贯整合。桁架拓扑优化的实际设计要求包括已经解决的关键问题，如局部或全局不稳定性 [ 38 , 39 ] [ 38 , 39 ] [38,39][38,39] 。桁架拓扑优化的最新进展表明，混合算法有效地处理拓扑、尺寸和形状优化，同时提高结构性能和效率。Dang 等人[40]提出了一种单步混合方法，增强了复杂桁架设计的计算效率，而 Lieu[41]则提出了一个统一框架，以解决桁架的静态、动态和瞬态行为。类似地，

Oguz 等人[42]提出了一种针对增材制造的全球优化策略，强调在拓扑优化过程中整合设计约束的重要性。进一步的研究应关注可制造性、简易性、工作流程集成、材料效率和条件约束，以确保拓扑优化在实际应用中的可行性。遵循这些原则可以确保开发出适合实际使用的高效、高性能桁架结构。

(1) 可制造性：优化设计必须与现有的制造技术相适应，考虑诸如构件长度、最小截面积和特定制造工艺等限制。解决可制造性问题可以确保设计到施工的过渡更加顺畅，减少定制制造的需求。

(2) 简单性和优雅性：在追求最佳性能的同时，桁架 TO 的结果应保持简单和优雅。过于复杂的解决方案可能会使制造、组装和维护变得复杂，而更简单的设计通常更实用且更具成本效益。

(3) 材料效率：采用具有不同拉伸和压缩性能的材料可以增强结构性能和可持续性。利用多种材料可以实现优化设计，充分发挥每种材料的优势，从而实现高性能、经济有效和可持续的结构。

(4) 条件约束：解决多个约束，如应力、构件截面面积和使用性，确保桁架设计满足现实工程应用的多样化和复杂需求。这些约束有助于结构的稳健性和可靠性。

(5) 与现有工作流程的集成：TO 过程应与当前的设计和工程工作流程无缝集成，包括与 CAD 软件和其他工程工具的兼容性。有效的集成有助于在更广泛的项目工作流程中采用优化设计。


2.2. 创建任意形态学 GS 的方法


所谓的 GS 方法涉及通过根据设计域的形状和边界条件连接离散节点来创建桁架网络。然而，使用传统的 MATLAB 代码生成 GS [43,44]，尽管系统化，但在处理复杂几何形状和实际应用于现实场景方面往往面临灵活性挑战。为了解决这个问题，我们开发了一种使用 Rhino 和 Grasshopper 创建 GS 的参数化方法。

Rhino 软件广泛应用于多个行业，擅长管理复杂和不规则的几何形状，提供比基于 MATLAB 的方法更大的灵活性。Rhino 的高级可视化工具使设计师能够在 3D 中可视化复杂的几何形状，帮助理解空间关系和潜在的设计问题。设计师可以轻松创建复杂的形状，并调整 GS 以适应非标准设计领域。Grasshopper 插件提供强大的参数化设计能力，使用户能够在可视化编程环境中快速修改和迭代设计。这种迭代过程和灵活性促进了高效的设计优化和各种设计替代方案的探索。此外，Rhino 与其他 CAD 软件无缝集成，增强了 GS 在更广泛设计过程中的实用性和适用性。

一种简单高效的方法，如图 1(a)所示，涉及使用“Grasshopper Python Script Editor”来创建 GS。用户在 Python 代码中输入相关参数，整个 GS 的 Abaqus 输入文件（.inp）信息会自动生成。例如，在创建规则平面网格 GS 时，如图 1(b)所示，过程首先从定义关键参数开始，例如网格尺寸、网格间距和连接级别。接下来，脚本生成一个点网格。

图 1. 使用“Grasshopper Python Script Editor”生成 GS 的示意图：(a) GS 的创建和 Abaqus 输入文件（.inp）的生成；(b) 一个常规平面 GS 的示例；(c) 其对应空间表面的基础网格；(d) 从常规平面 GS 的演变中得出的空间 GS。

在平面上，表示 GS 的潜在节点位置。最后，它根据指定的连接级别系统地在这些点之间创建连接（或杆），确保通过检查共线性并去除冗余杆，仅保留最短和最有效的连接。这种方法和原则与文献中建立的原则一致[43]。此外，用户可以将节点坐标投影到从规则平面网格到任意形状的空间表面，以生成空间形式的 GS，如图 1©和(d)所示。这种方法利用 Python 的脚本功能来自动化设计过程，生成一个清晰的 GS，可以在 Rhino 中轻松可视化和细化。

为了获得复杂和任意的 GS，另一种方法涉及使用 Rhino 和 Grasshopper 进行参数建模，通过“参数对象”和“关联逻辑”定义 3D 几何体。这种方法基于多个空间控制线创建具有任意形态的 GS，使得 GS 的 3D 几何体能够根据特定约束和参数变量进行调整，并自动生成 Abaqus 预处理输入数据。例如，在创建如图 2(a)和(b)所示的拱形 3D 空间 GS 时，首先生成顶部和底部的初始空间控制线，如图 2(c)所示；然后设置控制参数，例如沿控制线的节点数量和连接（或杆）之间的最大角度，如图 2(d)所示。这个基于参数的模型可以根据每个参数的值提供多种解决方案。

通过将建模过程中生成的参数文件定义为可重用的功能文件，可以显著提高建模 GS 的效率。

无论采用何种方法，初始的 GS 设计数据是在 Grasshopper 中生成的。有限元分析（FEA）的 GS 输入文件是严格按照 Abaqus 输入文件（.inp）的结构和格式创建和导出的。


2.3. TO 实施平台和工作流程


Abaqus 是一款多功能的通用有限元分析（FEA）软件包，因其在执行复杂模拟和分析方面的强大能力而被广泛使用。Abaqus 脚本接口由 Python 驱动，提供对软件高级有限元分析功能的全面访问，允许对建模、分析和后处理任务进行自动化。这种集成便于将自定义算法和工作流程无缝地纳入 Abaqus 环境，特别有利于迭代优化（TO）。我们的优化过程的实施遵循一个系统的工作流程，包括以下关键步骤：

(1) 导入初始设计：从 Grasshopper 导出的.inp 文件生成的初始设计作为优化过程的起点。

(a)


©

(d)

图 2. 使用“参数对象”和“关联逻辑”生成 GS 的示意图：(a) Grasshopper 参数和 Abaqus 输入文件；(b) 一个拱形的 3D 空间 GS；(c) 初始控制线；(d) 主要控制参数。

(2) 执行有限元分析：使用 Abaqus，初始 GS 进行有限元分析，以评估其在特定载荷情况下和边界条件下的结构性能。包括节点位移、单元应力和应变能的分析结果是优化算法的重要输入。

(3) 执行优化算法：TO 算法通过迭代修改结构，以最小化目标函数，同时满足设计约束。在每次迭代中，修订后的设计在 Abaqus 中进行分析，结果指导后续的设计演变。

(4) 更新设计并判断收敛：设计变量根据优化标准进行更新。这个迭代过程持续进行，直到满足收敛标准。

（5）导出最终设计：在收敛后，优化的桁架结构从 Abaqus 导出，准备进行进一步分析、验证或其他应用。

图 3 显示了用于拓扑优化实施平台的软件和插件。通过利用 Abaqus、Grasshopper 和 Python 的能力，这一拓扑优化过程有效地整合了现实世界的设计考虑，以满足现代工程应用的复杂需求。

图 3. 使用的软件程序和插件。

  3. 方法论


3.1. BESO 的理论背景


作为拓扑优化的重要分支，连续 BESO 方法允许材料的同时去除和添加。以最小结构应变能为目标函数，以单元密度为设计变量，以结构优化体积比和平衡方程为约束的数学模型可以分别表示如下。
  目标函数：
C ( x i ) = 1 2 P T U C x i = 1 2 P T U C(x_(i))=(1)/(2)P^(T)UC\left(x_{i}\right)=\frac{1}{2} \mathbf{P}^{\mathrm{T}} \mathbf{U}
  设计变量：
x = { x 1 , x 2 , … , x n } T , x ∈ Ω R x = x 1 , x 2 , … , x n T , x ∈ Ω R x={x_(1),x_(2),dots,x_(n)}^(T),x inOmega^(R)\boldsymbol{x}=\left\{x_{1}, x_{2}, \ldots, x_{\mathrm{n}}\right\}^{\mathrm{T}}, \boldsymbol{x} \in \Omega^{R}
  约束条件：
{ K T U = P ∑ i = 1 n V i x i − γ V 0 = 0 E i = E i 0 x i p x i = x min or 1 K T U = P ∑ i = 1 n   V i x i − γ V 0 = 0 E i = E i 0 x i p x i = x min  or  1 {[K^(T)U=P],[sum_(i=1)^(n)V_(i)x_(i)-gammaV_(0)=0],[E_(i)=E_(i0)x_(i)^(p)],[x_(i)=x_(min)" or "1]:}\left\{\begin{array}{l}\mathbf{K}^{\mathrm{T}} \mathbf{U}=\mathbf{P} \\ \sum_{i=1}^{n} V_{i} x_{i}-\gamma V_{0}=0 \\ E_{i}=E_{i 0} x_{i}^{p} \\ x_{i}=x_{\min } \text { or } 1\end{array}\right.

其中 i i ii 表示元素编号， i = 1 , 2 , … , n i = 1 , 2 , … , n i=1,2,dots,ni=1,2, \ldots, n 和 x i x i x_(i)x_{i} 是二进制设计变量，只能取值为 x min x min x_(min)x_{\min } 或 1。这里， x min x min x_(min)x_{\min } 是一个大于零的较小数，通常为 0.001。如果 x i = x min x i = x min x_(i)=x_(min)x_{i}=x_{\min } ，则该元素被视为已删除；如果 x i = 1 x i = 1 x_(i)=1x_{i}=1 ，则该元素保留在优化结构中。 C ( x i ) C x i C(x_(i))C\left(x_{i}\right) 表示结构的应变能，作为设计变量 x i . P , U x i . P , U x_(i).P,Ux_{i} . \mathbf{P}, \mathbf{U} 的函数， K K K\mathbf{K} 分别是外部载荷向量、位移向量和总刚度矩阵。 V i V i V_(i)V_{i} 是第 i i ii 个元素的体积， V 0 V 0 V_(0)V_{0} 是优化前的原始体积。 γ γ gamma\gamma 是目标体积优化系数，表示优化体积与原始体积的比率。 E i 0 E i 0 E_(i0)E_{i 0} 和 E i E i E_(i)E_{i} 分别是第 i i ii 个元素优化前后的弹性模量， p p pp 是元素密度的惩罚系数，通常根据文献 [45] 取为 3。

BESO 是一种基于梯度的方法，通过对目标函数进行微分来计算元素灵敏度。每个元素的灵敏度表示由于元素状态的变化而导致的目标函数的变化，通常在 x min x min x_(min)x_{\min } 到 1 之间变化。优化过程从一个包含所有元素的初始完整设计域开始。接下来是合规性最小化的目标。

在此过程中，灵敏度最高的元素从 1 过渡到 x min x min x_(min)x_{\min } ，而灵敏度最低的元素则从 x min x min x_(min)x_{\min } 过渡到 1。这一过程是寻找最优解的典型“最陡下降”迭代方法。关于连续 BESO 方法的更多细节可以在文献[46,47]中找到。

受到连续 BESO 方法优化策略的启发，蔡和冯提出了一种桁架-BESO（T-BESO）方法。该方法在 MATLAB 平台上实现，能够使用单一材料实现二维平面桁架优化。数值示例表明，该方法的计算结果与解析解的偏差小于 1 % 1 % 1%1 \% ，证明了其准确性和有效性。有关 T-BESO 方法的更多细节，请参见参考文献[37]。


3.2. 双材料 T-BESO


考虑到在桁架 TO 中多种材料在拉伸和压缩性能上的差异，对于优化材料使用和提高成本效益至关重要。根据力的特性分配材料不仅可以增强结构性能，还能提供创新和适应性强的解决方案。因此，我们在 TBESO 算法框架的基础上进一步开发了一种双材料 T-BESO (DMT-BESO) 方法，允许分别对拉伸和压缩杆使用具有显著差异的允许应力和弹性模量的两种材料。

这种方法在实际工程中尤其有价值，因为像钢和混凝土这样的材料具有本质上不同的拉伸和压缩强度，通常会一起使用以实现更高效和更具韧性的结构。通过结合具有不同特性的材料，桁架 TO 可以最小化材料消耗，降低施工成本，并提高最终设计的耐久性，使其在土木工程和基础设施开发等行业中具有很高的适用性。


3.2.1. 前提假设


在本研究中，我们简化了 DMT-BESO 过程，以关注整体结构性能和材料分布。关于界面区域的挑战在以下假设下得到解决：

(1) 线性弹性材料行为：假设这两种材料在整个优化过程中都表现出线性弹性行为。这一简化使我们能够专注于结构性能，避免与非线性相关的复杂性，例如塑性变形或材料在界面处的失效。

(2) 最小结构变形：假设结构在外力作用下发生最小变形。这个小变形假设简化了分析，确保应力和位移保持在线性弹性范围内。

(3) 界面的连续性：一个关键假设是材料在连接界面保持完全的连续性。位移和应力在界面上是连续的，确保完全耦合的机械行为。因此，界面上不会发生滑移、分离或应力集中，连接的行为就像这两种材料是完美结合的一样。

（4）局部节点连接：不同材料之间的连接被视为节点的局部问题。本研究不考虑特定的界面行为，如粘附、结合或局部应力集中。相反，连接被简化，以确保材料在全局优化过程中的粘合行为。处理界面行为的详细力学超出了本研究的范围，但可以在未来的研究中进行探索。


3.2.2. 优化公式


该优化问题的数学模型定义如下：

(1) 目标函数：在实际应用中，最小化应变能确保结构在施加荷载下实现最佳刚度，这对机械效率和可服务性至关重要。最小化应变能的结构本质上对变形更具抵抗力，并且在循环加载条件下往往表现出更优越的性能。DMT-BESO 方法通过将能量原理与全应力设计标准相结合，使用应变能最小化作为目标函数。当满足全应力条件时，桁架中的每个构件都被充分利用到其容量，从而有效地导致最佳材料分布。
C ( A i ) = 1 2 F T U C A i = 1 2 F T U C(A_(i))=(1)/(2)F^(T)UC\left(A_{i}\right)=\frac{1}{2} \mathbf{F}^{\mathrm{T}} \mathbf{U}

(2) 设计变量：与连续 BESO 方法不同，在该方法中，设计变量 x i x i x_(i)x_{i} 只能取值 x min x min x_(min)x_{\min } 或 1，DMTBESO 方法使用杆件的截面积作为设计变量。利用连续设计变量在优化过程中提供了更大的灵活性，以探索设计空间。虽然在实际应用中，使用从标准化目录值中选择的离散截面积是必要的，但过早施加此类限制可能会限制识别最佳拓扑的能力，并完全探索设计空间。因此，建议在拓扑优化阶段应用连续设计变量，以实现最佳或接近最佳的解决方案。然后，可以在后处理阶段执行聚类算法，以满足制造标准，如文献[36]中所建议的。
A = { A 1 , A 2 , … , A n } T , A ∈ Ω R A = A 1 , A 2 , … , A n T , A ∈ Ω R A={A_(1),A_(2),dots,A_(n)}^(T),A inOmega^(R)\boldsymbol{A}=\left\{A_{1}, A_{2}, \ldots, A_{\mathrm{n}}\right\}^{\mathrm{T}}, \boldsymbol{A} \in \Omega^{R}

DMT-BESO 方法定义了两种材料：材料 1 和材料 2，其属性列在表 1 中。
  表 1

材料的相关属性。
  材料
弹性模量 E E EE 		  允许的应力 [ σ ] [ σ ] [sigma][\sigma]
  材料 1 E t E t E_(t)E_{t}    σ 0 t σ 0 t sigma_(0t)\sigma_{0 t} （在张力中）
  材料 2 E c E c E_(c)E_{c}    σ 0 c σ 0 c sigma_(0c)\sigma_{0 c} （在压缩中）
Materials Modulus of elasticity E Allowable stress [sigma] material 1 E_(t) sigma_(0t) (in tension) material 2 E_(c) sigma_(0c) (in compression)| Materials | Modulus of elasticity $E$ | Allowable stress $[\sigma]$ | | :--- | :--- | :--- | | material 1 | $E_{t}$ | $\sigma_{0 t}$ (in tension) | | material 2 | $E_{c}$ | $\sigma_{0 c}$ (in compression) |

两种材料的泊松比为 μ μ mu\mu ，弹性模量分别为 E t E t E_(t)E_{t} 和 E c E c E_(c)E_{c} 。材料 1 更适合拉伸力，而材料 2 更适合压缩力。材料 1 在拉伸时的允许应力为 σ 0 t σ 0 t sigma_(0t)\sigma_{0 t} ，材料 2 在压缩时的允许应力为 σ 0 c σ 0 c sigma_(0c)\sigma_{0 c} 。根据公式(5)，引入一个新变量 m t i m t i mt_(i)m t_{i} 来描述材料的分布和性质：
m t i = { 1 ( material1 ) 2 ( material ) , i = 1 , 2 , … , n m t i = 1 (  material1  ) 2 (  material  ) , i = 1 , 2 , … , n mt_(i)={[1(" material1 ")],[2(" material ")],i=1,2,dots,n:}m t_{i}=\left\{\begin{array}{l}1(\text { material1 }) \\ 2(\text { material })\end{array}, \mathrm{i}=1,2, \ldots, \mathrm{n}\right.

(3) 约束：在工程设计中，强度是一个关键标准。应力约束确保优化设计符合行业标准和规范，保持安全性和可靠性。在桁架优化朝向实际设计应用的主要目标是确保所有杆件能够承受施加的荷载而不发生屈服。这需要包括应力约束，以保证每个杆件在安全应力限制内运行，从而防止结构失效。因此，我们纳入了两个材料的相关面积和应力约束，如公式（7）所定义：
{ K T U = F A t , min ⩽ A t ⩽ A t , max A c , min ⩽ A c ⩽ A c , max σ 0 , t − ⩽ σ t ⩽ σ 0 , t + σ 0 , c − ⩽ σ c ⩽ σ 0 , c + K T U = F A t ,  min  ⩽ A t ⩽ A t ,  max  A c ,  min  ⩽ A c ⩽ A c ,  max  σ 0 , t − ⩽ σ t ⩽ σ 0 , t + σ 0 , c − ⩽ σ c ⩽ σ 0 , c + {[K^(T)U=F],[A_(t," min ") <= A_(t) <= A_(t," max ")],[A_(c," min ") <= A_(c) <= A_(c," max ")],[sigma_(0,t)^(-) <= sigma_(t) <= sigma_(0,t)^(+)],[sigma_(0,c)^(-) <= sigma_(c) <= sigma_(0,c)^(+)]:}\left\{\begin{array}{l}\mathbf{K}^{\mathrm{T}} \mathbf{U}=\mathbf{F} \\ A_{t, \text { min }} \leqslant A_{t} \leqslant A_{t, \text { max }} \\ A_{c, \text { min }} \leqslant A_{c} \leqslant A_{c, \text { max }} \\ \sigma_{0, t}^{-} \leqslant \sigma_{t} \leqslant \sigma_{0, t}^{+} \\ \sigma_{0, c}^{-} \leqslant \sigma_{c} \leqslant \sigma_{0, c}^{+}\end{array}\right.

在公式（7）中，首先必须确保设计满足平衡方程。在识别拉伸和压缩杆并分配材料后，第 i i ii 根杆的截面积为 A t A t A_(t)A_{t} 或 A c A c A_(c)A_{c} 。拉伸杆的最小和最大截面积分别为 A t , min A t ,  min  A_(t," min ")A_{t, \text { min }} 和 A t , max A t ,  max  A_(t," max ")A_{t, \text { max }} ，而压缩杆的最小和最大截面积为 A c , min A c ,  min  A_(c," min ")A_{c, \text { min }} 和 A c , max A c ,  max  A_(c," max ")A_{c, \text { max }} 。为了优化桁架结构的尺寸， A t A t A_(t)A_{t} 和 A c A c A_(c)A_{c} 在范围 [ A t , min , A t , max A t ,  min  , A t ,  max  [A_(t," min "),A_(t," max "):}\left[A_{t, \text { min }}, A_{t, \text { max }}\right. ]和[ A c , min , A c , max A c , min , A c , max A_(c,min),A_(c,max)A_{c, \min }, A_{c, \max } ]内取连续值。此外， A t , min A t , min A_(t,min)A_{t, \min } 和 A c , min A c , min A_(c,min)A_{c, \min } 被设置为大于零的非常小的值，以防止结构刚度矩阵的奇异性。第 i i ii 根杆的拉伸和压缩应力分别用 σ t σ t sigma_(t)\sigma_{t} 和 σ c σ c sigma_(c)\sigma_{c} 表示。材料在拉伸下的允许拉伸和压缩应力为 σ 0 , t − σ 0 , t − sigma_(0,t)^(-)\sigma_{0, t}^{-} 和 σ 0 , t + σ 0 , t + sigma_(0,t)^(+)\sigma_{0, t}^{+} ，而在压缩下的材料允许拉伸和压缩应力为 σ 0 , c − σ 0 , c − sigma_(0,c)^(-)\sigma_{0, c}^{-} 和 σ 0 , c + σ 0 , c + sigma_(0,c)^(+)\sigma_{0, c}^{+} 。


3.2.3. 材料分配


在桁架结构的双材料分配策略中，具有优越拉伸性能的材料被分配给拉伸杆，而具有优越压缩性能的材料则被分配给压缩杆。这种直观的方法与多材料结构的实际工程设计原则相一致。因此，具有显著拉伸和压缩性能差异的材料可以根据其关键机械性能进行分配，例如弹性模量和允许应力：
E i = { E t ( σ i ≥ 0 ) E c ( σ i < 0 ) E i = E t σ i ≥ 0 E c σ i < 0 E_(i)={[E_(t),(sigma_(i) >= 0)],[E_(c),(sigma_(i) < 0)]:}E_{i}= \begin{cases}E_{t} & \left(\sigma_{i} \geq 0\right) \\ E_{c} & \left(\sigma_{i}<0\right)\end{cases}
σ 0 , i = { σ 0 , t ( σ i ≥ 0 ) σ 0 , c ( σ i < 0 ) σ 0 , i = σ 0 , t σ i ≥ 0 σ 0 , c σ i < 0 sigma_(0,i)={[sigma_(0,t),(sigma_(i) >= 0)],[sigma_(0,c),(sigma_(i) < 0)]:}\sigma_{0, i}= \begin{cases}\sigma_{0, t} & \left(\sigma_{i} \geq 0\right) \\ \sigma_{0, c} & \left(\sigma_{i}<0\right)\end{cases}

其中 E i , σ 0 , i E i , σ 0 , i E_(i),sigma_(0,i)E_{i}, \sigma_{0, i} 和 σ i σ i sigma_(i)\sigma_{i} 分别表示弹性模量、允许应力，以及每次优化迭代步骤中第 i i ii 根杆的正常应力。

根据拉伸和压缩应力状态分配材料发生在优化的早期阶段，因此，两种材料之间的机械性能显著差异可能导致初始优化过程中的剧烈和不连续变化。为此，我们采用以下策略，如下所示。

图 4. 迭代过程中材料机械性能变化的策略：(a) 弹性模量；(b) 允许应力。

图 4：最初，两个材料的初始弹性模量和允许应力被设置为它们的中位值。这确保了在优化开始时两个材料的机械性能是相同的。然后，随着优化过程的发展，在优化结构的拉伸和压缩材料分配稳定后，逐渐过渡到显著不同的性能。

对应于图 4，材料的弹性模量和允许应力随迭代过程变化的数学表达式如下：
E t = { ( E c + E t ) 2 + ( E t − E c ) 2 ( 1 1 + e − ρ ( m − k 2 ) ) ( m ≤ k ) E t = E c + E t 2 + E t − E c 2 1 1 + e − ρ m − k 2 ( m ≤ k ) E_(t)={((E_(c)+E_(t)))/(2)+((E_(t)-E_(c)))/(2)((1)/(1+e^(-rho(m-(k)/(2)))))quad(m <= k):}E_{t}=\left\{\frac{\left(E_{c}+E_{t}\right)}{2}+\frac{\left(E_{t}-E_{c}\right)}{2}\left(\frac{1}{1+e^{-\rho\left(m-\frac{k}{2}\right)}}\right) \quad(m \leq k)\right.
E t E t E_(t)E_{t}
( m > k ) ( m > k ) (m > k)(m>k)
E c = { ( E c + E t ) 2 + ( E c − E t ) 2 ( 1 1 + e − ρ ( m − k 2 ) ) ( m ≤ k ) E c = E c + E t 2 + E c − E t 2 1 1 + e − ρ m − k 2 ( m ≤ k ) E_(c)={((E_(c)+E_(t)))/(2)+((E_(c)-E_(t)))/(2)((1)/(1+e^(-rho(m-(k)/(2)))))quad(m <= k):}E_{c}=\left\{\frac{\left(E_{c}+E_{t}\right)}{2}+\frac{\left(E_{c}-E_{t}\right)}{2}\left(\frac{1}{1+e^{-\rho\left(m-\frac{k}{2}\right)}}\right) \quad(m \leq k)\right.
E c E c E_(c)E_{c}
( m > k ) ( m > k ) (m > k)(m>k)
σ 0 , t = { ( σ 0 , c + σ 0 , t ) 2 + ( σ 0 , t − σ 0 , c ) 2 ( 1 1 + e − ρ ( m − k 2 ) ) σ 0 , t ( m ≤ k ) ( m > k ) σ 0 , t = σ 0 , c + σ 0 , t 2 + σ 0 , t − σ 0 , c 2 1 1 + e − ρ m − k 2 σ 0 , t ( m ≤ k ) ( m > k ) {:[sigma_(0,t)={[((sigma_(0,c)+sigma_(0,t)))/(2)+((sigma_(0,t)-sigma_(0,c)))/(2)((1)/(1+e^(-rho(m-(k)/(2)))))],[sigma_(0,t)](m <= k):}],[(m > k)]:}\begin{aligned} \sigma_{0, t}=\left\{\begin{array}{l} \frac{\left(\sigma_{0, c}+\sigma_{0, t}\right)}{2}+\frac{\left(\sigma_{0, t}-\sigma_{0, c}\right)}{2}\left(\frac{1}{1+e^{-\rho\left(m-\frac{k}{2}\right)}}\right) \\ \sigma_{0, t} \end{array}(m \leq k)\right. \\ (m>k) \end{aligned}
σ 0 , c = { ( σ 0 , c + σ 0 , t ) 2 + ( σ 0 , t − σ 0 , c ) 2 ( 1 1 + e − ρ ( m − k 2 ) ) σ 0 , c ( m ≤ k ) ( m > k ) σ 0 , c = σ 0 , c + σ 0 , t 2 + σ 0 , t − σ 0 , c 2 1 1 + e − ρ m − k 2 σ 0 , c ( m ≤ k ) ( m > k ) {:[sigma_(0,c)={[((sigma_(0,c)+sigma_(0,t)))/(2)+((sigma_(0,t)-sigma_(0,c)))/(2)((1)/(1+e^(-rho(m-(k)/(2))))):}],[sigma_(0,c)](m <= k):}],[(m > k)]:}\begin{aligned} \sigma_{0, c}=\left\{\begin{array}{l} \frac{\left(\sigma_{0, c}+\sigma_{0, t}\right)}{2}+\frac{\left(\sigma_{0, t}-\sigma_{0, c}\right)}{2}\left(\frac{1}{\left.1+e^{-\rho\left(m-\frac{k}{2}\right)}\right)}\right. \\ \sigma_{0, c} \end{array}(m \leq k)\right. \\ (m>k) \end{aligned}

其中 m m mm 表示迭代步骤的数量， k k kk 表示在优化过程中完全完成材料属性变化所需的总步骤数。根据我们的测试，当 k k kk 是不小于 400 的整数且 ρ ρ rho\rho 是 0.025 时，可以实现平滑过渡。


3.2.4. 敏感性分析


进行灵敏度分析是评估在桁架元素的截面面积变化时应变能如何变化。这涉及计算应变能相对于每个设计变量的导数，从而允许进行有根据的更新。

设目标函数为 Eq. (4) 为 C ( A ) = 1 2 F T U = W e x C ( A ) = 1 2 F T U = W e x C(A)=(1)/(2)F^(T)U=W_(ex)C(A)=\frac{1}{2} \mathbf{F}^{\mathrm{T}} \mathbf{U}=W_{e x} ，其中 W e x W e x W_(ex)W_{e x} 表示外部载荷所做的功。由于杆件的轴向变形，桁架结构的内部功为：
W i n = ∑ i = 1 n 1 2 F i Δ i = ∑ i = 1 n F i 2 L i 2 E i A i W i n = ∑ i = 1 n   1 2 F i Δ i = ∑ i = 1 n   F i 2 L i 2 E i A i W_(in)=sum_(i=1)^(n)(1)/(2)F_(i)Delta_(i)=sum_(i=1)^(n)(F_(i)^(2)L_(i))/(2E_(i)A_(i))W_{i n}=\sum_{i=1}^{n} \frac{1}{2} F_{i} \Delta_{i}=\sum_{i=1}^{n} \frac{F_{i}^{2} L_{i}}{2 E_{i} A_{i}}

其中 F i , Δ i F i , Δ i F_(i),Delta_(i)F_{i}, \Delta_{i} 和 L i L i L_(i)L_{i} 分别表示第 i i ii 根杆的轴向力、伸长和长度， E i E i E_(i)E_{i} 表示其材料的弹性模量。

根据能量原理，外部功 W e x W e x W_(ex)W_{e x} 等于内部功 W i n W i n W_(in)W_{i n} 。因此，目标函数可以修订为：
C ( A i ) = ∑ i = 1 n F i 2 L i 2 E i A i C A i = ∑ i = 1 n   F i 2 L i 2 E i A i C(A_(i))=sum_(i=1)^(n)(F_(i)^(2)L_(i))/(2E_(i)A_(i))C\left(A_{i}\right)=\sum_{i=1}^{n} \frac{F_{i}^{2} L_{i}}{2 E_{i} A_{i}}

为了建立目标函数与设计变量之间的关系，方程（15）是关于设计变量 A i A i A_(i)A_{i} 推导的，得到：
α i = ∂ C ( A i ) ∂ A i = − F i 2 L i 2 E i A i 2 = − σ i 2 L i 2 E i α i = ∂ C A i ∂ A i = − F i 2 L i 2 E i A i 2 = − σ i 2 L i 2 E i alpha_(i)=(del C(A_(i)))/(delA_(i))=-(F_(i)^(2)L_(i))/(2E_(i)A_(i)^(2))=-(sigma_(i)^(2)L_(i))/(2E_(i))\alpha_{i}=\frac{\partial C\left(A_{i}\right)}{\partial A_{i}}=-\frac{F_{i}^{2} L_{i}}{2 E_{i} A_{i}^{2}}=-\frac{\sigma_{i}^{2} L_{i}}{2 E_{i}}

其中 α i α i alpha_(i)\alpha_{i} 是元素 i i ii 的灵敏度，反映了截面 A i A i A_(i)A_{i} 的变化对结构整体应变能的影响。由于 GS 中大多数杆件长度不一致，为消除杆件长度对灵敏度的影响，可以将公式 (16) 修订为：
α ~ i = − α i L i = σ i 2 2 E i α ~ i = − α i L i = σ i 2 2 E i widetilde(alpha)_(i)=-(alpha_(i))/(L_(i))=(sigma_(i)^(2))/(2E_(i))\widetilde{\alpha}_{i}=-\frac{\alpha_{i}}{L_{i}}=\frac{\sigma_{i}^{2}}{2 E_{i}}

结合方程（17）和（8），杆的元素灵敏度：
α ~ i = { σ i 2 2 E t ( σ i ≥ 0 ) σ i 2 2 E c ( σ i < 0 ) α ~ i = σ i 2 2 E t σ i ≥ 0 σ i 2 2 E c σ i < 0 widetilde(alpha)_(i)={[(sigma_(i)^(2))/(2E_(t)),(sigma_(i) >= 0)],[(sigma_(i)^(2))/(2E_(c)),(sigma_(i) < 0)]:}\widetilde{\alpha}_{i}= \begin{cases}\frac{\sigma_{i}^{2}}{2 E_{t}} & \left(\sigma_{i} \geq 0\right) \\ \frac{\sigma_{i}^{2}}{2 E_{c}} & \left(\sigma_{i}<0\right)\end{cases}

最佳桁架设计旨在实现完全应力状态，在该状态下，所有杆件在施加荷载下同时达到其允许的应力极限。这种方法最大化了材料效率和结构性能。

性能。Ahrari 和 Atai [48] 以及 He 等人 [49] 的研究表明，最佳桁架设计实现了完全应力条件。因此，公式 (18) 可以修订如下：
α ~ 0 = { σ 0 , t 2 2 E t ( σ i ≥ 0 ) σ 0 , c 2 2 E c ( σ i < 0 ) α ~ 0 = σ 0 , t 2 2 E t σ i ≥ 0 σ 0 , c 2 2 E c σ i < 0 widetilde(alpha)_(0)={[(sigma_(0,t)^(2))/(2E_(t)),(sigma_(i) >= 0)],[(sigma_(0,c)^(2))/(2E_(c)),(sigma_(i) < 0)]:}\widetilde{\alpha}_{0}= \begin{cases}\frac{\sigma_{0, t}^{2}}{2 E_{t}} & \left(\sigma_{i} \geq 0\right) \\ \frac{\sigma_{0, c}^{2}}{2 E_{c}} & \left(\sigma_{i}<0\right)\end{cases}

其中 α ~ 0 α ~ 0 widetilde(alpha)_(0)\widetilde{\alpha}_{0} 是灵敏度阈值。

DMT-BESO 算法基于全应力标准更新杆件的截面积，动态调整材料分布以优化结构性能。当元素灵敏度 α ~ i α ~ i widetilde(alpha)_(i)\widetilde{\alpha}_{i} 小于 α ~ 0 α ~ 0 widetilde(alpha)_(0)\widetilde{\alpha}_{0} 时，通过减少杆件的截面积可以提高材料利用效率。相反，当 α ~ i α ~ i widetilde(alpha)_(i)\widetilde{\alpha}_{i} 超过 α ~ 0 α ~ 0 widetilde(alpha)_(0)\widetilde{\alpha}_{0} 时，需要增加截面积以确保结构安全。该策略通过调整杆件的截面积有效优化了材料利用。若截面积达到最小值 A t , min A t , min A_(t,min)A_{t, \min } 或 A c , min A c , min A_(c,min)A_{c, \min } ，则由于刚度不足，杆件被视为从结构中移除，从而实现结构的拓扑变化。

值得注意的是，连续 BESO 方法在根据灵敏度阈值添加或删除元素之前需要对元素灵敏度进行排序。然而，DMT-BESO 方法通过采用不同的进化优化策略绕过了这一排序要求。


3.2.5. 设计变量更新策略


DMT-BESO 方法遵循典型的进化优化程序，需要对设计变量进行迭代更新以达到最优解。虽然足够的迭代步骤可以确保平滑和稳定的优化过程，但过多的步骤会降低计算效率。为了解决这个问题，我们提出以下设计变量更新策略：在优化过程开始时，较大的步长更新有利于快速找到可行解。随着过程的进行并接近收敛，较小的步长更新有助于避免振荡，实现稳定和准确的收敛。为了通过最小化迭代步骤的数量来提高计算效率，同时确保平滑性和稳定性，本研究引入了一种动态自适应更新机制，以控制设计变量的更新速率，如下所示：
ω i = { ( 1.0 + 0.1 e − μ m ) ( α ~ i ≥ α ~ 0 ) ( 1.0 − 0.1 e − μ m ) ( α ~ i < α ~ 0 ) ω i = 1.0 + 0.1 e − μ m α ~ i ≥ α ~ 0 1.0 − 0.1 e − μ m α ~ i < α ~ 0 omega_(i)={[(1.0+0.1e^(-mu m)),( widetilde(alpha)_(i) >= widetilde(alpha)_(0))],[(1.0-0.1e^(-mu m)),( widetilde(alpha)_(i) < widetilde(alpha)_(0))]:}\omega_{i}= \begin{cases}\left(1.0+0.1 \mathrm{e}^{-\mu m}\right) & \left(\widetilde{\alpha}_{i} \geq \widetilde{\alpha}_{0}\right) \\ \left(1.0-0.1 \mathrm{e}^{-\mu m}\right) & \left(\widetilde{\alpha}_{i}<\widetilde{\alpha}_{0}\right)\end{cases}
μ = { 0.005 e − 15 λ ( λ < 0.1 ) 0.001 ( λ ≥ 0.1 ) μ = 0.005 e − 15 λ ( λ < 0.1 ) 0.001 ( λ ≥ 0.1 ) mu={[0.005e^(-15 lambda)(lambda < 0.1)],[0.001(lambda >= 0.1)]:}\mu=\left\{\begin{array}{l}0.005 e^{-15 \lambda}(\lambda<0.1) \\ 0.001(\lambda \geq 0.1)\end{array}\right.

(a)

在公式 (20) 中， ω i ω i omega_(i)\omega_{i} 代表第 i i ii 个条形的更新因子。指数衰减项 e − μ m e − μ m e^(-mum)\mathrm{e}^{-\mu \mathrm{m}} 确保 ω i ω i omega_(i)\omega_{i} 随着迭代步骤 m m mm 的增加而减少，从而在解接近最优状态时允许进行更细致的调整。在公式 (21) 中， λ λ lambda\lambda 表示收敛度量。通过区分 λ λ lambda\lambda 小于 0.1 和不小于 0.1 的情况，该函数应用不同的衰减速率。这种区分使得在解远离收敛时进行更激进的更新，而在接近时则进行更谨慎的更新。

图 5 展示了由方程(20)和(21)表示的自适应更新策略，显示了更新因子在迭代步骤中的变化。使用不同的指数速率（蓝线和红线）为更新策略提供了灵活性。图 5(a)和(b)都强调了一个自适应范围，在该范围内更新因子动态调整，表示为阴影区域。在每个迭代步骤中，更新因子可以根据收敛指标 λ λ lambda\lambda 在蓝色和红色范围内自适应选择。这个自适应机制对于实现快速收敛与精确微调之间的平衡至关重要，最终导致最佳和稳健的设计。

对于增加设计变量：
A i , m = { max ( A t , max , ω i A i , m − 1 ) ( σ i ≥ 0 ) max ( A c , max , ω i A i , m − 1 ) ( σ i < 0 ) A i , m = max A t , max , ω i A i , m − 1 σ i ≥ 0 max A c ,  max  , ω i A i , m − 1 σ i < 0 A_(i,m)={[max(A_(t,max),omega_(i)A_(i,m-1)),(sigma_(i) >= 0)],[max(A_(c," max "),omega_(i)A_(i,m-1)),(sigma_(i) < 0)]:}A_{i, m}= \begin{cases}\max \left(A_{t, \max }, \omega_{i} A_{i, m-1}\right) & \left(\sigma_{i} \geq 0\right) \\ \max \left(A_{c, \text { max }}, \omega_{i} A_{i, m-1}\right) & \left(\sigma_{i}<0\right)\end{cases}

对于减少设计变量：
A i , m = { min ( A t , min , ω i A i , m − 1 ) ( σ i ≥ 0 ) min ( A c , min , ω i A i , m − 1 ) ( σ i < 0 ) A i , m = min A t , min , ω i A i , m − 1 σ i ≥ 0 min A c , min , ω i A i , m − 1 σ i < 0 A_(i,m)={[min(A_(t,min),omega_(i)A_(i,m-1)),(sigma_(i) >= 0)],[min(A_(c,min),omega_(i)A_(i,m-1)),(sigma_(i) < 0)]:}A_{i, m}= \begin{cases}\min \left(A_{t, \min }, \omega_{i} A_{i, m-1}\right) & \left(\sigma_{i} \geq 0\right) \\ \min \left(A_{c, \min }, \omega_{i} A_{i, m-1}\right) & \left(\sigma_{i}<0\right)\end{cases}

其中 A i , m A i , m A_(i,m)A_{i, m} 是第 i i ii 个条在第 m m mm 次迭代中的更新区域，而 A i , m − 1 A i , m − 1 A_(i,m-1)A_{i, m-1} 是来自上一次迭代的区域。

方程（22）和（23）对迭代优化过程至关重要。它们平衡了根据性能标准调整设计变量的需要，同时遵循方程（7）中的面积约束。使用 ω i A i , m − 1 ω i A i , m − 1 omega_(i)A_(i,m-1)\omega_{i} A_{i, m-1} 允许逐步调整，确保优化过程的稳定性。方程（22）确保 A i , m A i , m A_(i,m)A_{i, m} 不超过拉伸或压缩杆的最大允许面积，防止杆面积过大。同样，方程（23）确保 A i , m A i , m A_(i,m)A_{i, m} 不低于拉伸或压缩杆的最小允许面积，防止杆面积过小。

值得注意的是，与通常每次迭代仅改变 2 % 2 % 2%2 \% 个元素的连续 BESO 方法不同，DMT-BESO 方法可以在每次迭代中同时优化所有元素。


3.2.6. 多重载荷情况


考虑到桁架拓扑优化中的多种荷载情况，确保设计的稳健性，并在各种实际场景中表现良好。这种方法增强了优化结构的安全性和可靠性。

(b)

图 5. 设计变量自适应更新策略示意图：(a) 增加和 (b) 减少。

设计。为此，我们使用应力包络法计算每个元素在所有载荷条件下所经历的最大和最小应力。设 σ i s σ i s sigma_(i)^(s)\sigma_{i}^{s} 表示在载荷条件 S S SS 下第 i i ii 个元素的正应力。最大和最小应力 σ max , i σ max  , i sigma_("max ",i)\sigma_{\text {max }, i} 和 σ min , i σ min , i sigma_(min,i)\sigma_{\min , i} 由以下公式给出：
σ max , i = max S ∈ V ( σ i S ) σ max , i = max S ∈ V   σ i S sigma_(max,i)=max_(S inV)(sigma_(i)^(S))\sigma_{\max , i}=\max _{S \in \mathscr{\mathscr { V }}}\left(\sigma_{i}^{S}\right)
σ min , i = min S ∈ Y ( σ i S ) σ min , i = min S ∈ Y   σ i S sigma_(min,i)=min_(S inY)(sigma_(i)^(S))\sigma_{\min , i}=\min _{S \in \mathscr{\mathscr { Y }}}\left(\sigma_{i}^{S}\right)

元素 i i ii 的临界应力 σ i σ i sigma_(i)\sigma_{i} 随后根据最大和最小应力的绝对值确定：
σ i = { σ max , i if | σ max , i | ≥ | σ min , i | σ min , i otherwise σ i = σ max , i  if  σ max , i ≥ σ min , i σ min , i  otherwise  sigma_(i)={[sigma_(max,i)" if "|sigma_(max,i)| >= |sigma_(min,i)|],[sigma_(min,i)" otherwise "]:}\sigma_{i}=\left\{\begin{array}{l}\sigma_{\max , i} \text { if }\left|\sigma_{\max , i}\right| \geq\left|\sigma_{\min , i}\right| \\ \sigma_{\min , i} \text { otherwise }\end{array}\right.

同样，让 α i ( S ) α i ( S ) alpha_(i)^((S))\alpha_{i}^{(S)} 表示在载荷条件 S S SS 下第 i i ii 个元素的灵敏度。元素 i i ii 在所有载荷情况下的最大灵敏度计算为：
α max , i = max S ( α i ( S ) ) α max , i = max S   α i ( S ) alpha_(max,i)=max_(S)(alpha_(i)^((S)))\alpha_{\max , i}=\max _{S}\left(\alpha_{i}^{(S)}\right)

通过识别所有负载条件下的最大灵敏度，可以捕捉到灵敏度的最坏情况，从而确保设计在最敏感的条件下保持稳健。

使用每个元素的最大绝对应力和最大灵敏度是一种保守的方法，确保在最恶劣的加载条件下的稳健性。然而，这种方法在实际工程应用中防止了潜在的失效模式，使结构更加可靠。


3.2.7. 收敛标准


随着有限元分析（FEA）和设计变量更新的循环继续，优化结构的总体积趋于稳定。为了确定程序是否收敛，应计算每次迭代的结构总体积。总体积的相对变化衡量程序的收敛性，表示为：
λ = | ∑ i = 1 N V k − i + 1 − ∑ i = 1 N V k − N − i + 1 | ∑ i = 1 N V k − i + 1 ≤ τ λ = ∑ i = 1 N   V k − i + 1 − ∑ i = 1 N   V k − N − i + 1 ∑ i = 1 N   V k − i + 1 ≤ τ lambda=(|sum_(i=1)^(N)V_(k-i+1)-sum_(i=1)^(N)V_(k-N-i+1)|)/(sum_(i=1)^(N)V_(k-i+1)) <= tau\lambda=\frac{\left|\sum_{i=1}^{N} V_{k-i+1}-\sum_{i=1}^{N} V_{k-N-i+1}\right|}{\sum_{i=1}^{N} V_{k-i+1}} \leq \tau

在优化过程中， k k kk 代表当前迭代次数，而 τ τ tau\tau 是确定优化何时稳定的收敛容忍度。 N N NN 是用于计算平均体积变化的迭代次数。通常， N N NN 设置为 5 ，这意味着监测前 10 次迭代的平均总量变化，以确保其足够小并在指定的容忍度范围内。这确保在结束优化过程之前，设计已经稳定。在所呈现的数值示例中，使用了收敛容忍度 τ = 0.0001 τ = 0.0001 tau=0.0001\tau=0.0001 。

  3.2.8. 程序


DMT-BESO 方法的进化迭代过程包括以下步骤：

(1) 建立 GS，定义载荷、边界条件，并初始化设计变量。

(2) 执行有限元分析以计算每个元素的正常应力。

(3) 根据迭代次数计算优化过程中拉伸和压缩材料性能的变化。

(4) 计算元素灵敏度和灵敏度阈值。

(5) 将拉伸材料属性分配给受拉的杆，将压缩材料属性分配给受压的杆。

(6) 比较元素灵敏度与阈值，并相应调整设计变量。

(7) 使用指定的标准检查收敛性。如果未达到收敛，重复从步骤 (2) 开始的过程。

图 6 显示了 DMT-BESO 方法的流程图。该算法用 Python 编程，并与 Abaqus 软件集成，后者用作有限元分析求解器。


3.3. 结构复杂性控制策略


在桁架 TO 的背景下，控制结构复杂性对于实现实用和可制造的设计至关重要。过于复杂的桁架结构难以制造、组装和维护，导致成本增加。通过减少细杆或非必要杆件的数量来最小化不必要的复杂性，有助于设计简单而优雅的结构，同时保持承载能力和稳健性。

为了实现这一目标，探索一种高效的方法来生成近似最优的桁架设计，而不是严格的理论最优解是至关重要的。因此，我们提出了一种简单有效的方法：对桁架结构中的细杆进行惩罚，并在一定次数的迭代后根据应力限制迭代调整元素面积。“惩罚”是指在优化结构中人为降低某些杆（或元素）的重要性。

惩罚函数旨在减少细杆，这些细杆对整体结构刚度的贡献最小。该过程首先计算所有杆中最大的面积 A max A max A_(max)A_{\max } ，代表结构中最大的横截面积：
A max = max ( A i ) A max  = max A i A_("max ")=max(A_(i))A_{\text {max }}=\max \left(A_{i}\right)

为了确保横截面积较小的杆件受到更严格的惩罚，我们将惩罚阈值定义为 5 % 5 % 5%5 \% 的 A max A max  A_("max ")A_{\text {max }} 。面积低于该阈值的杆件将受到惩罚函数的影响。惩罚函数 P y ( A i , A min ) P y A i , A min P_(y)(A_(i),A_(min))P_{y}\left(A_{i}, A_{\min }\right) 定义如下：
P y ( A i , A min ) = min ( 0.9 , 1 1 + e − 0.01 ( A i / A min ) ) P y A i , A min = min 0.9 , 1 1 + e − 0.01 A i / A min P_(y)(A_(i),A_(min))=min(0.9,(1)/(1+e^(-0.01(A_(i)//A_(min)))))P_{y}\left(A_{i}, A_{\min }\right)=\min \left(0.9, \frac{1}{1+e^{-0.01\left(A_{i} / A_{\min }\right)}}\right)

公式 (30) 计算基于面积 A i A i A_(i)A_{i} 与最小允许面积 A min A min A_(min)A_{\min } 的比率的惩罚因子。它使用类似于 sigmoid 的函数来确保惩罚因子保持在 0.5 到 0.9 的范围内。当 A i A i A_(i)A_{i} 更接近 A min A min A_(min)A_{\min } 时，惩罚因子更接近 0.5，从而更快地惩罚对结构刚度贡献较小的细杆。

如果 A i < 0.05 A max A i < 0.05 A max A_(i) < 0.05A_(max)A_{i}<0.05 A_{\max } 和 α ~ i < α ~ 0 α ~ i < α ~ 0 widetilde(alpha)_(i) < widetilde(alpha)_(0)\widetilde{\alpha}_{i}<\widetilde{\alpha}_{0} 都满足，则条形的面积将使用公式(30)计算的惩罚因子进行惩罚：
A i = { A i Py ( A i , A t , min ) if σ i ≥ 0 A i Py ( A i , A c , min ) otherwise A i = A i Py A i , A t ,  min   if  σ i ≥ 0 A i Py A i , A c ,  min   otherwise  A_(i)={[A_(i)Py(A_(i),A_(t," min "))," if "sigma_(i) >= 0],[A_(i)Py(A_(i),A_(c," min "))," otherwise "]:}A_{i}= \begin{cases}A_{i} \operatorname{Py}\left(A_{i}, A_{t, \text { min }}\right) & \text { if } \sigma_{i} \geq 0 \\ A_{i} \operatorname{Py}\left(A_{i}, A_{c, \text { min }}\right) & \text { otherwise }\end{cases}

然而，对细杆的惩罚可能会导致其余杆件的应力超限。因此，在优化过程中稍后对应力超限的杆件进行迭代调整。应力超限杆件的累积面积增加更新方法如下：
φ i ← φ i ( 1.0 + 0.01 e − 0.05 ( m − k ) ) φ i ← φ i 1.0 + 0.01 e − 0.05 ( m − k ) varphi_(i)larrvarphi_(i)(1.0+0.01e^(-0.05(m-k)))\varphi_{i} \leftarrow \varphi_{i}\left(1.0+0.01 \mathrm{e}^{-0.05(m-k)}\right)

其中 φ i φ i varphi_(i)\varphi_{i} 是累积面积增加因子， e − 0.05 ( m − k ) e − 0.05 ( m − k ) e^(-0.05(m-k))\mathrm{e}^{-0.05(m-k)} 确保随着 m m mm 超过 k k kk 调整变得不那么激进，从而允许在优化过程中逐渐收敛。

基于 φ i φ i varphi_(i)\varphi_{i} 的区域更新为：
A i = { min ( A t , max , A i φ i ) if σ i > σ 0 , t min ( A c , max , A i φ i ) if σ i < σ 0 , c A i = min A t ,  max  , A i φ i  if  σ i > σ 0 , t min A c ,  max  , A i φ i  if  σ i < σ 0 , c A_(i)={[min(A_(t," max "),A_(i)varphi_(i))," if "sigma_(i) > sigma_(0,t)],[min(A_(c," max "),A_(i)varphi_(i))," if "sigma_(i) < sigma_(0,c)]:}A_{i}=\left\{\begin{array}{lll}\min \left(A_{t, \text { max }}, A_{i} \varphi_{i}\right) & \text { if } \sigma_{i}>\sigma_{0, t} \\ \min \left(A_{\mathrm{c}, \text { max }}, A_{i} \varphi_{i}\right) & \text { if } \sigma_{i}<\sigma_{0, c}\end{array}\right.

这种方法在保持结构可行性与优化目标之间取得平衡。通过根据应力条件应用累积面积增加更新，确保超应力的构件以受控方式进行调整，防止剧烈变化，同时稳步引导设计朝向最佳应力分布。指数衰减

图 6. DMT-BESO 方法的流程图。

该函数用于在优化过程中逐渐减少调整因子。这确保了早期和后期阶段之间的平滑过渡——早期阶段需要较大的更新以快速探索，而后期阶段则需要更精细的调整以防止振荡并提高稳定性。指数衰减提供的平滑减少使优化过程既稳健又可靠，通过平衡探索和收敛。

图 7 显示了结构复杂性控制策略的工作流程。该工作流程整合了惩罚和迭代调整。对低灵敏度的细杆施加适当的惩罚可以减少过于纤细的杆件数量，从而导致更简单、更优雅的优化结构。在优化的后期阶段，对超应力元素区域的迭代调整确保了遵守应力约束，从而在拓扑优化结果中保持结构安全。

图 7. 结构复杂性控制策略的工作流程。

  4. 数值示例


本节的主要目的是通过一系列数值示例来证明上述方法的有效性、稳健性和准确性。这些示例旨在突出 DMT-BESO 方法通过利用两种不同材料的独特性质来优化桁架结构的能力。此外，这些示例展示了该方法处理复杂几何形状的能力，确保符合行业标准，并通过实施的复杂性控制策略保持结构的简单性和可制造性。


4.1. 操作环境和参数设置


为了验证 DMT-BESO 方法的有效性，我们在以下研究案例中进行了单材料和双材料基础桁架拓扑优化的比较研究。所有数值示例均在一台配备 12 代 Intel® Core™ i7-12700 CPU @ 2.10 GHz 和 32 GB RAM 的台式计算机上进行，运行 64 位 Windows 10 专业版操作系统。拓扑优化过程使用 Abaqus 2023 版本实现。

对于双材料桁架 TO，采用 Q235 钢材作为拉伸材料，C50 混凝土作为压缩材料。对于单材料桁架 TO，使用两种具有相同机械性能的材料 Q235 钢材，以便通过 DMT-BESO 方法轻松便捷地实现优化。表 2 列出了根据[50,51]的钢材和混凝土材料的性质。在单材料和双材料情况下，每根杆件的初始截面积均设定为 100 mm 2 100 mm 2 100mm^(2)100 \mathrm{~mm}^{2} ，外部载荷为 P = 10 kN P = 10 kN P=10kNP=10 \mathrm{kN} 。截面积约束的最小值和最大值分别设定为 A t , min = A c , min = 1.0 × 10 − 6 A t , min = A c ,  min  = 1.0 × 10 − 6 A_(t,min)=A_(c," min ")=1.0 xx10^(-6)A_{t, \min }=A_{c, \text { min }}=1.0 \times 10^{-6} mm 2 mm 2 mm^(2)\mathrm{mm}^{2} 和 A t , max = A c , max = 1.0 × 10 6 mm 2 A t ,  max  = A c ,  max  = 1.0 × 10 6 mm 2 A_(t," max ")=A_(c," max ")=1.0 xx10^(6)mm^(2)A_{t, \text { max }}=A_{c, \text { max }}=1.0 \times 10^{6} \mathrm{~mm}^{2} 。
  表 2

材料的相关属性。
  材料
弹性模量 E ( MPa ) E ( MPa ) E(MPa)E(\mathrm{MPa})
允许应力 [ σ [ σ [sigma[\sigma ] (MPa)
  Q235 钢 2.06 × 10 5 2.06 × 10 5 2.06 xx10^(5)2.06 \times 10^{5} 215
  C60 混凝土 3.6 × 10 4 3.6 × 10 4 3.6 xx10^(4)3.6 \times 10^{4} 27.5
Materials Modulus of elasticity E(MPa) Allowable stress [sigma ] (MPa) Q235 steel 2.06 xx10^(5) 215 C60 concrete 3.6 xx10^(4) 27.5| Materials | Modulus of elasticity $E(\mathrm{MPa})$ | Allowable stress $[\sigma$ ] (MPa) | | :--- | :--- | :--- | | Q235 steel | $2.06 \times 10^{5}$ | 215 | | C60 concrete | $3.6 \times 10^{4}$ | 27.5 |

注意：允许应力取决于材料的设计强度。


4.2. 2 D 2 D 2D2 D 简单支撑桁架


在第一个数值示例中，优化了一个二维简支桁架。使用第 2.2 节中介绍的第一种方法生成的 GS 如图 8(a)所示。GS 的基础网格尺寸为 400 mm × 400 mm 400 mm × 400 mm 400mmxx400mm400 \mathrm{~mm} \times 400 \mathrm{~mm} ，由 99 根杆件组成。图 8(b)和图 8(c)显示，优化后的结构与传统桁架结构非常相似。此外，优化结构中杆件的正常应力接近各自的应力约束。图 8(b)中优化结构的体积和应变能与文献[37]中相同案例的结果误差在 0.5%以内。图 8(b)显示，拉伸和压缩杆件的厚度相当，而图 8(c)中的压缩杆件由于 C50 混凝土的允许应力显著低于 Q235 钢，因此厚度远大于拉伸杆件。单材料和双材料案例的拓扑优化过程分别在 408 和 431 次迭代后收敛，总耗时分别为 99 分钟和 106 分钟。


4.3. 2 D 2 D 2D2 D 米歇尔桁架


作为一个成熟的最优解，Michell 桁架作为验证和比较新优化算法和方法的基准。根据文献 [28]，其最小体积的解析解为：
V anal = P ( L x 2 ) ( 1 2 + π 4 ) [ 1 σ + + 1 σ − ] V anal  = P L x 2 1 2 + π 4 1 σ + + 1 σ − V_("anal ")=P((L_(x))/(2))((1)/(2)+(pi)/(4))[(1)/(sigma^(+))+(1)/(sigma^(-))]V_{\text {anal }}=P\left(\frac{L_{x}}{2}\right)\left(\frac{1}{2}+\frac{\pi}{4}\right)\left[\frac{1}{\sigma^{+}}+\frac{1}{\sigma^{-}}\right]

其中 P P PP 表示集中力值， L x L x L_(x)L_{x} 和 L y L y L_(y)L_{y} 分别表示设计域的长度和宽度。

在第二个数值示例中，我们使用 Michell 桁架验证了 DMT-BESO 方法的准确性和效率。如图 9(a)和(b)所示，设计域的长度为 5 米，宽度为 2 米，基础网格大小为 250 mm × 250 mm 250 mm × 250 mm 250mmxx250mm250 \mathrm{~mm} \times 250 \mathrm{~mm} ，由 160 个基础网格和 189 个节点组成。在生成的 GS 中，每个节点与最多三个相邻基础网格中的所有节点相连， resulting in a total of 3116 bars.

图 9©至(f)展示了使用单材料和双材料方法进行 Michell 桁架优化的结果，包含和不包含结构复杂性控制。如图 9©和(e)所示，在没有复杂性控制的情况下，优化后的结构表现出密集的细杆。

图 8. 简支桁架的优化：(a) 带有尺寸、载荷和边界条件的 GS；(b) 单材料结果；(c) 双材料结果（单位：尺寸为毫米，应力为兆帕）。

围绕主要承载应力的构件，无论是单材料还是双材料情况。相反，如图 9(d)和(f)所示，通过复杂性控制，优化后的结构在保持效率的同时显著简化，更适合实际应用。

图 10 显示了具有复杂性控制的 Michell 桁架优化的演变历史。在图 10(a)和 10(b)中，应变能和体积最初迅速演变，表明较大的步长更新有效地促进了优化过程朝向可行解的方向发展。后续迭代中逐渐减小的步长确保了准确的收敛。在整个优化过程中，特别是在早期阶段，没有观察到严重的振荡，证明了第 3.3 节中描述的材料分配方法的可行性和有效性。在双材料情况下的应变能低于单材料情况，充分利用了两种材料的不同特性，以实现更优的整体结构性能。

在这个数值示例中，最小体积的解析解为 2.989 × 10 5 mm 3 2.989 × 10 5 mm 3 2.989 xx10^(5)mm^(3)2.989 \times 10^{5} \mathrm{~mm}^{3} ，在全应力下，结构应变能为 33.539 J，为评估优化性能提供了理论基准。表 3 比较了单材料和双材料情况下，优化体积 ( V opt V opt  V_("opt ")V_{\text {opt }} ) 和结构应变能 ( C opt C opt  C_("opt ")C_{\text {opt }} ) 在有无复杂性控制下的情况。单材料情况下，无论是否有复杂性控制， V opt V opt  V_("opt ")V_{\text {opt }} 和 C opt C opt  C_("opt ")C_{\text {opt }} 的值都相对接近解析解。在没有复杂性控制的情况下， V o p t V o p t V_(opt)V_{o p t} 和 C o p t C o p t C_(opt)C_{o p t} 的值分别比解析解高出约 1.19 % 1.19 % 1.19%1.19 \% 和 1.05 % 1.05 % 1.05%1.05 \% 。在有复杂性控制的情况下， V o p t V o p t V_(opt)V_{o p t} 和 C o p t C o p t C_(opt)C_{o p t} 的值分别比解析解高出约 2.04 % 2.04 % 2.04%2.04 \% 和 0.19%。将双材料情况与单材料情况进行比较，前者的 V o p t V o p t V_(opt)V_{o p t} 和 C o p t C o p t C_(opt)C_{o p t} 值差异显著。这一差异归因于混凝土的弹性模量和允许应力明显低于钢材。 然而， V opt V opt  V_("opt ")V_{\text {opt }} 和 C o p t C o p t C_(opt)C_{o p t} 的值在有复杂性控制和没有复杂性控制的情况下是可比的，这表明复杂性控制略微增加了材料使用，但确保了更可行和可制造的设计。


4.4. 2 D 2 D 2D2 D 大麻悬臂


如图 11(a)和(b)所示，Hemp 悬臂桁架在右侧中部承受集中载荷 P P PP ，并由左侧的两个固定铰链支撑。设计域为 3.6 m × 3.6 m 3.6 m × 3.6 m 3.6mxx3.6m3.6 \mathrm{~m} \times 3.6 \mathrm{~m} ，基础网格大小为 300 mm × 300 mm 300 mm × 300 mm 300mmxx300mm300 \mathrm{~mm} \times 300 \mathrm{~mm} ，由 144 个基础网格和 169 个节点组成。生成的 GS 中的每个节点与最多四个相邻基础网格中的所有节点连接，形成总共 2840 根杆件。根据 Hemp [52]，悬臂的最佳分析体积为 4.32168 P L / σ 4.32168 P L / σ 4.32168 PL//sigma4.32168 P L / \sigma 。因此，本数值示例中的分析解为 7.236 × 10 5 mm 3 7.236 × 10 5 mm 3 7.236 xx10^(5)mm^(3)7.236 \times 10^{5} \mathrm{~mm}^{3} ，对应的应变能为 81.119 J。对于复杂性控制的单材料案例，优化结构的体积和应变能分别为 7.467 × 10 5 mm 3 7.467 × 10 5 mm 3 7.467 xx10^(5)mm^(3)7.467 \times 10^{5} \mathrm{~mm}^{3} 和 82.526 J。与分析解相比，增加量仅为 3.2 % 3.2 % 3.2%3.2 \% 和 1.6 % 1.6 % 1.6%1.6 \% 。这再次证明了 DMT-BESO 方法和结构复杂性控制策略在生成近似最佳桁架设计方面的准确性和有效性。

图 11©和(d)展示了在结构复杂性控制下，单材料和双材料情况下的优化结果。优化后的结构从载荷施加点到固定支撑展现出清晰、高效的载荷路径。此外，两种优化场景中的应力分布控制良好，杆件应力接近各自的约束，表明材料特性的有效利用。在单材料情况下，拉伸（红色）和压缩（蓝色）杆件的尺寸相同，拓扑形态对称。相反，图 11(d)中的双材料结果显示出显著不同的尺寸和不对称性。

图 9. Michell 桁架的优化：(a) 基础网格尺寸；(b) 带载荷和边界条件的 GS；(c) 和(d) 无/有结构复杂性控制的单材料结果；(e) 和(f) 无/有结构复杂性控制的双材料结果（单位：尺寸为毫米，压力为兆帕）。

图 10. Michell 桁架的演化历史比较：(a) 单材料与结构复杂性控制；(b) 双材料与结构复杂性控制。
  表 3

单材料和双材料案例的比较，带有/不带复杂性控制。
  研究案例
没有复杂性控制		   带有复杂性控制
V opt ( mm 3 ) V opt  mm 3 V_("opt ")(mm^(3))V_{\text {opt }}\left(\mathrm{mm}^{3}\right) C opt ( J ) C opt  ( J ) C_("opt ")(J)C_{\text {opt }}(\mathrm{J})   总时间（分钟）   总迭代次数 V opt ( mm 3 ) V opt  mm 3 V_("opt ")(mm^(3))V_{\text {opt }}\left(\mathrm{mm}^{3}\right) C opt C opt  C_("opt ")C_{\text {opt }} (J)   总时间（分钟）   总迭代次数
  单一材料案例 3.024 × 10 5 3.024 × 10 5 3.024 xx10^(5)3.024 \times 10^{5} 33.892 270 406 3.05 × 10 5 3.05 × 10 5 3.05 xx10^(5)3.05 \times 10^{5} 33.603 273 410
  双材料外壳 16.843 × 10 5 16.843 × 10 5 16.843 xx10^(5)16.843 \times 10^{5} 27.751 279 409 16.974 × 10 5 16.974 × 10 5 16.974 xx10^(5)16.974 \times 10^{5} 27.560 295 432
Study case Without complexity control With complexity control V_("opt ")(mm^(3)) C_("opt ")(J) Total time (minute) Total iterations V_("opt ")(mm^(3)) C_("opt ") (J) Total time (minute) Total iterations Single-material case 3.024 xx10^(5) 33.892 270 406 3.05 xx10^(5) 33.603 273 410 Dual- material case 16.843 xx10^(5) 27.751 279 409 16.974 xx10^(5) 27.560 295 432| Study case | Without complexity control | | | | With complexity control | | | | | :---: | :---: | :---: | :---: | :---: | :---: | :---: | :---: | :---: | | | $V_{\text {opt }}\left(\mathrm{mm}^{3}\right)$ | $C_{\text {opt }}(\mathrm{J})$ | Total time (minute) | Total iterations | $V_{\text {opt }}\left(\mathrm{mm}^{3}\right)$ | $C_{\text {opt }}$ (J) | Total time (minute) | Total iterations | | Single-material case | $3.024 \times 10^{5}$ | 33.892 | 270 | 406 | $3.05 \times 10^{5}$ | 33.603 | 273 | 410 | | Dual- material case | $16.843 \times 10^{5}$ | 27.751 | 279 | 409 | $16.974 \times 10^{5}$ | 27.560 | 295 | 432 |

图 11. Hamp 悬臂梁的优化：(a) 基础网格尺寸；(b) 带载荷和边界条件的 GS；(c) 具有结构复杂性控制的单材料结果；(d) 具有结构复杂性控制的双材料结果（单位：尺寸为毫米，压力为 MPa）。

与图 11©中的单材料结果相比，拉伸和压缩杆的拓扑形态。图 11(d)中的压缩杆比图 11©中的更粗壮。这与使用多材料连续体 BESO 方法获得的结果非常相似。双材料设计使得压缩杆能够保持较低的应力水平并使用更大的截面积，从而增加了弯矩惯性并减少了长细比。这大大增强了结构的整体稳定性。图 12 显示了复杂性控制下的麻质悬臂优化的演变历史，规律性与之前的例子一致。TO 过程在单材料和双材料情况下分别经过 431 和 452 次迭代后收敛，总耗时分别为 256 和 271 分钟。


4.5. 3 D 3 D 3D3 D 空间穹顶


第四个数值例子验证了 DMT-BESO 方法在空间桁架中的有效性。通过将图 11(a)所示的平面网格投影到 z = 1.5 ( x 2 + y 2 ) z = 1.5 x 2 + y 2 z=1.5(x^(2)+y^(2))z=1.5\left(x^{2}+y^{2}\right) 空间表面，可以获得图 13(a)所示的空间网格。因此，使用第 2.2 节中介绍的方法，生成了与图 11(b)对应的这个空间穹顶的 GS，如图 13(b)所示。在杆件数量和连接层次方面，它与第三个数值例子相匹配。不同之处在于，集中力施加在最高点的节点上，而四个固定铰支撑设置在最低点。

图 13©和(f)中的单材料和双材料优化结果展示了从复杂的 3D GS 中确定最有效力传递路径的能力。结构复杂性控制方法有效地简化了设计，同时保持了

图 12. Hamp 悬臂梁的演化历史比较：(a) 具有结构复杂性控制的单材料；(b) 具有结构复杂性控制的双材料。

图 13. 3D 空间穹顶的优化：(a) 基础网格及其对应的空间表面；(b) 带有载荷和边界条件的 GS；(c) 到(e) 单材料结果与结构复杂性控制；(f) 到(h) 双材料结果与结构复杂性控制（单位：位移为 mm，压力为 MPa）。

优化结构的完整性。这种简化不仅保留了结构性能，还促进了更简单和更具成本效益的施工实践。图 13(d)和(g)中的应力颜色映射清楚地显示出杆件应力在安全范围内，证明了优化设计的可靠性和安全性。尽管钢的弹性模量高于混凝土，但图 13(e)和(h)之间的比较表明，双材料情况下的结构刚度大于单材料情况下的刚度，Z 方向的位移证明了这一点。单材料和双材料情况下的拓扑优化过程分别在 406 和 408 次迭代后收敛，总耗时分别为 239 和 243 分钟。


5. 实际应用


在本节中，我们应用了 DMT-BESO 方法对单跨拱桥进行结构优化，重点优化吊杆布置以提高结构性能。该桥设计主跨为 120 米，宽度为 10 米，采用内倾的篮把拱形设计，以增强稳定性和美观性。复合结构采用超高性能混凝土（UHPC），因其优异的抗压强度，以及高强度钢材，因其优越的抗拉性能，具体见表 4。这些材料的弹性模量和允许应力参考自 [ 53 , 54 ] [ 53 , 54 ] [53,54][53,54] 。

为了对这座桥进行结构优化，考虑的主要荷载条件是恒载和车辆荷载，这些荷载被转换为施加在拉杆节点上的集中荷载。此外，为了简化计算，同时提供安全可靠的结果，
  表 4

主要建筑材料的参数值。
  材料

弹性模量 E E EE ( MPa ) ( MPa ) (MPa)(\mathrm{MPa})
Modulus of elasticity E (MPa)| Modulus of elasticity $E$ | | :--- | | $(\mathrm{MPa})$ |

允许应力 [ σ ] [ σ ] [sigma][\sigma] ( MPa ) ( MPa ) (MPa)(\mathrm{MPa})
Allowable stress [sigma] (MPa)| Allowable stress $[\sigma]$ | | :--- | | $(\mathrm{MPa})$ |

高强度钢 (Q550)
High-strength steel (Q550)| High-strength steel | | :--- | | (Q550) |
 2.06 × 10 5 2.06 × 10 5 2.06 xx10^(5)2.06 \times 10^{5} 500
UHPC (UC140) 4.43 × 10 4 4.43 × 10 4 4.43 xx10^(4)4.43 \times 10^{4} 68
Materials "Modulus of elasticity E (MPa)" "Allowable stress [sigma] (MPa)" "High-strength steel (Q550)" 2.06 xx10^(5) 500 UHPC (UC140) 4.43 xx10^(4) 68| Materials | Modulus of elasticity $E$ <br> $(\mathrm{MPa})$ | Allowable stress $[\sigma]$ <br> $(\mathrm{MPa})$ | | :--- | :--- | :--- | | High-strength steel <br> (Q550) | $2.06 \times 10^{5}$ | 500 | | UHPC (UC140) | $4.43 \times 10^{4}$ | 68 |

车辆静载荷近似和适当的冲击系数用于捕捉动态载荷的影响。由于桥梁结构在使用过程中会经历可变载荷，导致不同的应力分布和潜在的失效模式，因此针对单一载荷条件的优化可能导致在不同或组合载荷下的设计失效。因此，在优化过程中，除了应用结构复杂性控制以简化结构外，我们还采用了多种载荷情况，即将恒载和车辆移动载荷叠加为五个主要载荷情况，以模拟车辆在桥上行驶的情况。TO 过程在 419 次迭代后收敛，总耗时为 216 分钟。

图 14(a)所示的弓形 3D 空间 GS 由 814 根杆件组成，采用第 2.2 节中介绍的第二种方法生成。如图 14(b)所示，节点仅分布在 GS 的顶部和底部边缘，而杆件相交但不形成额外的交点。每根杆件的初始截面积以及截面积约束的最小值和最大值与第 4.1 节中描述的相符。图 14(c)至(g)展示了在五个主要载荷情况下的边界条件下的 GS 及其相应的优化设计。应力颜色图表明，所有杆件的应力水平均保持在可接受的范围内，确保桥梁能够安全承受施加的载荷，同时有效利用材料。优化后的结构符合通过系杆拱桥的基本力学特性，超高性能混凝土（UHPC）分配给拱肋，高强度钢分配给系梁和吊杆，确保每个组件适合其各自的力学状态。

通过将优化后的结构从 Abaqus 导出为*.odb 文件格式，并将 CAD 模型纳入 3D 环境，可以进一步分析和评估设计。图 15 所示的渲染突显了这座桥的优雅和实用性，展示了新型结构优化技术与实用设计美学的无缝结合。

  6. 结论


本文提出了一种新颖的进化拓扑优化方法

图 14. 单跨拱桥的优化：(a) GS 的透视图和立面图；(b) GS 中节点和杆件的分布；(c) 到(g) 在多种荷载情况下的 GS 及其边界条件，以及相应的优化设计（单位：尺寸为毫米，压力为兆帕）。

图 15. 优化的单跨拱桥的渲染。

方法，双材料桁架 BESO（DMT-BESO），通过同时使用具有不同拉伸和压缩特性的两种材料来优化桁架结构。可以得出以下结论：

（1）桁架拓扑优化的实际设计考虑因素包括可制造性、简洁性、材料效率、特定约束和与现有工作流程的平滑集成。在理论优化与实际实施之间取得平衡需要一种全面的方法。

(2) 一种新的生成任意形状地面结构的方法已经开发，利用了 Grasshopper 和 Rhino 的参数化能力。Abaqus、Grasshopper 和 Python 的集成使得处理复杂几何形状和解决实际工程问题变得更加方便。

通过采用最小能量原理和全应力设计标准，DMT-BESO 方法实现了拓扑和尺寸的同时优化，而应力约束确保符合行业标准，促进了安全性和材料效率。

(4) DMT-BESO 方法及其结构复杂性控制策略的有效性已通过多个数值示例得到验证，包括复合材料拱桥的设计，展示了其提供高效且实用的桁架结构的能力。

DMT-BESO 方法为实际桁架设计提供了一种强大的工具，结合了理论基础和实际设计考虑。未来的深入研究可以集中在将欧拉屈曲约束纳入优化过程，以提高计算效率，并将应用扩展到更大、更复杂的结构系统。


CRediT 作者贡献声明


赖雅平：撰写 - 原始草稿，可视化，监督，软件，方法论，概念化。蔡琦：验证，概念化。李宇：验证，概念化。陈家勇：软件。谢怡敏：撰写 - 审阅与编辑，监督，方法论，概念化。


利益冲突声明


作者声明他们没有已知的竞争性财务利益或个人关系，这些关系可能会影响本文所报告的工作。

  致谢


该项目得到了澳大利亚研究委员会（FL190100014）和中国博士后科学基金（GZC20241226）博士后奖学金项目的支持。作者感谢匿名评审人，他们的建设性意见无疑有助于提高手稿的整体质量。


附录 A. 支持信息


与本文相关的补充数据可以在在线版本中找到，网址为 doi:10.1016/j.engstruct.2024.119326。

  数据可用性


数据将根据请求提供。

  参考文献


[1] Bendsøe MP, Kikuchi N. 使用均匀化方法在结构设计中生成最优拓扑。计算方法与应用机械工程 1988;71:197-224. https://doi.org/10.1016/0045-7825(88)90086-2.

[2] Bendsøe MP, Sigmund O. 拓扑优化中的材料插值方案。应用力学档案 1999;69:635-54。https://doi.org/10.1007 s004190050248。

[3] Rietz A. 在 SIMP（幂律）方法中有限指数的充分性。结构多学科优化 2001;21:159-63。https://doi.org/10.1007/s001580050180。

[4] Andreassen E, Clausen A, Schevenels M, Lazarov BS, Sigmund O. 在 MATLAB 中使用 88 行代码进行高效拓扑优化。Struct Multidiscip Optim 2011;43: 1-16. https://doi.org/10.1007/s00158-010-0594-7.

[5] Xie YM, Steven GP. 一种简单的进化程序用于结构优化。计算结构 1993;49:885-96。https://doi.org/10.1016/0045-7949(93)90035-C。

[6] Xie YM, Steven GP. 进化结构优化. 施普林格; 1997.

[7] Huang X, Xie Y-M. 对拓扑优化的 ESO 类型方法的进一步回顾。Struct Multidiscip Optim 2010;41:671-83. https://doi.org/10.1007/ s00158-010-0487-9.

[8] Wang MY, Wang X, Guo D. 一种用于结构拓扑优化的水平集方法。计算方法与应用机械工程 2003;192:227-46。https://doi.org/10.1016/ S0045-7825(02)00559-5。

[9] van Dijk NP, Maute K, Langelaar M, van Keulen F. 结构拓扑优化的水平集方法：综述。Struct Multidisc Optim 2013;48:437-72. https:// doi.org/10.1007/s00158-013-0912-y.

[10] Bourdin B, Chambolle A. 最优设计中的相场法。IUTAM 研讨会。拓扑设计。优化。结构。机械。材料。施普林格；2006 年。第 207-15 页。

[11] Mukherjee S, Lu D, Raghavan B, Breitkopf P, Dutta S, Xiao M, 等。加速大规模拓扑优化：现状与挑战。Arch Comput Methods Eng 2021;28:4549-71。https://doi.org/10.1007/s11831-021-09544-3。

[12] Das S, Sutradhar A. 功能梯度可控多孔结构的多物理场拓扑优化：应用于散热问题。Mater Des 2020;193:108775. https://doi.org/10.1016/j.matdes.2020.108775.

[13] Wu J, Sigmund O, Groen JP. 多尺度结构的拓扑优化：综述。结构多学科优化 2021;63:1455-80。https://doi.org/10.1007/ s00158-021-02881-8。

[14] Zhang S, Li H, Huang Y. 一种基于 SIMP 方法的改进多目标拓扑优化模型，用于包括自重的连续体结构。Struct Multidiscip Optim 2021;63:211-30. https://doi.org/10.1007/s00158-020-026852.

[15] Zheng J, Luo Z, Jiang C, Gao J. 针对混合不确定性下动态结构的并行设计的鲁棒拓扑优化。机械系统与信号处理 2019; 120:540-59. https://doi.org/10.1016/j.ymssp.2018.10.026.

[16] Gao K, Do DM, Chu S, Wu G, Kim HA, Featherston CA. 在不确定的模糊随机场传播下，结构的稳健拓扑优化。薄壁结构 2022;175:109238。https://doi.org/ 10.1016/j.tws.2022.109238。

[17] Zhu J-H, Zhang W-H, Xia L. 飞机和航空航天结构设计中的拓扑优化。Arch Comput Methods Eng 2016;23:595-622. https://doi.org/ 10.1007/s11831-015-9151-2.

[18] Munk DJ, Auld DJ, Steven GP, Vio GA. 关于将拓扑优化应用于飞机组件结构设计的好处。Struct Multidiscip Optim 2019;60:1245-66. https://doi.org/10.1007/s00158-019-02250-6.

[19] Aulig N, Nutwell E, Menzel S, Detwiler D. 基于偏好的拓扑优化用于车辆概念设计，考虑静态和碰撞载荷情况。Struct Multidiscip Optim 2018;57:251-66. https://doi.org/10.1007/s00158-017-1751-z.

[20] Sotiropoulos S, Kazakis G, Lagaros ND. 基于拓扑优化和预制构件的结构系统概念设计。计算结构 2020;226:106136。https://doi.org/10.1016/j.compstruc.2019.106136。

[21] Li Y, Ding J, Zhang Z, Zhou X, Makvandi M, Yuan PF, 等. 多材料拓扑优化在标志性建筑性能基础设计中的实际应用. 复合材料结构 2023;325:117603. https://doi.org/10.1016/j. compstruct.2023.117603.

[22] Lai Y, Li Y, Liu Y, Chen P, Zhao L, Li J, 等. 双向进化结构优化在创新人行桥设计中的应用. AI Civ Eng 2024;3:9. https://doi.org/10.1007/s43503-024-00027-5.

[23] Yan X, Bao D, Zhou Y, Xie Y, Cui T. 建筑设计与开发中的拓扑优化详细控制策略。前沿建筑研究 2022;11: 340-56. https://doi.org/10.1016/j.foar.2021.11.001.

[24] Mizobuti V, Vieira Junior LCM. 基于结构拓扑优化的生物启发建筑设计。前沿建筑研究 2020;9:264-76。https://doi.org/ 10.1016/j.foar.2019.12.002。

[25] Li Y, Lai Y, Lu G, Yan F, Wei P, Xie YM. 使用多材料拓扑优化的长跨钢混凝土组合桥的创新设计。工程结构 2022;269:114838. https://doi.org/10.1016/j.engstruct.2022.114838.

[26] 麦克斯韦 JC. I. 关于互惠图形、框架和力的图示。转自爱丁堡皇家学会 1870;26:1-40. https://doi.org/10.1017/S0080456800026351.

[27] Cilley FH. 静不确定框架的精确设计。对其可能性但无用性的阐述。美国土木工程师学会会刊 1900;43:353-407。https://doi. org/10.1061/TACEAT.0001397。

[28] Michell AGM. LVIII. 框架结构中材料经济性的限制。伦敦爱丁堡都柏林哲学杂志科学 1904;8:589-97。https://doi.org/10.1080/ 14786440409463229。

[29] Dorn W., Gomory R.E., Grenberg H. 最优结构的自动设计 1964。

[30] Martínez P, Martí P, Querin OM. 用于桁架结构的尺寸、拓扑和几何优化的生长方法。结构多学科优化 2007;33:13-26. https:// doi.org/10.1007/s00158-006-0043-9.

[31] McKeown JJ. 生长最佳铰接框架。结构优化 1998;15:92-100。https://doi.org/10.1007/BF01278495。

[32] 规则 WK。通过优化生长进行自动桁架设计。J Struct Eng 1994;120: 3063-70。https://doi.org/10.1061/(ASCE)0733-9445(1994)120:10(3063)。

[33] He L, Li Q, Gilbert M, Shepherd P, Rankine C, Pritchard T, 等。基于优化驱动的桁架结构概念设计在参数化建模环境中的应用。Structures 2022;37:469-82。https://doi.org/10.1016/j. istruc.2021.12.048。

[34] Gilbert M, Tyas A. 大规模销接框架的布局优化。工程计算 2003;20:1044-64。https://doi.org/10.1108/02644400310503017。

[35] Wang Q, Arora JS. 使用稀疏 SAND 公式优化大规模桁架结构。国际数值方法工程 2007;69:390-407。

[36] Lu H, Xie YM. 在桁架布局优化中减少不同成员的数量。Struct Multidiscip Optim 2023;66:52. https://doi.org/10.1007/ s00158-023-03514-y.

[37] 蔡琦, 冯锐. 基于改进的双向进化结构优化方法的桁架结构拓扑优化. 建筑结构学报 2022;43: 68-76. https://doi.org/10.14006/j.jzjgxb.2020.0704.

[38] 蔡琦, 冯睿, 张志. 考虑实际局部屈曲稳定性的桁架拓扑优化. 结构 2022;41:1710-8. https://doi. org/10.1016/j. istruc.2022.05.109.

[39] Weldeyesus AG, Gondzio J, He L, Gilbert M, Shepherd P, Tyas A. 桁架几何和拓扑优化与全球稳定性约束。结构多学科优化 2020;62:1721-37。https://doi.org/10.1007/s00158-020-02634-z。

[40] Dang KD, Nguyen-Van S, Thai S, Lee S, Luong VH, Lieu QX. 一种使用混合差分进化和共生生物搜索的桁架拓扑、尺寸和形状的单步优化方法。计算结构 2022;270: 106846. https://doi.org/10.1016/j.compstruc.2022.106846.

[41] Lieu QX. 一种新颖的拓扑框架，用于在静态、自由振动和瞬态行为下对桁架进行拓扑、尺寸和形状的同时优化。工程计算 2022;38:1-25。https://doi.org/10.1007/s00366-022-01599-5。

[42] Toragay O, Silva DF, Vinel A, Shamsaei N. 精确的框架结构全局优化用于增材制造。Struct Multidiscip Optim 2022;65:97. https://doi.org/10.1007/s00158-022-03178-0.

[43] Zegard T, Paulino GH. GRAND - 基于地面结构的任意二维域的拓扑优化，使用 MATLAB。结构多学科优化 2014;50:861-82。https://doi.org/10.1007/s00158-014-1085-z。

[44] Zegard T, Paulino GH. GRAND3 - 基于地面结构的拓扑优化，用于任意 3D 域，使用 MATLAB。结构多学科优化 2015;52: 1161-84。https://doi.org/10.1007/s00158-015-1284-2。

[45] 黄旭, 谢勇明. 连续体结构的进化拓扑优化：方法与应用. Wiley; 2010.

[46] 黄旭, 谢勇明. 双向进化结构优化方法的收敛性和网格无关解. 有限元分析与设计 2007;43:1039-49. https://doi.org/10.1016/j.finel.2007.06.006.

[47] Xia L, Xia Q, Huang X, Xie YM. 关于先进结构和材料的双向进化结构优化：一项综合评述。Arch Comput Methods Eng 2018;25:437-78. https://doi.org/10.1007/s11831-016-9203-2.

[48] Ahrari A, Atai AA. 钢架结构形状和尺寸优化的完全应力设计演变策略。计算结构 2013;123:58-67。https://doi.org/ 10.1016/j. compstruc.2013.04.013。

[49] He L, Gilbert M, Song X. 一种用于桁架自适应布局优化的 Python 脚本。结构多学科优化 2019;60:835-47。https://doi.org/10.1007/ s00158-019-02226-6。

[50] 中华人民共和国住房和城乡建设部。钢结构设计标准（GB 50017-2017）。北京：中国建筑工业出版社；2017。

[51] 中华人民共和国住房和城乡建设部。混凝土结构设计规范（GB 50010-2010）。北京：中国建筑工业出版社；2010。

[52] Hemp WS. Michell 框架在固定支撑之间的均匀载荷。工程优化 1974;1:61-9。https://doi.org/10.1080/03052157408960577。

[53] 中华人民共和国住房和城乡建设部. 高强度钢结构设计标准 (JGJ/T483-2020). 北京: 中国建筑工业出版社; 2020.

[54] 中国建筑材料联合会。《超高性能混凝土结构设计规范》（T/CCPA 35-2022）。北京：中国建筑工业出版社；2022。

  2. 注意： V opt V opt  V_("opt ")V_{\text {opt }} 和 C opt C opt  C_("opt ")C_{\text {opt }} 分别代表优化体积和结构应变能。