CN111709097B - 一种基于零亏格网格曲面连续变形的柔顺机构生成方法 - Google Patents

一种基于零亏格网格曲面连续变形的柔顺机构生成方法 Download PDF

Info

Publication number
CN111709097B
CN111709097B CN202010514409.0A CN202010514409A CN111709097B CN 111709097 B CN111709097 B CN 111709097B CN 202010514409 A CN202010514409 A CN 202010514409A CN 111709097 B CN111709097 B CN 111709097B
Authority
CN
China
Prior art keywords
mesh
isogeometric
matrix
iteration
unit
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010514409.0A
Other languages
English (en)
Other versions
CN111709097A (zh
Inventor
李宝童
丁森卯
苏文杰
程阿康
刘宏磊
洪军
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xian Jiaotong University
Original Assignee
Xian Jiaotong University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xian Jiaotong University filed Critical Xian Jiaotong University
Priority to CN202010514409.0A priority Critical patent/CN111709097B/zh
Publication of CN111709097A publication Critical patent/CN111709097A/zh
Application granted granted Critical
Publication of CN111709097B publication Critical patent/CN111709097B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/17—Mechanical parametric or variational design
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling

Landscapes

  • Physics & Mathematics (AREA)
  • Geometry (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Image Generation (AREA)

Abstract

一种基于零亏格网格曲面连续变形的柔顺机构生成方法，先配置初始组件布局，初始化基于等几何分析的计算层，在计算层的等几何单元控制点上施加载荷和约束，初始化移动渐近线法的迭代参数，优化迭代；计算三角形网格体积，将上层三角形网格描述的拓扑映射至下层计算层，计算控制点投影值，组装整体刚度矩阵；计算结构响应的目标函数值，求解体积的灵敏度，输入至优化器获得更新坐标值；根据更新制顶点坐标值，对网格细化、网格分裂和网格重划分操作的触发条件进行判定，完成网格再次更新，计算出非控制顶点更新坐标，更新三角形网格曲面形状，完成一次迭代；重复迭代直至满足收敛性条件；本发明可靠实现拓扑优化过程中的结构响应计算和敏度分析。

Description

一种基于零亏格网格曲面连续变形的柔顺机构生成方法
技术领域
本发明属于柔顺机构的设计技术领域，具体涉及一种基于零亏格网格曲面连续变形的柔顺机构生成方法。
背景技术
柔顺机构是一种整体的、连续的结构，能够通过自身弹性变形将运动，将力和能量从输入端口传递到输出端口，该机构与由刚性铰链和连杆组成的刚性机构存在很大不同；由于无需组装和润滑，没有摩擦，具有更长的使用寿命和更小的空间等优点，柔顺机构已广泛用于光纤对准、微机电系统、生物细胞操纵、精密工程等领域当中，以实现超高精度。
柔顺机构的设计方法主要包括两类，第一类是基于运动学的方法，通常利用伪刚体模型来设计由刚性连杆和柔性铰链构成的集中式柔顺机构；在该方法中，通过柔性铰链代替已知刚性连杆机构中的传统铰链来获得柔顺机构，因此可以直接将刚性机构的分析方法应用于柔顺机构中；这种便利性使得柔顺机构可以很容易地从现有的刚性机构发展获得；然而，基于运动学的方法高度依赖于已知的刚性连杆机构，这使得设计过程严重依赖于设计者的经验。
柔顺机构设计的另一类方法是拓扑优化方法，通常用于设计分布式柔顺机构；迄今为止，拓扑优化已经得到了全面发展，多种拓扑优化方法相继被提出；在柔顺机构的设计领域中，以带有惩罚的固体各向同性材料(SIMP)方法和水平集方法最常用。在SIMP方法中，设计域由合理分辨率的像素离散化，利用黑白像素的分布表示结构，由于其简单的概念和易用性，该方法已广泛应用于各种商业软件。然而SIMP方法由像素形成的边界缺乏明确的几何信息，甚至包含没有明确物理意义的灰色元素，因此无法准确地描述结构的几何形状和拓扑，这将在后处理中造成很大困难；其次，缺乏参数化几何信息会导致在拓扑优化过程中难以控制结构的特征尺寸，优化产生的柔性铰链会在柔顺机构工作时造成严重的应力集中；此外，在考虑与制造相关的要求(例如，增材制造中圆角的最小半径和结构部件的悬垂角)时，仍然缺乏明确的边界几何信息。水平集方法不同于SIMP方法，后者侧重于微观级别的材料密度，而水平集方法侧重边界演化；水平集方法仍将设计域离散化为有限元，并通过数值分析来计算结构响应，并且可以通过提取水平集函数的零轮廓来获得优化的拓扑。然而通过水平集函数获得的几何信息(例如边界曲率)同样不是参数化信息，因此在后处理过程中仍然存在一些SIMP中常见的问题，这是SIMP和水平集方法作为隐式方法存在的固有缺陷。
在柔顺机构的拓扑优化设计中，为了克服隐式方法存在的问题已经发展了多种显式拓扑优化方法；代表性的显式方法包括几何投影方法和进一步发展的移动可变形组件方法。在显式方法中，通常使用几何组件(图元)来描述拓扑，并通过直接控制组件形状参数来获得最佳拓扑；通过变形，移动，重叠和隐藏这些组件，可以描述各种组合的几何形状；显式拓扑优化框架被视为两层结构，其中上层用于几何描述，而由有限元离散的下层执行结构响应计算；在每次优化迭代中，将上层中的组合几何体投影到下层，投影边界内的有限元节点设置为实体材料，反之设置为弱材料，从而将所有有限元划分为实体，空隙和相交元素，并进一步利用插值算法计算相交元素的物理属性；在下层计算了结构响应和灵敏度后，将迭代信息反馈到上层以进一步控制下一次变形，并重复进行迭代以获得最佳结构。
在基于组件(图元)的显式拓扑优化方法中，柔顺机构的拓扑是由多个组件的变形、移动和组合形成的，因此需要足够数量的组件来满足拓扑中孔洞的形成，在没有足够组件的情况下显然无法形成最佳结构；因此优化结果高度依赖于组件的数量和初始布局，这种依赖性使得优化变得复杂；为了减少初始组件布局对柔顺机构最优结构的影响，有必要考虑如何仅使用简单的初始组件布局来有效，稳健地形成全局最优结构；此外，目前大多数拓扑优化方法都采用有限元进行设计域离散和结构响应计算，由于有限元本身的低稳定性和低阶连续，限制了拓扑优化最优解的获得，因此需要找到更可靠的分析方法。
发明内容
为了克服上述现有技术的缺点，本发明的目的在于提供了一种基于零亏格网格曲面连续变形的柔顺机构生成方法，使用可变三角形网格描述拓扑并投影到等几何单元进行结构响应计算和敏度分析，在克服初始布局依赖性的同时减少计算量，通过等几何分析更可靠地实现柔顺机构拓扑优化设计。
为了达到上述目的，本发明采取的技术方案为：
一种基于零亏格网格曲面连续变形的柔顺机构生成方法，包括以下步骤：
1)初始化零亏格网格曲面，配置初始组件布局；设置一系列顶点坐标并连接顶点形成几何，参数化几何进一步被划分为更小的三角形网格，形成网格曲面描述的初始组件布局；给定三角形网格细化操作和分裂操作的阈值LS和DS；设置材料用量上限V，设置最大迭代步数Loop以及针对变量变化的收敛标准tv；
2)初始化基于等几何分析的计算层；在计算层中使用等几何单元进行离散，利用等参变换计算等几何单元统一的单元刚度矩阵；
3)根据给定的边界条件在计算层的等几何单元控制点上施加载荷和约束；
4)迭代准备，初始化作为优化器的移动渐近线法的迭代参数；
5)启动优化迭代过程，根据上层三角形网格的半边数据结构，通过遍历所有三角单元获得三角形网格曲面的几何边界，将边界上的顶点用作控制顶点，即优化变量；
6)计算三角形网格体积V；
7)根据投影规则将上层三角形网格描述的拓扑映射至下层计算层，计算背景等几何单元上所有控制点的投影值；
8)根据控制点的投影值计算各等几何单元的单元刚度矩阵，边界通过的交叉单元利用插值算法求得，并依据所有单元刚度矩阵组装整体刚度矩阵；
9)利用等几何分析计算结构响应的目标函数值；
10)求解目标函数和体积的灵敏度；
11)将目标函数和体积的灵敏度输入至优化器中分析计算，获得优化变量更新后的坐标值；
12)根据更新后的制顶点坐标值，依次对网格细化、网格分裂和网格重划分操作的触发条件进行判定，满足触发条件则进行相应网格操作，完成网格再次更新；
13)根据网格操作后的控制顶点的坐标，采用保刚性算法计算出其他非控制顶点的更新坐标，以此更新三角形网格曲面的形状，完成一次迭代过程；
14)收敛性判断，当迭代步数k>Loop时，强制停止迭代过程，结束迭代设计；
15)收敛性判断，当相邻两个迭代步骤中变量的最大变化值小于tv时，即max(|ck +1-ck|)＜tv，满足收敛条件，迭代停止，优化设计过程完成；
16)返回步骤5)重复迭代过程，直至满足收敛性条件。
所述的步骤1)中利用三角形网格曲面显式地描述拓扑时，拓扑优化问题公式表达为：
Figure BDA0002529451390000051
式中，d＝(d1,...,dnc)T表示用作设计变量的所有控制顶点的坐标向量，nt表示优化问题中涉及的控制顶点的总数，I是优化问题的目标函数，Ωi，i＝1，2，...，n表示第i个三角形单元占据的区域，ε(·)表示线性应变矢量，Di代表第i个三角形单元的弹性矩阵，bi和
Figure BDA0002529451390000052
分别表示体力密度和在诺依曼边界上的表面牵引力，u和v分别表示位移场和相应的测试函数，
Figure BDA0002529451390000053
是狄利克雷边界上的给定位移，
Figure BDA0002529451390000054
表示体积V(d)的上限，Ud是d所属的可允许集合。
所述的步骤2)中等几何单元的单刚阵的计算方法如下：
Figure BDA0002529451390000055
式中，ξi和ηj表示第e号等几何单元ξ方向和η方向上的高斯积分点，m和n表示高斯积分点个数；雅克比矩阵J1和J2用于物理空间Ωe，参数空间
Figure BDA0002529451390000061
和母空间
Figure BDA0002529451390000062
之间的映射；B为应变位移矩阵，ωi,j表示等几何分析中的NURBS基函数，每个等几何单元的弹性矩阵D会随着三角形网格曲面的变形而变化。
所述的步骤6)中三角形网格体积V的计算方法如下：
Figure BDA0002529451390000063
式中ng是三角形单元的个数，lng是三角单元的边长；p的计算方法如下：
Figure BDA0002529451390000064
所述的步骤7)的投影过程中，对于背景等几何网格上任意一个控制点N，根据该控制点的坐标和三角单元投影ΔABC之间的位置关系，控制点的投影值计算如下：
Figure BDA0002529451390000065
式中，
Figure BDA0002529451390000066
以及
Figure BDA0002529451390000067
的具体计算方法如下：
Figure BDA0002529451390000068
式中，其中ng是三角形网格曲面中三角形单元的数量，e是欧拉数，b是一个很大的常数，此处取b＝100000。
所述的步骤8)中各等几何单元的单刚阵需要在步骤2)中所求的统一单刚阵基础上代入各单元的杨氏模量，底层等几何单元的杨氏模量的计算方法如下：
Figure BDA0002529451390000071
式中，ngp表示每个等几何单元中高斯点的个数，ωigp表示高斯点的权重，Ricp是NURBS基函数，E0是每个控制点处固体材料的杨氏模量，pe是惩罚系数，
Figure BDA0002529451390000072
icp＝1，2，...，ncp表示第e个单元中第icp个控制点的海维赛德函数值，
Figure BDA0002529451390000073
为控制点的投影值。
所述的步骤9)的目标函数值的具体计算方法如下：
Figure BDA0002529451390000074
式中，J*是人为设置的期望雅可比矩阵，wj和ws是权重因子；kin和kout是机构的输入和输出刚度；f和s为机构的运动要求和强度要求，目标函数将设计视为多目标问题。
所述的步骤10)中目标函数和体积的灵敏度的具体计算方法如下：
机构柔度矩阵的敏度：
Figure BDA0002529451390000075
式中，Cm为机构柔度矩阵，目标函数的敏度由
Figure BDA0002529451390000076
决定，L1和L2是常数矩阵，用于将K22和K21从机构刚度矩阵Km中提取出来；NE是等几何单元的个数，Uj(Uk)是虚拟单位力Fj(Fk)下的位移；
体积的灵敏度：
Figure BDA0002529451390000081
所述的步骤13)中已知控制顶点的更新坐标后，非控制顶点的移动根据最小化能量函数获得，能量函数的具体计算方法如下：
Figure BDA0002529451390000082
式中，M和M′表示变形前后的三角形网格，v′i和v′j表示vi和vj移动后的笛卡尔坐标，Ri是一个2×2的旋转矩阵；更新控制顶点坐标后需要求出未知的v′i，v′j以及Ri；首先每个顶点的平移设为零以求得近似的旋转矩阵，然后固定旋转矩阵计算新的位置v′i和v′j；因此问题转换为线性最小二乘问题，此时求解方程(14)的偏导数的线性方程组即获得新的变形。
本发明具有如下有益的技术效果：
1)本发明使用三角形网格显式地描述几何，解决了传统隐式方法中优化结果缺乏明确几何信息、优化过程无法有效控制柔顺机构特征尺寸的问题，同时只配置一个初始组件，通过网格变形和网格操作技术实现连续变拓扑，解决了一般显式方法中优化结果对初始布局存在严重依赖性的问题。
2)本发明采用等几何分析计算拓扑的结构响应和敏度值，相比较于传统的有限元方法能提高网格单元之间的连续性，从而提高计算精度，获得最佳柔顺机构拓扑。
3)本发明提出一种几何层和计算层之间的投影方法，通过控制点和投影三角形顶点之间的矢量积求得控制点的投影值，从而求出各等几何单元的单元刚度矩阵，完成从上层的几何到下层等几何单元的映射，为计算结构响应计算做准备。
4)本发明采用基于雅克比矩阵的方法构造柔顺机构的设计公式，将多目标优化问题转换为单目标优化问题，最小化目标函数时同时优化柔顺机构的雅克比矩阵和出入输出刚度，使设计结果同时满足运动要求和强度要求。
附图说明
图1是本发明的流程图。
图2是显式三角形网格曲面投影到等几何层的两层模型。
图3是柔顺机构物理模型。
图4是角点控制顶点驱动的网格变形。
图5是基于所提方法设计柔顺位移换向器的优化结果。
具体实施方式
下面结合附图对本发明做进一步详细描述。
如图1所示，一种基于零亏格网格曲面连续变形的柔顺机构生成方法，包括以下步骤：
1)初始化零亏格网格曲面，配置初始组件布局；设置一系列顶点坐标并连接顶点形成几何，参数化几何进一步被划分为更小的三角形网格，形成网格曲面描述的初始组件布局；给定三角形网格细化操作和分裂操作的阈值，分别为设计域最短边的25％和30％；设置材料用量上限为20％，设置最大迭代步数为1000以及针对变量变化的收敛标准为0.001；
利用三角形网格曲面进行拓扑描述时，柔顺机构设计问题的公式可以表达为：
Figure BDA0002529451390000101
式中，d＝(d1,...,dnc)T表示用作设计变量的所有控制顶点的坐标向量，nt表示优化问题中涉及的控制顶点的总数，I是优化问题的目标函数，Ωi，i＝1，2，...，n表示第i个三角形单元占据的区域，ε(·)表示线性应变矢量，Di代表第i个三角形单元的弹性矩阵，bi和
Figure BDA0002529451390000102
分别表示体力密度和在诺依曼边界上的表面牵引力，u和v分别表示位移场和相应的测试函数，
Figure BDA0002529451390000103
是狄利克雷边界上的给定位移，
Figure BDA0002529451390000104
表示体积V(d)的上限，Ud是d所属的可允许集合；
2)初始化基于等几何分析的计算层；在计算层中使用等几何单元代替传统有限元进行离散，利用等参变换计算等几何单元统一的单元刚度矩阵，此时等几何单元的单刚阵由公式(2)计算得到：
Figure BDA0002529451390000111
式中，ξi和ηj表示第e号等几何单元ξ方向和η方向上的高斯积分点，m和n表示高斯积分点个数；雅克比矩阵J1和J2用于物理空间Ωe，参数空间
Figure BDA0002529451390000112
和母空间
Figure BDA0002529451390000113
之间的映射；B为应变位移矩阵，ωi,j表示等几何分析中的NURBS基函数，每个等几何单元的弹性矩阵D会随着三角形网格曲面的变形而变化；
在此2D问题中，弹性矩阵D可以计算为：
Figure BDA0002529451390000114
3)根据给定的边界条件在计算层的等几何单元控制点上施加载荷和约束；
4)迭代准备，初始化作为优化器的移动渐近线法的迭代参数；
5)启动优化迭代过程，顶层三角形网格的几何信息存储在半边结构中，任意半边结构都存储了每个半边的起点、终点及其上一下、下一个和相对的半边，因此可以在不发生任何拓扑错误下遍历整个三角形网格，从而轻松识别三角形网格几何边界，并将几何边界上的顶点用作控制顶点，即优化变量；
6.计算三角形网格的体积V，具体计算公式如下：
Figure BDA0002529451390000115
式中ng是三角形单元的个数，lng是三角单元的边长；p由海伦公式(5)计算获得：
Figure BDA0002529451390000121
7)根据投影规则将上层的几何映射到底层，投影过程如图2所示，当控制点被上层几何覆盖时，设为实体材料，否则设为弱材料，通过投影值φ记录控制点状态；背景网格上控制点的投影值φ的计算公式如下：
Figure BDA0002529451390000122
式中，
Figure BDA0002529451390000123
以及
Figure BDA0002529451390000124
的具体计算方法如下：
Figure BDA0002529451390000125
式中，其中ng是三角形网格曲面中三角形单元的数量，e是欧拉数，b是一个很大的常数，此处取b＝100000；
8)根据控制点的投影值计算各等几何单元的单元刚度矩阵，下层等几何单元的物理场是由控制点插值得到，因此下层单元经过投影后的状态也由控制点的投影值φ决定，经过投影后各等几何单元的刚度矩阵需要在统一的单刚阵基础上乘投影系数，该投影系数是一个0-1的值，0表示该单元没有被覆盖，1表示该单元完全被覆盖，如图2所示，投影后的单元分为实体单元、弱单元以及边界单元；投影系数影响下各单元的杨氏模量计算如下：
Figure BDA0002529451390000126
式中，ngp表示每个等几何单元中高斯点的个数，ωigp表示高斯点的权重，Ricp是NURBS基函数，E0是每个控制点处固体材料的杨氏模量，pe是惩罚系数，
Figure BDA0002529451390000131
icp＝1，2，...，ncp表示第e个单元中第icp个控制点的海维赛德函数值，
Figure BDA0002529451390000132
为控制点的投影值；
获得各个单刚阵后，按照自由度标号组装整体刚度矩阵K；
9)利用等几何分析计算结构响应的目标函数值；图3为柔顺机构的物理模型，柔顺机构优化问题的目标函数值可以计算为：
Figure BDA0002529451390000133
式中，J*是人为设置的期望雅可比矩阵，wj和ws是权重因子；kin和kout是机构的输入和输出刚度；f和s为机构的运动要求和强度要求，目标函数将设计视为多目标问题；
在此优化设计过程中，设计目标为在各约束下最小化目标函数值；
10)求解目标函数和体积的灵敏度；本发明采用梯度法进行迭代优化，因此需要求解敏度值，目标函数和体积的灵敏度的具体计算方法如下：
雅克比矩阵的灵敏度：
Figure BDA0002529451390000134
输入刚度的灵敏度：
Figure BDA0002529451390000135
由于存在关系
Figure BDA0002529451390000136
Cm的灵敏度的计算如下：
Figure BDA0002529451390000141
式中，Cm为机构柔度矩阵，目标函数的敏度由
Figure BDA0002529451390000142
决定，L1和L2是常数矩阵，用于将K22和K21从机构刚度矩阵Km中提取出来；NE是等几何单元的个数，Uj(Uk)是虚拟单位力Fj(Fk)下的位移；
体积的灵敏度：
Figure BDA0002529451390000143
11)根据步骤10)计算的目标函数值和敏度值，将计算结果输入至优化器中进行分析计算，获得优化变量新的坐标值；
12)根据上述控制顶点更新后的坐标值进行网格操作判定；三角形网格单元边长Lmax＞LS时启动网格细分操作；相邻三角形网格单元之间未共享两顶点之间距离Dij＞DS时启动网格分裂操作；通过Greiner-Hormann算法进行自相交判定，采用ConstrainedDelaunayTriangulation算法进行网格重划分操作；执行网格操作后完成网格再次更新；
13)三角形网格的运动规则如图4所示，假定边界上的控制顶点发生移动后，内部的其他顶点需要按变形规则发生运动；因此，根据上述获得更新后的控制顶点坐标后，采用保刚性算法计算出其他非控制顶点的更新坐标，非控制顶点的移动可根据最小化能量函数获得，能量函数的具体计算方法如下：
Figure BDA0002529451390000151
式中，M和M′表示变形前后的三角形网格，v′i和v′j表示vi和vj移动后的笛卡尔坐标，Ri是一个2×2的旋转矩阵；更新控制顶点坐标后需要求出未知的v′i，v′j以及Ri；
首先每个顶点的平移设为零以求得近似的旋转矩阵，然后固定旋转矩阵计算新的位置v′i和v′j；因此问题转换为线性最小二乘问题，此时求解方程(14)的偏导数的线性方程组即可获得新的变形，以此更新三角形网格曲面的形状，完成一次迭代过程；
14)收敛性判断，当迭代步数k>1000时，强制停止迭代过程，结束迭代设计；
15.收敛性判断，当相邻迭代步骤中变量的最大变化值小于0.001时，即max(|ck+1-ck|)＜tv，满足收敛准则，迭代停止，优化设计过程完成；
16)返回步骤5)重复迭代过程，直至满足收敛性条件；图5为位移换向器的设计迭代设计结果，该拓扑由网格曲面显示表达，设计获得的机构满足柔顺机构的运动要求和强度要求。
综上所述，本发明提供的柔顺机构拓扑优化方法结合了新开发的几何变形技术与等几何分析，并提出了一种新的两层集成模型；在拓扑优化两层模型的上层，柔顺机构的拓扑由三角形网格曲面显式描述，采用半边结构存储和识别三角形网格几何信息；通过移动、分裂和细化这些三角形网格单元，将生成的形状投影至使用等几何单元离散的下层；在下层中使用等几何分析计算目标函数以及敏度，并将敏度值输入至优化器进行变量更新，并将更新信息反馈至上层，以驱动新的迭代；迭代至一定步数满足收敛条件后完成柔顺机构的整个设计过程，最终获得最优的柔顺机构拓扑。

