CN103020361B - 一种柔顺机构无棋盘格拓扑图提取方法 - Google Patents

Abstract

本发明公开了一种柔顺机构无棋盘格拓扑图提取方法，首先建立具有体积约束的柔顺机构拓扑优化模型；然后对拓扑图进行棋盘格分析，在此基础上，建立消除棋盘格的约束条件；之后，建立具有消除棋盘格的设计变量更新最佳准则；直至提取出柔顺机构无棋盘格拓扑图。本发明可有效克服基于SIMP方法松弛设计变量时所产生的棋盘格问题，提取出柔顺机构无棋盘格拓扑图；无需平滑设计域内的单元密度，能极大地改善拓扑优化结果中的中间单元现象；消除棋盘格的约束条件只施加在棋盘格区域的较低密度单元与其邻域内密度大于或等于给定阈值的单元之间，因而计算复杂度低。

Description

一种柔顺机构无棋盘格拓扑图提取方法

技术领域

本发明涉及柔顺机构拓扑优化设计中的拓扑图提取技术领域，特别涉及一种柔顺机构无棋盘格拓扑图提取方法。

背景技术

随着微机电系统、微加工和微操作以及新材料等领域的迅速发展，柔顺机构的设计已经成为目前国内外机构学领域的研究热点。

采用拓扑优化方法研究柔顺机构的设计只需给定设计域和指定输入输出位置，无须从一个已知的刚性机构出发，且所得的机构具有分布式柔性的优越性能，因而引起了极大的重视。这种方法通常是以有限元分析为基础，在拓扑优化的初始阶段，首先将设计域离散成一定数量的有限元网格，再利用优化方法确定单元材料的保留与删除，以满足预定的目标和约束条件。在优化结果中，单元密度的理想取值为0或1，当单元密度取值为0时，表示该单元被删除，单元密度取值为1时，该单元被保留，于是，由高密度单元连接构成机构的拓扑图。为了求解这种整数规划问题，人们常采用形如ρ^P,(P＞1)的SIMP方法松弛设计变量。

然而，这种方法容易导致棋盘格现象。它是指材料高低密度分布呈周期性交替的一种数值不稳定现象，给机构的拓扑图提取及后续的设计、加工都带来极大的困难。因此，在柔顺机构拓扑优化设计中，必须发展一种有效的方法，促使有材料的单元能够紧密连接，从而提取出具有材料连续分布属性的柔顺机构无棋盘格拓扑图。

目前，在以柔顺机构设计为代表拓扑优化领域中，主要包括以下几种解决棋盘格问题的方法。一种是高阶单元法，其实质是利用八节点或九节点等参元模型化设计问题。这种方法能在一定程度上改善棋盘格现象，但计算量太大，故不是很适用。

另一种常用的方法是施加几何约束，如周长法、局部密度斜率约束、密度过滤、灵敏度过滤、基于密度单调性变化的方法，等等。施加几何约束的方法的实质是平滑单元密度，因而有很好的消除棋盘格的效果，但是容易导致中间单元问题。

之后，又相继发展了基于小波的方法、非协调元和杂交元相结合的方法、拓扑分析和应力分析相结合的方法以及基于拓扑描述函数的方法等，这些方法也具有一定的消除棋盘格的效果，但由于它们仍需平滑单元密度，因而也无法避免中间单元的出现。

基于水平集的拓扑优化方法也是解决棋盘格问题的一种有效方法，但是该方法也具有初始敏感性、不能生成新孔、计算效率低和难以收敛到不光滑的角点等缺陷，虽然目前已提出一些改进的方法，但这些问题尚未完全得到很好的解决。

发明内容

本发明的发明目的是针对现有柔顺机构拓扑图提取方法的技术不足，提供一种柔顺机构无棋盘格拓扑图提取方法。

本发明首先建立具有体积约束的柔顺机构拓扑优化模型；然后对拓扑图进行棋盘格分析，在此基础上，建立消除棋盘格的约束条件；之后，建立具有消除棋盘格的设计变量更新最佳准则；最后再通过有限元分析、灵敏度计算及基于最佳准则的设计变量更新，获得一次迭代优化的拓扑图；重复消除棋盘格约束条件的建立、具有消除棋盘格的设计变量更新最佳准则的建立及一次迭代优化的拓扑图的获得过程，直至提取出柔顺机构无棋盘格拓扑图。

为实现上述发明目的，本发明采用的技术方案为：

提供一种柔顺机构无棋盘格拓扑图提取方法，其特征在于包括如下步骤：

1)：建立具有体积约束的柔顺机构拓扑优化模型；

2)：建立消除柔顺机构拓扑图中的棋盘格的约束条件，步骤如下：

2-1)：分析柔顺机构拓扑图中的棋盘格；

2-2)：施加消除棋盘格的约束条件；

3)建立具有消除棋盘格的设计变量更新最佳准则；

4)：通过一次迭代优化，得到优化后的拓扑图；

5)重复步骤2)至4)，直至达到最大迭代数或单元密度变化最大值小于阈值为止；提取到柔顺机构无棋盘格拓扑图。

优选地，步骤1)中，

以Ω为设计域，其为柔顺机构拓扑优化模型可利用材料域，P_i和P_o分别为柔顺机构拓扑优化模型载荷输入点和位移输出点，F_in和F_d分别为柔顺机构拓扑优化模型输入的载荷和沿输出位移方向的虚拟单位载荷，k_in和k_out分别为柔顺机构拓扑优化模型输入和输出弹簧刚度，柔顺机构的目标体积为V^*，将设计域离散成N个单元；柔顺机构拓扑优化模型的应变能和互应变能如下：

E_s＝∫_Ωε(u)^TDε(u)dΩ＝U^TKU

$E_{\mathrm{ms}}={\∫}_{\mathrm{\Ω}}\mathrm{\ϵ}{\left(u_d\right)}^T\mathrm{D\ϵ}\left(u\right)\mathrm{d\Ω}=U_d^T\mathrm{KU}$