Claims (9)

1.一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，包括以下步骤：
1)初始化零亏格网格曲面，配置初始组件布局；设置一系列顶点坐标并连接顶点形成几何，参数化几何进一步被划分为更小的三角形网格，形成网格曲面描述的初始组件布局；给定三角形网格细化操作和分裂操作的阈值LS和DS；设置材料用量上限
Figure FDA0002529451380000011
设置最大迭代步数Loop以及针对变量变化的收敛标准tv；
2)初始化基于等几何分析的计算层；在计算层中使用等几何单元进行离散，利用等参变换计算等几何单元统一的单元刚度矩阵；
3)根据给定的边界条件在计算层的等几何单元控制点上施加载荷和约束；
4)迭代准备，初始化作为优化器的移动渐近线法的迭代参数；
5)启动优化迭代过程，根据上层三角形网格的半边数据结构，通过遍历所有三角单元获得三角形网格曲面的几何边界，将边界上的顶点用作控制顶点，即优化变量；
6)计算三角形网格体积V；
7)根据投影规则将上层三角形网格描述的拓扑映射至下层计算层，计算背景等几何单元上所有控制点的投影值；
8)根据控制点的投影值计算各等几何单元的单元刚度矩阵，边界通过的交叉单元利用插值算法求得，并依据所有单元刚度矩阵组装整体刚度矩阵；
9)利用等几何分析计算结构响应的目标函数值；
10)求解目标函数和体积的灵敏度；
11)将目标函数和体积的灵敏度输入至优化器中分析计算，获得优化变量更新后的坐标值；
12)根据更新后的制顶点坐标值，依次对网格细化、网格分裂和网格重划分操作的触发条件进行判定，满足触发条件则进行相应网格操作，完成网格再次更新；
13)根据网格操作后的控制顶点的坐标，采用保刚性算法计算出其他非控制顶点的更新坐标，以此更新三角形网格曲面的形状，完成一次迭代过程；
14)收敛性判断，当迭代步数k>Loop时，强制停止迭代过程，结束迭代设计；
15)收敛性判断，当相邻两个迭代步骤中变量的最大变化值小于tv时，即max(|ck+1-ck|)＜tv，满足收敛条件，迭代停止，优化设计过程完成；
16)返回步骤5)重复迭代过程，直至满足收敛性条件。
2.根据权利要求1所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤1)中利用三角形网格曲面显式地描述拓扑时，拓扑优化问题公式表达为：
Figure FDA0002529451380000021
式中，d＝(d1,...,dnc)T表示用作设计变量的所有控制顶点的坐标向量，nt表示优化问题中涉及的控制顶点的总数，I是优化问题的目标函数，Ωi，i＝1，2，...，n表示第i个三角形单元占据的区域，ε(·)表示线性应变矢量，Di代表第i个三角形单元的弹性矩阵，bi和
Figure FDA0002529451380000038
分别表示体力密度和在诺依曼边界上的表面牵引力，u和v分别表示位移场和相应的测试函数，
Figure FDA0002529451380000031
是狄利克雷边界上的给定位移，
Figure FDA0002529451380000032
表示体积V(d)的上限，Ud是d所属的可允许集合。
3.根据权利要求2述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤2)中等几何单元的单刚阵的计算方法如下：
Figure FDA0002529451380000033
式中，ξi和ηj表示第e号等几何单元ξ方向和η方向上的高斯积分点，m和n表示高斯积分点个数；雅克比矩阵J1和J2用于物理空间Ωe，参数空间
Figure FDA0002529451380000034
和母空间
Figure FDA0002529451380000035
之间的映射；B为应变位移矩阵，ωi,j表示等几何分析中的NURBS基函数，每个等几何单元的弹性矩阵D会随着三角形网格曲面的变形而变化。
4.根据权利要求3所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤6)中三角形网格体积V的计算方法如下：
Figure FDA0002529451380000036
式中ng是三角形单元的个数，lng是三角单元的边长；p的计算方法如下：
Figure FDA0002529451380000037
5.根据权利要求4所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤7)的投影过程中，对于背景等几何网格上任意一个控制点N，根据该控制点的坐标和三角单元投影ΔABC之间的位置关系，控制点的投影值计算如下：
Figure FDA0002529451380000041
式中，
Figure FDA0002529451380000042
以及
Figure FDA0002529451380000043
的具体计算方法如下：
Figure FDA0002529451380000044
式中，其中ng是三角形网格曲面中三角形单元的数量，e是欧拉数，b是一个很大的常数，此处取b＝100000。
6.根据权利要求5所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤8)中各等几何单元的单刚阵需要在步骤2)中所求的统一单刚阵基础上代入各单元的杨氏模量，底层等几何单元的杨氏模量的计算方法如下：
Figure FDA0002529451380000045
式中，ngp表示每个等几何单元中高斯点的个数，ωigp表示高斯点的权重，Ricp是NURBS基函数，E0是每个控制点处固体材料的杨氏模量，pe是惩罚系数，
Figure FDA0002529451380000046
表示第e个单元中第icp个控制点的海维赛德函数值，
Figure FDA0002529451380000047
为控制点的投影值。
7.根据权利要求6所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤9)的目标函数值的具体计算方法如下：
Figure FDA0002529451380000051
式中，J*是人为设置的期望雅可比矩阵，wf和ws是权重因子；kin和kout是机构的输入和输出刚度；f和s为机构的运动要求和强度要求，目标函数将设计视为多目标问题。
8.根据权利要求7所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤10)中目标函数和体积的灵敏度的具体计算方法如下：
机构柔度矩阵的敏度：
Figure FDA0002529451380000052
式中，Cm为机构柔度矩阵，目标函数的敏度由
Figure FDA0002529451380000053
决定，L1和L2是常数矩阵，用于将K22和K21从机构刚度矩阵Km中提取出来；NE是等几何单元的个数，Uj(Uk)是虚拟单位力Fj(Fk)下的位移；
体积的灵敏度：
Figure FDA0002529451380000061
9.根据权利要求8所述的一种基于零亏格网格曲面连续变形的柔顺机构生成方法，其特征在于，所述的步骤13)中已知控制顶点的更新坐标后，非控制顶点的移动根据最小化能量函数获得，能量函数的具体计算方法如下：
Figure FDA0002529451380000062
式中，M和M′表示变形前后的三角形网格，v′i和v′j表示vi和vj移动后的笛卡尔坐标，Ri是一个2×2的旋转矩阵；更新控制顶点坐标后需要求出未知的v′i，v′j以及Ri；首先每个顶点的平移设为零以求得近似的旋转矩阵，然后固定旋转矩阵计算新的位置v′i和v′j；因此问题转换为线性最小二乘问题，此时求解方程14)的偏导数的线性方程组即获得新的变形。
CN202010514409.0A 2020-06-08 2020-06-08 一种基于零亏格网格曲面连续变形的柔顺机构生成方法 Active CN111709097B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010514409.0A CN111709097B (zh) 2020-06-08 2020-06-08 一种基于零亏格网格曲面连续变形的柔顺机构生成方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010514409.0A CN111709097B (zh) 2020-06-08 2020-06-08 一种基于零亏格网格曲面连续变形的柔顺机构生成方法