式中，E_s是系统的应变能，应变能越小则表明系统的刚度越大；E_ms是系统的互应变能，互应变能越大则表明系统的柔性越大；D为弹性矩阵，K是系统整体刚度矩阵，U是F_in作用下的节点位移向量，U_d是F_d作用下的节点位移向量，ε(u)和u是设计域内任一点在载荷F作用下的应变和弹性变形，ε(u_d)和u_d是设计域内任一点在载荷F_d作用下的应变和弹性变形；

为使柔顺机构既有足够大的刚度又有足够大的柔性，通过多目标优化而得到柔顺机构的应变能和互应变能的关系如下：式中，符号Min代表最小值；

采用相对密度法松弛设计变量，使柔顺机构拓扑优化模型的单元密度可在0-1范围内取值，柔顺机构拓扑优化模型的单元密度如下：

0＜ρ_min≤ρ_i≤ρ_max＝1,i＝1,2,…,N；

式中，ρ_i是单元i的密度，ρ_min是单元密度下限，ρ_max是单元密度上限；

该柔顺机构拓扑优化模型的整体刚度矩阵如下：

$\begin{array}{c}K=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{K}_i\\ =\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\∫}_{V_e}{\mathrm{\ρ}}_i^P{B}^T\mathrm{DBdV}\\ =\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\ρ}}_i^P{K}^0\end{array};$

式中，K_i是单元i的刚度矩阵，V_e是任一实心单元的材料体积，P为密度ρ_i的指数，且P∈Z,P＞1，K⁰为任一实心单元的单元刚度矩阵，且B是任一实心单元的应变矩阵，该柔顺机构拓扑优化模型的体积约束如下：

$V\left(\mathrm{\ρ}\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{V}_e{\mathrm{\ρ}}_i\≤{\mathrm{\θ}}^{*}{V}_0;$

式中，V₀是柔顺机构拓扑优化模型可利用材料的体积，ρ是由ρ_i所构成的列向量，i＝1,2,…,N；

综合以上目标、整体刚度矩阵和体积约束，得到具有体积约束的柔顺机构拓扑优化模型，如下：

$\left\{\begin{array}{c}\begin{array}{cc}\mathrm{Min}:& f\left(\mathrm{\&rho;}\right)=-\frac{E_{\mathrm{ms}}}{E_s}\end{array}\\ \begin{array}{cc}s.t.& \mathrm{KU}=F\\ & {\mathrm{KU}}_d=F_d\\ & \underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{V}_e{\mathrm{\&rho;}}_i\&le;{\mathrm{\&theta;}}^{*}{V}_0\\ & 0<{\mathrm{\&rho;}}_{\min }\&le;{\mathrm{\&rho;}}_i\&le;{\mathrm{\&rho;}}_{\max }=1,i=\mathrm{1,2},...,N\end{array}.\end{array}\right.$

优选地，步骤2-1)中，步骤如下：

首先，选取[α₁,α₂]为棋盘格区域的低密度单元的密度搜索区间，且ρ_min＝α₁＜α₂＜1；

然后，在一次优化迭代初始阶段，先检测拓扑图中各单元的密度，当单元m的密度满足α₁≤ρ_m≤α₂时，再判断以下棋盘格生成条件是否成立：ρ_s-ρ_m≥α₃

式中，ρ_m是单元m的密度，ρ_s是与单元m具有公共边的四个邻接单元中的任一单元s的密度；参数α₃如下：α₃＝max(α₀-[10ρ_m]Δα,α_lim)

式中，符号[]表示取整运算，α₀为给定的初值，且ρ_min＜α₀≤1，Δα为增量，且0＜Δα＜α₀，α_lim为给定的α₃最小值，且ρ_min＜α_lim≤1；

如果棋盘格生成条件对单元m恒成立，则由单元m及其具有公共边的四个邻接单元构成棋盘格区域。

优选地，步骤2-2)中，步骤如下：

首先，以棋盘格区域的低密度单元m为中心，确定施加约束条件的邻域大小，并选取该邻域内单元密度阈值α₄,ρ_min＜α₄≤1，然后，施加消除棋盘格的约束条件，如下：

$\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})=\frac{1}{\underset{n\&Element;N_m}{\mathrm{\&Sigma;}}\frac{1}{r_{\mathrm{mn}}^2}}\underset{n\&Element;N_m}{\mathrm{\&Sigma;}}\frac{1}{r_{\mathrm{mn}}^2}{\mathrm{\&rho;}}_m^3{\mathrm{\&rho;}}_n^3\&GreaterEqual;G;$

式中，参数G为消除棋盘格约束的下限值，且0＜G≤1，N_m是施加约束条件的邻域内密度大于或等于α₄且不包括m在内的所有单元的集合，ρ_n为集合N_m中的单元n的密度，是集合N_m中的所有单元密度构成的列向量，r_mn是单元m和单元n之间的单元中心点距离。

优选地，步骤3)中，具有消除棋盘格的拉格朗日函数如下：

$\begin{array}{c}L=f+{\mathrm{\&lambda;}}_1(V\left(\mathrm{\&rho;}\right)-V^{*})+{\mathrm{\&lambda;}}_2^T(\mathrm{KU}-F_{\mathrm{in}})+{\mathrm{\&lambda;}}_3^T({\mathrm{KU}}_d-F_d)+\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_4^i({\mathrm{\&rho;}}_{\min }-{\mathrm{\&rho;}}_i+b_i^2)\\ +\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_5^i({\mathrm{\&rho;}}_i-{\mathrm{\&rho;}}_{\max }+c_i^2)+\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m(G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m}))\end{array};$

式中，λ₁，λ₂，λ₃，和β_m为拉格朗日乘子；V^*为柔顺机构的目标体积，和是松弛因子；

继续得到具有消除棋盘格的Kuhn-Tuck必要条件如下：

$\left\{\begin{array}{c}\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}+{\mathrm{\&lambda;}}_1{V}_i+{\mathrm{\&lambda;}}_2^T\frac{\&PartialD;\left(\mathrm{KU}\right)}{{\&PartialD;\mathrm{\&rho;}}_i}+{\mathrm{\&lambda;}}_3^T\frac{\&PartialD;\left({\mathrm{KU}}_d\right)}{{\&PartialD;\mathrm{\&rho;}}_i}-\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m\frac{\&PartialD;\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})}{{\&PartialD;\mathrm{\&rho;}}_i}\left\{\begin{array}{ccc}=0& \mathrm{if}& {\mathrm{\&rho;}}_{\min }<{\mathrm{\&rho;}}_i<{\mathrm{\&rho;}}_{\max }\\ >0& \mathrm{if}& {\mathrm{\&rho;}}_i={\mathrm{\&rho;}}_{\min }\\ <0& \mathrm{if}& {\mathrm{\&rho;}}_i={\mathrm{\&rho;}}_{\max }\end{array}\right.\\ V\left(\mathrm{\&rho;}\right)=V^{*}\\ \mathrm{KU}-F=0\\ {\mathrm{KU}}_d-F_d=0\\ G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})=0,m\&Element;N_{\mathrm{ck}}\\ i=\mathrm{1,2},...,N\end{array}\right..$

式中，N_ck是本次优化迭代中所有对棋盘格生成条件恒成立的低密度单元的集合；

首先，得到优化目标、体积约束和消除棋盘格的灵敏度分别如下：

优化目标的灵敏度： $\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}=\frac{E_s\left(P{\left({\mathrm{\&rho;}}_i\right)}^{P-1}{\left(u_i\right)}^T{K}^0{u}_{\mathrm{di}}\right)-E_{\mathrm{ms}}\left(P{\left({\mathrm{\&rho;}}_i\right)}^{P-1}{\left(u_i\right)}^T{K}^0{u}_i\right)}{{\left(E_s\right)}^2},i=1,...,N;$

体积约束的灵敏度： $\frac{\&PartialD;V}{{\&PartialD;\mathrm{\&rho;}}_i}=V_e,i=1,...,N;$

消除棋盘格约束的灵敏度：

然后，代入Kuhn-Tuck必要条件，得到灵敏度关系如下：

$Q_i^t=\frac{-\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}+\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m\frac{\&PartialD;\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})}{{\&PartialD;\mathrm{\&rho;}}_i}}{{\mathrm{\&lambda;}}_1\frac{\&PartialD;V}{{\&PartialD;\mathrm{\&rho;}}_i}}=1;$

于是，得到具有消除棋盘格的设计变量更新最佳准则，如下：

${\mathrm{\&rho;}}_i^{t+1}=\left\{\begin{array}{ccc}{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}& \mathrm{if}& \max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })<{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}<\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })\\ \min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })& \mathrm{if}& {\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}\&GreaterEqual;\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })\\ \max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })& \mathrm{if}& {\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}\&le;\max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })\end{array}\right.,i=\mathrm{1,2},...,N;$

式中，t为迭代数,η为松弛因子,且0＜η＜1,ζ为一个较小的移动极限,M_i是一个非负数，

$M_i^t=\max (0,Q_i^t);$

由于是λ₁和β_m的函数，i＝1,2,…,N，因此，在设计变量更新过程中，须采用二分法确定λ₁和β_m的值，以使更新后的密度满足体积约束和消除棋盘格约束，即：

${\mathrm{\&theta;}}^{*}{V}_0-\{\underset{\mathrm{icas}1}{\mathrm{\&Sigma;}}{V}_e{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}+\underset{\mathrm{icas}2}{\mathrm{\&Sigma;}}{V}_e\min \{(1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max }\}+\underset{\mathrm{icas}3}{\mathrm{\&Sigma;}}{V}_e\max \{(1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min }\left\}\right\}=0;$

$\begin{array}{c}G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m^t,{\mathrm{\&rho;}}_{N_m}^t)-\underset{\mathrm{icas}1}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t({\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}-{\mathrm{\&rho;}}_i^t)-\underset{\mathrm{icas}2}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t(\max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })-{\mathrm{\&rho;}}_i^t)\\ -\underset{\mathrm{icas}3}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t(\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })-{\mathrm{\&rho;}}_i^t)=0,m\&Element;N_{\mathrm{ck}}\end{array};$

式中，和分别是基于最佳准则所更新的三类设计变量的和。

优选地，步骤4)中，首先对柔顺机构拓扑优化模型进行有限元分析，并得出优化目标、体积约束和消除棋盘格约束的灵敏度；然后计算所有消除棋盘格约束值；最后基于具有消除棋盘格的最佳准则更新设计变量，得到一次迭代优化的拓扑图。

本发明相对于现有技术，具有以下有益效果：

1、本发明可有效克服基于SIMP方法的柔顺机构拓扑优化设计中所产生的棋盘格问题，提取出柔顺机构无棋盘格拓扑图。

2.本发明无需平滑设计域内的单元密度，因而能极大地改善拓扑优化结果中的中间单元现象。

3.本发明中的消除棋盘格的约束条件只施加在棋盘格区域的低密度单元与其邻域内密度大于或等于给定阈值的单元之间，因而计算复杂度低。

4.本发明具有较好的收敛稳定性。

附图说明

图1为本发明所提出的一种柔顺机构无棋盘格拓扑图提取方法流程图；

图2为本发明采用的一种典型的柔顺力-位移反向机构设计域和边界条件示意图；

图3为棋盘格现象示意图；

图4为α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.5，网格划分为60×30，消除棋盘格约束的邻域大小为3×3时，实施例柔顺机构无棋盘格拓扑图提取结果；

图5为α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.8，网格划分为60×30，消除棋盘格约束的邻域大小为3×3时，实施例柔顺机构无棋盘格拓扑图提取结果；

图6为α₀＝0.8，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.5，网格划分为60×30，消除棋盘格约束的邻域大小为3×3时，实施例柔顺机构无棋盘格拓扑图提取结果；

图7为α₀＝0.6，α₂＝0.9，α₄＝0.3，Δα＝0.1，α_lim＝0.3，G＝0.5，网格划分为60×30，消除棋盘格约束的邻域大小为3×3时，实施例柔顺机构无棋盘格拓扑图提取结果；

图8为α₀＝0.6，α₂＝0.9，α₄＝0.3，Δα＝0.1，α_lim＝0.5，G＝0.5，网格划分为60×30，消除棋盘格约束的邻域大小为3×3时，实施例柔顺机构无棋盘格拓扑图提取结果；

图9为α₀＝0.6，α₂＝0.9，α₄＝0.2，Δα＝0.1，α_lim＝0.3，G＝0.5，网格划分为60×30，消除棋盘格约束的邻域大小为3×3时，实施例柔顺机构无棋盘格拓扑图提取结果；

图10为α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.5，网格划分为70×35，消除棋盘格约束的邻域大小为6×6时，实施例柔顺机构无棋盘格拓扑图提取结果；

图11为α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝1，网格划分为70×35，消除棋盘格约束的邻域大小为5.1×5.1时，实施例柔顺机构无棋盘格拓扑图提取结果；

图12为α₀＝0.7，α₂＝0.5，α₄＝0.35，Δα＝0.18，α_lim＝0.3，G＝0.5，网格划分为70×35，消除棋盘格约束的邻域大小为4.3×4.3时，实施例柔顺机构无棋盘格拓扑图提取结果。

具体实施方式

下面结合附图和具体实施例对本发明的发明目的作进一步详细地描述，实施例不能在此一一赘述，但本发明的实施方式并不因此限定于以下实施例。除非特别说明，本发明采用的材料和加工方法为本技术领域常规材料和加工方法。

如图1所示，为本发明所提出的一种柔顺机构无棋盘格拓扑图提取方法流程图。

本发明实施例是一种典型的柔顺力-位移反向机构，其设计域和边界条件如图2所示。其中，设计域大小为60×60，材料的弹性模量和泊松比分别为E＝1和v＝0.3，P_i和P_o分别为载荷输入点和位移输出点，输入载荷F_in＝10，同时，沿输出位移方向作用有一虚拟单位载荷F_d，输入和输出弹簧刚度分别为k_in＝0.1和k_out＝0.1，设计问题的目标体积比为θ^*＝0.3。由于设计域和边界条件具有对称性，因此仅取设计域的下半区域进行优化。首先建立具有体积约束的柔顺机构拓扑优化模型为：

$\left\{\begin{array}{c}\begin{array}{cc}\mathrm{Min}:& f\left(\mathrm{\&rho;}\right)=-\frac{E_{\mathrm{ms}}}{E_s}\end{array}\\ \begin{array}{cc}s.t.& \mathrm{KU}=F\\ & {\mathrm{KU}}_d=F_d\\ & \underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{V}_e{\mathrm{\&rho;}}_i\&le;{\mathrm{\&theta;}}^{*}{V}_0\\ & 0<{\mathrm{\&rho;}}_{\min }\&le;{\mathrm{\&rho;}}_i\&le;{\mathrm{\&rho;}}_{\max }=1,i=\mathrm{1,2},...,N\end{array};\end{array}\right.$

式中，V_e是任一实心单元的材料体积,V₀是柔顺机构拓扑优化模型可利用材料的体积，ρ_min是单元密度下限，ρ_max是单元密度上限，N是设计域内的单元数，ρ_i是单元i的密度，U是载荷F_in作用下的节点位移向量，U_d是载荷F_d作用下的节点位移向量，ρ是由ρ_i所构成的列向量，i＝1,2,…,N；

，K是系统的整体刚度矩阵，公式如下：

$\begin{array}{c}K=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{K}_i\\ =\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\&Integral;}_{V_e}{\mathrm{\&rho;}}_i^P{B}^T\mathrm{DBdV}\\ =\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&rho;}}_i^P{K}^0\end{array};$

式中，K_i是单元i的刚度矩阵，K⁰为任一实心单元的单元刚度矩阵，且B是任一实心单元的应变矩阵，P为密度ρ_i的指数，且P∈Z,P＞1，本实施例中，为了加快收敛，取P＝3；

此外，E_s和E_ms分别是是系统的应变能和互应变能，公式如下：

E_s＝∫_Ωε(u)^TDε(u)dΩ＝U^TKU；

$E_{\mathrm{ms}}={\&Integral;}_{\mathrm{\&Omega;}}\mathrm{\&epsiv;}{\left(u_d\right)}^T\mathrm{D\&epsiv;}\left(u\right)\mathrm{d\&Omega;}=U_d^T\mathrm{KU};$

式中，D为弹性矩阵，U是F_in作用下的节点位移向量，U_d是F_d作用下的节点位移向量，ε(u)是设计域内任一点在载荷F作用下的应变，ε(u_d)是设计域内任一点在载荷F_d作用下的应变；

然后，选取[α₁,α₂]为棋盘格区域的低密度单元的密度搜索区间，且ρ_min＝α₁＜α₂＜1。

之后，在一次优化迭代初始阶段，先检测拓扑图中各单元的密度，当单元m的密度满足α₁≤ρ_m≤α₂时，再判断以下棋盘格生成条件是否成立：ρ_s-ρ_m≥α₃；

式中，ρ_m是单元m的密度，ρ_s是与单元m具有公共边的四个邻接单元中的任一单元s的密度。参数α₃的计算公式如下：

α₃＝max(α₀-[10ρ_m]Δα,α_lim)

式中，max代表最大值，符号[]表示取整运算，α₀为给定的初值，且ρ_min＜α₀≤1，Δα为增量，且0＜Δα＜α₀，α_lim为给定的α₃最小值，且ρ_min＜α_lim≤1。

如果棋盘格生成条件对单元m恒成立，则由单元m及其具有公共边的四个邻接单元构成棋盘格区域，如图3所示。

进一步地，以棋盘格区域的低密度单元m为中心，确定施加约束条件的邻域大小，并选取该邻域内单元密度阈值α₄,ρ_min＜α₄≤1，然后，施加消除棋盘格的约束条件，如下：

$\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})=\frac{1}{\underset{n\&Element;N_m}{\mathrm{\&Sigma;}}\frac{1}{r_{\mathrm{mn}}^2}}\underset{n\&Element;N_m}{\mathrm{\&Sigma;}}\frac{1}{r_{\mathrm{mn}}^2}{\mathrm{\&rho;}}_m^3{\mathrm{\&rho;}}_n^3\&GreaterEqual;G;$

式中，参数G为消除棋盘格约束的下限值，且0＜G≤1，N_m是施加约束条件的邻域内密度大于或等于α₄且不包括m在内的所有单元的集合，ρ_n为集合N_m中的单元n的密度，ρ_nei是集合N_m中的所有单元密度构成的列向量，r_mn是单元m和单元n之间的单元中心点距离。

因此，具有消除棋盘格的柔顺机构拓扑优化模型的拉格朗日函数是：

$\begin{array}{c}L=f+{\mathrm{\&lambda;}}_1(V\left(\mathrm{\&rho;}\right)-V^{*})+{\mathrm{\&lambda;}}_2^T(\mathrm{KU}-F_{\mathrm{in}})+{\mathrm{\&lambda;}}_3^T({\mathrm{KU}}_d-F_d)+\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_4^i({\mathrm{\&rho;}}_{\min }-{\mathrm{\&rho;}}_i+b_i^2)\\ +\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_5^i({\mathrm{\&rho;}}_i-{\mathrm{\&rho;}}_{\max }+c_i^2)+\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m(G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m}))\end{array};$

式中，λ₁，λ₂，λ₃，和β_m为拉格朗日乘子，V^*为柔顺机构的目标体积，和是松弛因子。继续得到具有消除棋盘格的Kuhn-Tucker必要条件是：

$\left\{\begin{array}{c}\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}+{\mathrm{\&lambda;}}_1{V}_i+{\mathrm{\&lambda;}}_2^T\frac{\&PartialD;\left(\mathrm{KU}\right)}{{\&PartialD;\mathrm{\&rho;}}_i}+{\mathrm{\&lambda;}}_3^T\frac{\&PartialD;\left({\mathrm{KU}}_d\right)}{{\&PartialD;\mathrm{\&rho;}}_i}-\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m\frac{\&PartialD;\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})}{{\&PartialD;\mathrm{\&rho;}}_i}\left\{\begin{array}{ccc}=0& \mathrm{if}& {\mathrm{\&rho;}}_{\min }<{\mathrm{\&rho;}}_i<{\mathrm{\&rho;}}_{\max }\\ >0& \mathrm{if}& {\mathrm{\&rho;}}_i={\mathrm{\&rho;}}_{\min }\\ <0& \mathrm{if}& {\mathrm{\&rho;}}_i={\mathrm{\&rho;}}_{\max }\end{array}\right.\\ V\left(\mathrm{\&rho;}\right)=V^{*}\\ \mathrm{KU}-F=0\\ {\mathrm{KU}}_d-F_d=0\\ G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})=0,m\&Element;N_{\mathrm{ck}}\\ i=\mathrm{1,2},...,N\end{array}\right..$

式中，N_ck是本次优化迭代中所有对棋盘格生成条件恒成立的低密度单元的集合；

进一步地，优化目标、体积约束和消除棋盘格约束的灵敏度计算公式是：

$\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}=\frac{E_s\left(P{\left({\mathrm{\&rho;}}_i\right)}^{P-1}{\left(u_i\right)}^T{K}^0{u}_{\mathrm{di}}\right)-E_{\mathrm{ms}}\left(P{\left({\mathrm{\&rho;}}_i\right)}^{P-1}{\left(u_i\right)}^T{K}^0{u}_i\right)}{{\left(E_s\right)}^2},i=1,...,N;$

$\frac{\&PartialD;V}{{\&PartialD;\mathrm{\&rho;}}_i}=V_e,i=1,...,N;$

然后，将以上灵敏度公式代入Kuhn-Tuck必要条件，得到公式：

$Q_i^t=\frac{-\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}+\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m\frac{\&PartialD;\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})}{{\&PartialD;\mathrm{\&rho;}}_i}}{{\mathrm{\&lambda;}}_1\frac{\&PartialD;V}{{\&PartialD;\mathrm{\&rho;}}_i}}=1$

因此，得到具有消除棋盘格的设计变量更新最佳准则，如下：

${\mathrm{\&rho;}}_i^{t+1}=\left\{\begin{array}{ccc}{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}& \mathrm{if}& \max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })<{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}<\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })\\ \min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })& \mathrm{if}& {\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}\&GreaterEqual;\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })\\ \max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })& \mathrm{if}& {\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}\&le;\max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })\end{array}\right.,i=\mathrm{1,2},...,N;$

式中，t为迭代数,η为松弛因子,且0＜η＜1,ζ为一个较小的移动极限，在本实施例中，为了保证密度更新过程的稳定性，η和ζ分别取值为0.3和0.05。

为一个非负数，即：

由于是λ₁和β_m的函数，i＝1,2,…,N，因此，在设计变量更新过程中，须采用二分法确定λ₁和β_m的值，以使更新后的密度满足体积约束和消除棋盘格约束，即：

${\mathrm{\&theta;}}^{*}{V}_0-\{\underset{\mathrm{icas}1}{\mathrm{\&Sigma;}}{V}_i{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}+\underset{\mathrm{icas}2}{\mathrm{\&Sigma;}}{V}_i\min \{(1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max }\}+\underset{\mathrm{icas}3}{\mathrm{\&Sigma;}}{V}_i\max \{(1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min }\left\}\right\}=0;$

$\begin{array}{c}G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m^t,{\mathrm{\&rho;}}_{N_m}^t)-\underset{\mathrm{icas}1}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t({\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}-{\mathrm{\&rho;}}_i^t)-\underset{\mathrm{icas}2}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t(\max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })-{\mathrm{\&rho;}}_i^t)\\ -\underset{\mathrm{icas}3}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t(\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })-{\mathrm{\&rho;}}_i^t)=0,m\&Element;N_{\mathrm{ck}}\end{array};$

式中，和分别是基于最佳准则所更新的三类设计变量的和。

之后，对柔顺机构拓扑优化模型进行有限元分析，得出优化目标、体积约束和消除棋盘格约束的灵敏度；然后计算所有消除棋盘格约束值；最后基于具有消除棋盘格的最佳准则更新设计变量，得到一次迭代优化的拓扑图。

其中，基于具有消除棋盘格的最佳准则更新设计变量的具体过程是：

以具有一个消除棋盘格约束的情形为例,首先给定β_m的最小值和最大值，分别是β_m1＝0和β_m2＝100000，计算其中点值：β_m＝0.5(β_m1+β_m2)；

然后保持β_m不变，确定λ₁的值并更新设计变量，其步骤是：

给定λ₁的最小值和最大值，分别是λ₁₁＝0和λ₁₂＝100000；然后计算中点值：λ₁＝0.5(λ₁₁+λ₁₂)；将β_m和λ₁代入设计变量更新公式，得到更新后的单元密度；计算体积约束余量再根据值修改λ₁₁或λ₁₂的值：如果则使λ₁₁＝λ₁，否则如果则使λ₁₂＝λ₁；重复λ₁₁和λ₁₂的中点值的计算、设计变量更新以及λ₁₁或λ₁₂的修改过程，直至λ₁₂与λ₁₁的差小于一个阈值为止。

再计算消除棋盘格约束的余量： $\&dtri;\mathrm{\&Phi;}=G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m^t,{\mathrm{\&rho;}}_{N_m}^t)-\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t({\mathrm{\&rho;}}_i^{t+1}-{\mathrm{\&rho;}}_i^t);$ 并根据值修改β_m1和β_m2的值：如果则使β_m1＝β_m，否则如果则使β_m2＝β_m，计算中点值β_m。

重复λ₁的确定、设计变量更新以及β_m1和β_m2的修改及中点值计算过程，直至β_m2与β_m1的差小于一个阈值为止。

重复消除棋盘格约束条件的建立、具有消除棋盘格的设计变量更新最佳准则的建立及一次迭代的最佳拓扑结果获得过程，直至达到最大迭代数或前后的单元密度变化最大值小于一个阈值为止，从而提取出柔顺机构无棋盘格拓扑图。

在本实施例中，为了验证不同的参数设置对本发明所提出的方法的性能的影响，对参数α₀、α₂、α₄、Δα、α_lim和G选取了九组数值，分别是：

(1)α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.5；

(2)α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.8；

(3)α₀＝0.8，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.5；

(4)α₀＝0.6，α₂＝0.9，α₄＝0.3，Δα＝0.1，α_lim＝0.3，G＝0.5；

(5)α₀＝0.6，α₂＝0.9，α₄＝0.3，Δα＝0.1，α_lim＝0.5，G＝0.5；

(6)α₀＝0.6，α₂＝0.9，α₄＝0.2，Δα＝0.1，α_lim＝0.3，G＝0.5；

(7)α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝0.5；

(8)α₀＝0.7，α₂＝0.5，α₄＝0.3，Δα＝0.2，α_lim＝0.3，G＝1；

(9)α₀＝0.7，α₂＝0.5，α₄＝0.35，Δα＝0.18，α_lim＝0.3，G＝0.5；

其中，(1)至(6)组参数取值对应的网格离散方式为60×30，消除棋盘格约束的邻域大小为3×3；

(7)至(10)组参数取值对应的网格离散方式为70×35，消除棋盘格约束的邻域大小分别为6×6、5.1×5.1和4.3×4.3。

图4至图12分别为对应于以上九组参数取值的实施例柔顺机构无棋盘格拓扑图提取结果。这些结果表明本发明所提出的柔顺机构无棋盘格拓扑图提取方法具有较好的收敛稳定性。

上述实施例仅为本发明的较佳实施例，并非用来限定本发明的实施范围。即凡依本发明内容所作的均等变化与修饰，都为本发明权利要求所要求保护的范围所涵盖。

Claims (2)

1.一种柔顺机构无棋盘格拓扑图提取方法，其特征在于包括如下步骤：

1)：建立具有体积约束的柔顺机构拓扑优化模型；

2)：建立消除柔顺机构拓扑图中的棋盘格的约束条件，步骤如下：

2-1)：分析柔顺机构拓扑图中的棋盘格；

2-2)：施加消除棋盘格的约束条件；

3)建立具有消除棋盘格的设计变量更新最佳准则；

4)：通过一次迭代优化，得到优化后的拓扑图；

5)重复步骤2)至4)，直至达到最大迭代数或单元密度变化最大值小于阈值为止；提取到柔顺机构无棋盘格拓扑图；

步骤1)中，

以Ω为设计域，其为柔顺机构拓扑优化模型可利用材料域，P_i和P_o分别为柔顺机构拓扑优化模型载荷输入点和位移输出点，F_in和F_d分别为柔顺机构拓扑优化模型的输入载荷和沿输出位移方向的虚拟单位载荷，k_in和k_out分别为柔顺机构拓扑优化模型输入和输出弹簧刚度，柔顺机构的目标体积比为θ^*，将设计域离散成N个单元；柔顺机构拓扑优化模型的应变能和互应变能如下：

E_s＝∫_Ωε(u)^TDε(u)dΩ_＝U^TKU

$E_{\mathrm{ms}}={\&Integral;}_{\mathrm{\&Omega;}}\mathrm{\&epsiv;}{\left(u_d\right)}^T\mathrm{D\&epsiv;}\left(u\right)\mathrm{d\&Omega;}=U_d^T\mathrm{KU}$

式中，E_s是系统的应变能，应变能越小则表明系统的刚度越大；E_ms是系统的互应变能，互应变能越大则表明系统的柔性越大；D为弹性矩阵，K是系统整体刚度矩阵，U是F_in作用下的节点位移向量，U_d是F_d作用下的节点位移向量，ε(u)和u是设计域内任一点在载荷F作用下的应变和弹性变形，ε(u_d)和u_d是设计域内任一点在载荷F_d作用下的应变和弹性变形；

为使柔顺机构既有足够大的刚度又有足够大的柔性，通过多目标优化而得到柔顺机构的应变能和互应变能的关系如下：式中，符号Min代表最小值；

采用相对密度法松弛设计变量，使柔顺机构拓扑优化模型的单元密度可在0-1范围内取值，柔顺机构拓扑优化模型的单元密度如下：

0＜ρ_min≤ρ_i≤ρ_max＝1,i＝1,2,…,N；

式中，ρ_i是单元i的密度，ρ_min是单元密度下限，ρ_max是单元密度上限；

该柔顺机构拓扑优化模型的整体刚度矩阵如下：

$\begin{array}{c}K=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{K}_i\\ =\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\&Integral;}_{V_e}{\mathrm{\&rho;}}_i^P{B}^T\mathrm{DBdV}\\ =\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&rho;}}_i^P{K}^0\end{array};$

式中，K_i是单元i的刚度矩阵，V_e是任一实心单元的材料体积，P为密度ρ_i的指数，且P∈Z,P＞1，K⁰为任一实心单元的单元刚度矩阵，且B是任一实心单元的应变矩阵，该柔顺机构拓扑优化模型的体积约束如下：

$V\left(\mathrm{\&rho;}\right)=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{V}_e{\mathrm{\&rho;}}_i\&le;{\mathrm{\&theta;}}^{*}{V}_0;$

式中，V₀是柔顺机构拓扑优化模型可利用材料的体积，ρ是由ρ_i所构成的列向量，i＝1,2,…,N_；

综合以上目标、整体刚度矩阵和体积约束，得到具有体积约束的柔顺机构拓扑优化模型，如下：

$\left\{\begin{array}{cc}\mathrm{Min}:& f\left(\mathrm{\&rho;}\right)=-\frac{E_{\mathrm{ms}}}{E_s}\\ s.t.& \mathrm{KU}=F\\ & {\mathrm{KU}}_d=F_d\\ & \underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{V}_e{\mathrm{\&rho;}}_i\&le;{\mathrm{\&theta;}}^{*}{V}_0\\ & 0<{\mathrm{\&rho;}}_{\min }\&le;{\mathrm{\&rho;}}_i\&le;{\mathrm{\&rho;}}_{\max }=1,i=\mathrm{1,2},...,N\end{array}\right.;$

步骤2-1)中，步骤如下：

首先，选取[α₁,α₂]为棋盘格区域的低密度单元的密度搜索区间，且ρ_min＝α₁＜α₂＜1_；

然后，在一次优化迭代初始阶段，先检测拓扑图中各单元的密度，当单元m的密度满足α₁≤ρ_m≤α₂时，再判断以下棋盘格生成条件是否成立：ρ_s-ρ_m≥α₃；

式中，ρ_m是单元m的密度，ρ_s是与单元m具有公共边的四个邻接单元中的任一单元s的密度；参数α₃如下：α₃＝max(α₀-[10ρ_m]Δα,α_lim)

式中，max代表最大值，符号[]表示取整运算，α₀为给定的初值，且ρ_min＜α₀≤1，Δα为增量，且0＜Δα＜α₀，α_lim为给定的α₃最小值，且ρ_min＜α_lim≤1；

如果棋盘格生成条件对单元m恒成立，则由单元m及其具有公共边的四个邻接单元构成棋盘格区域；

步骤2-2)中，步骤如下：

首先，以棋盘格区域的低密度单元m为中心，确定施加约束条件的邻域大小，并选取该邻域内单元密度阈值α₄,ρ_min＜α₄≤1，然后，施加消除棋盘格的约束条件，如下：

$\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})=\frac{1}{\underset{n\&Element;N_m}{\mathrm{\&Sigma;}}\frac{1}{r_{\mathrm{mn}}^2}}\underset{n\&Element;N_m}{\mathrm{\&Sigma;}}\frac{1}{r_{\mathrm{mn}}^2}{\mathrm{\&rho;}}_m^3{\mathrm{\&rho;}}_n^3\&GreaterEqual;G;$

式中，参数G为消除棋盘格约束的下限值，且0＜G≤1，N_m是施加约束条件的邻域内密度大于或等于α₄且不包括m在内的所有单元的集合，ρ_n为集合N_m中的单元n的密度，是集合N_m中的所有单元密度构成的列向量，r_mn是单元m和单元n之间的单元中心点距离；

步骤3)中，具有消除棋盘格的拉格朗日函数如下：

$\begin{array}{c}L=f+{\mathrm{\&lambda;}}_1(V\left(\mathrm{\&rho;}\right)-V^{*})+{\mathrm{\&lambda;}}_2^T(\mathrm{KU}-F_{\mathrm{in}})+{\mathrm{\&lambda;}}_3^T({\mathrm{KU}}_d-F_d)+\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_4^i({\mathrm{\&rho;}}_{\min }-{\mathrm{\&rho;}}_i+b_i^2)\\ +\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_5^i({\mathrm{\&rho;}}_i-{\mathrm{\&rho;}}_{\max }+c_i^2)+\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m(G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m}))\end{array};$

式中，λ₁，λ₂，λ₃，和β_m为拉格朗日乘子；V^*为柔顺机构的目标体积，和是松弛因子；

继续得到具有消除棋盘格的Kuhn-Tuck必要条件如下：

$\left\{\begin{array}{c}\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}+{\mathrm{\&lambda;}}_1{V}_i+{\mathrm{\&lambda;}}_2^T\frac{\&PartialD;\left(\mathrm{KU}\right)}{{\&PartialD;\mathrm{\&rho;}}_i}+{\mathrm{\&lambda;}}_3^T\frac{\&PartialD;\left({\mathrm{KU}}_d\right)}{{\&PartialD;\mathrm{\&rho;}}_i}-\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m\frac{\&PartialD;\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})}{{\&PartialD;\mathrm{\&rho;}}_i}\left\{\begin{array}{ccc}=0& \mathrm{if}& {\mathrm{\&rho;}}_{\min }<{\mathrm{\&rho;}}_i<{\mathrm{\&rho;}}_{\max }\\ >0& \mathrm{if}& {\mathrm{\&rho;}}_i={\mathrm{\&rho;}}_{\max }\\ <0& \mathrm{if}& {\mathrm{\&rho;}}_i={\mathrm{\&rho;}}_{\max }\end{array}\right.\\ V\left(\mathrm{\&rho;}\right)=V^{*}\\ \mathrm{KU}-F=0\\ {\mathrm{KU}}_d-F_d=0\\ G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})=0,m\&Element;N_{\mathrm{ck}}\\ i=\mathrm{1,2},...,N\end{array}\right.;$

式中，N_ck是本次优化迭代中所有对棋盘格生成条件恒成立的低密度单元的集合；

首先，得到优化目标、体积约束和消除棋盘格约束的灵敏度分别如下：

优化目标的灵敏度为： $\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}=\frac{E_s\left(P{\left({\mathrm{\&rho;}}_i\right)}^{P-1}{\left(u_i\right)}^T{K}^0{u}_{\mathrm{di}}\right)-E_{\mathrm{ms}}\left(P{\left({\mathrm{\&rho;}}_i\right)}^{P-1}{\left(u_i\right)}^T{K}^0{u}_i\right)}{{\left(E_s\right)}^2},i=1,...,N;$

体积约束的灵敏度为： $\frac{\&PartialD;V}{{\&PartialD;\mathrm{\&rho;}}_i}=V_e,i=1,...,N;$

消除棋盘格约束的灵敏度为：

然后，代入Kuhn-Tuck必要条件，得到灵敏度关系如下：

$Q_i^t=\frac{-\frac{\&PartialD;f}{{\&PartialD;\mathrm{\&rho;}}_i}+\underset{m}{\mathrm{\&Sigma;}}{\mathrm{\&beta;}}_m\frac{\&PartialD;\mathrm{\&Phi;}({\mathrm{\&rho;}}_m,{\mathrm{\&rho;}}_{N_m})}{{\&PartialD;\mathrm{\&rho;}}_i}}{{\mathrm{\&lambda;}}_1\frac{\&PartialD;V}{{\&PartialD;\mathrm{\&rho;}}_i}}=1;$

于是，得到具有消除棋盘格的设计变量更新最佳准则，如下：

${\mathrm{\&rho;}}_i^{t+1}=\left\{\begin{array}{ccc}{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}& \mathrm{if}& \max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })<{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}<\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })\\ \min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })& \mathrm{if}& {\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}\&GreaterEqual;\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })\\ \max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })& \mathrm{if}& {\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}\&le;\max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })\end{array}\right.,i=\mathrm{1,2},...,N;$

式中，t为迭代数,η为松弛因子,且0＜η＜1,ζ为一个较小的移动极限,是一个非负数，

$M_i^t=\max (0,Q_i^t);$

由于是λ₁和β_m的函数，i＝1,2,…,N，因此，在设计变量更新过程中，须采用二分法确定λ₁和β_m的值，以使更新后的密度满足体积约束和消除棋盘格约束，即：

$\begin{array}{c}{\mathrm{\&theta;}}^{*}{V}_0-\{\underset{\mathrm{icas}1}{\mathrm{\&Sigma;}}{V}_e{\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}+\underset{\mathrm{icas}2}{\mathrm{\&Sigma;}}{V}_e\min \{(1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max }\}+\underset{\mathrm{icas}3}{\mathrm{\&Sigma;}}{V}_e\max \{(1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min }\left\}\right\}=0;\\ G-\mathrm{\&Phi;}({\mathrm{\&rho;}}_m^t,{\mathrm{\&rho;}}_{N_m}^t)-\underset{\mathrm{icas}1}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t({\mathrm{\&rho;}}_i^t{\left(M_i^t\right)}^{\mathrm{\&eta;}}-{\mathrm{\&rho;}}_i^t)-\underset{\mathrm{icas}2}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t(\max ((1-\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\min })-{\mathrm{\&rho;}}_i^t)\\ -\underset{\mathrm{icas}3}{\mathrm{\&Sigma;}}{\left(\frac{\&PartialD;\mathrm{\&Phi;}}{{\&PartialD;\mathrm{\&rho;}}_i}\right)}^t(\min ((1+\mathrm{\&zeta;}){\mathrm{\&rho;}}_i^t,{\mathrm{\&rho;}}_{\max })-{\mathrm{\&rho;}}_i^t)=0,\&Element;N_{\mathrm{ck}}\end{array};$

式中，和分别是基于最佳准则所更新的三类设计变量的和。

2.根据权利要求1所述的柔顺机构无棋盘格拓扑图提取方法，其特征在于：步骤4)中，

首先对柔顺机构拓扑优化模型进行有限元分析，得出优化目标、体积约束和消除棋盘格约束的灵敏度；然后计算所有消除棋盘格约束值；最后基于具有消除棋盘格的最佳准则更新设计变量，得到一次迭代优化的拓扑图。