Publications (2)

Publication Number Publication Date
CN111709097A CN111709097A (zh) 2020-09-25
CN111709097B true CN111709097B (zh) 2022-12-09

Family

ID=72539385

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010514409.0A Active CN111709097B (zh) 2020-06-08 2020-06-08 一种基于零亏格网格曲面连续变形的柔顺机构生成方法

Country Status (1)

Country Link
CN (1) CN111709097B (zh)

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112084697B (zh) * 2020-10-26 2022-12-09 西安交通大学 一种可实现变形能传递解耦的柔顺机构设计方法
CN112949111B (zh) * 2021-01-29 2022-08-23 东风银轮(十堰)非金属部件有限公司 一种板簧设计优化方法
CN113591206B (zh) * 2021-07-09 2023-04-28 武汉理工大学 一种基于几何变形空间特征的船型优化设计方法及系统
CN115688497B (zh) * 2023-01-05 2023-05-05 安世亚太科技股份有限公司 载荷图元的显示方法、装置、电子设备及存储介质
CN116306088B (zh) * 2023-01-13 2024-02-06 华中科技大学 基于共形参数化的多尺度扫掠曲面结构设计方法及设备
CN117521234B (zh) * 2023-02-08 2024-03-12 中国空气动力研究与发展中心计算空气动力研究所 飞翼布局内外流一体化的参数化建模方法、设备及介质
CN116563490B (zh) * 2023-03-27 2024-01-23 浙江大学 一种基于刚性映射的全自动四边形附面层网格生成方法
CN116843656B (zh) * 2023-07-06 2024-03-15 安徽正汇汽配股份有限公司 钢带管的涂塑控制方法及其系统
CN116644620B (zh) * 2023-07-27 2023-10-17 华东交通大学 一种基于等几何分析的柔顺机构应力约束拓扑优化方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108897962A (zh) * 2018-07-05 2018-11-27 西北工业大学 一种基于非线性无网格法的柔性机构拓扑优化设计方法
CN109670200B (zh) * 2018-11-13 2022-04-22 华中科技大学 一种等几何材料密度场结构拓扑优化方法
CN109583091B (zh) * 2018-11-30 2022-11-29 长沙理工大学 基于自适应约束的柔性机构拓扑优化设计方法

Also Published As

Publication number Publication date
CN111709097A (zh) 2020-09-25

Similar Documents

Publication Publication Date Title
CN111709097B (zh) 一种基于零亏格网格曲面连续变形的柔顺机构生成方法
CN111709171B (zh) 一种热流强耦合问题的等几何求解及散热拓扑生成方法
JP4662959B2 (ja) デザイン最適化のための自由形状変形（ＦｒｅｅＦｏｒｍＤｅｆｏｒｍａｔｉｏｎ）表現の進化的直接操作（ＥｖｏｌｕｔｉｏｎａｒｙＤｉｒｅｃｔＭａｎｉｐｕｌａｔｉｏｎ）
CN106384384B (zh) 一种三维产品模型的形状优化方法
Kim et al. Eulerian shape design sensitivity analysis and optimization with a fixed grid
CN114741753B (zh) 一种薄壁加筋结构优化方法、装置、电子设备和存储介质
Zuo et al. Isogeometric analysis for CSG models
Sieger et al. A comprehensive comparison of shape deformation methods in evolutionary design optimization
CN115688276A (zh) 一种基于离散伴随方法的飞行器外形自动化优化方法、系统、设备、介质
Vucina et al. Computational procedure for optimum shape design based on chained Bezier surfaces parameterization
CN110689620B (zh) 一种多层次优化的网格曲面离散样条曲线设计方法
EP1840842A1 (en) Evolutionary design optimisation by means of extended direct manipulation of free form deformations
Schmidt et al. Extended isogeometric analysis of multi-material and multi-physics problems using hierarchical B-splines
Park et al. Boundary representation tolerance impacts on mesh generation and adaptation
CN116822038B (zh) 基于数据驱动的异型封闭加筋拓扑优化方法
CN115630542B (zh) 一种薄壁加筋结构的加筋布局优化方法
Ding et al. Isogeometric topology optimization of compliant mechanisms using transformable triangular mesh (TTM) algorithm
Schneider et al. Regularization and discretization error estimates for optimal control of ODEs with group sparsity
Chauhan et al. Wing shape optimization using FFD and twist parameterization
EP2442279B1 (en) Design optimization of physical bodies
Yang et al. Finite element mesh deformation with the skeleton-section template
Marco et al. Structural shape optimization using Cartesian grids and automatic h-adaptive mesh projection
Smith et al. A geometry projection method for the design exploration of wing-box structures
Ćurković et al. Novel re-parameterization for shape optimization and comparison with knot-based gradient fitting method
CN116049925A (zh) 面向建模-分析-优化一体化的组合薄壁元件结构自适应几何建模方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant