跳转至

牵引法

On the formulation and implementation of geometric and manufacturing constraints in node–based shape optimization

概述

论文将设计节点坐标分为优化节点和依赖节点,对于无法用隐式方法处理的约束(如脱模约束),通过不等式约束附加到优化问题中

论文分为三种节点坐标:

  • 设计节点坐标 \(x_d\): 设计更新位移
  • 受控节点坐标 \(x_c\): 不直接更新但可能受控位移以避免单元畸变
  • 固定节点坐标 \(x_f\): 优化过程中完全固定

受控节点 \(x_c\) 通过辅助问题与设计变量 \(x_d\)相关联, \(x_c=x_c(x_d)\)

位移场 U 通过在有限元分析步骤中求解力学平衡方程 \(KU=F\)得到

其中:

  • K是全局刚度矩阵
  • F是载荷向量

由于 \(x_c\) 是关于\(x_d\)的函数,故位移场 \(u\) 最终是设计变量 \(x_d\)的隐函数

灵敏度分析

首先需要明确以下变量含义

  • \(x_d\):主动节点坐标
  • \(x_c\): 被动节点坐标,\(x_c=x_c(x_d)\)
  • \(u\): 位移场,通过有限元方程\(Ku=F\)求得
  • \(f\): 目标函数,如柔度,\(f=f(u,x_c,x_d)\)
  • \(r\): 系统残差,\(r=KU-F\),在平衡状态下,\(r=0\)
  • \(R\): 定义几何约束的​​辅助问题的残差或最优性条件,辅助问题的解\(x_c^*\)满足\(R(x_c^{*},x_d)=0\),隐式表示\(x_c^*=x_c^*(x_d)\)

分为两个阶段处理

第一阶段

不考虑\(x_c\)\(x_d\)之间的关系,均看成\(x\)

应用链式法则:

\[ \frac{\partial f^{*}}{\partial \mathbf{x}} = \frac{\partial f}{\partial \mathbf{u}} \frac{\partial \mathbf{u}}{\partial \mathbf{x}} + \frac{\partial f}{\partial \mathbf{x}} \qquad{(1)} \]

难点在于计算 \(\frac{\partial \mathbf{u}}{\partial \mathbf{x}}\)。我们从系统方程 $ \mathbf{r}(\mathbf{u}, \mathbf{x}) = \mathbf{0} $ 出发。因为 \(\mathbf{r}\) 恒为 \(\mathbf{0}\),其全导数为 \(\mathbf{0}\)

\[ \frac{d \mathbf{r}}{d \mathbf{x}} = \frac{\partial \mathbf{r}}{\partial \mathbf{u}} \frac{\partial \mathbf{u}}{\partial \mathbf{x}} + \frac{\partial \mathbf{r}}{\partial \mathbf{x}} = \mathbf{0} \qquad{(2)} \]

对于线弹性问题,$ \frac{\partial \mathbf{r}}{\partial \mathbf{u}} = \mathbf{K} $(刚度矩阵)。因此:

\[ \mathbf{K} \frac{\partial \mathbf{u}}{\partial \mathbf{x}} = -\frac{\partial \mathbf{r}}{\partial \mathbf{x}} \qquad{(3)} \]

故(1)式可转化为

\[ \frac{\partial f^{*}}{\partial \mathbf{x}} = -\frac{\partial f}{\partial \mathbf{u}} \mathbf{K}^{-1} \frac{\partial \mathbf{r}}{\partial \mathbf{x}} + \frac{\partial f}{\partial \mathbf{x}} \qquad{(4)} \]

由于求解这个方程若直接求解,计算量巨大,故使用伴随法求解,引入\(\lambda^*\)

\[ \mathbf{K}^{T} \boldsymbol{\lambda}^{*} = -\left( \frac{\partial f}{\partial \mathbf{u}} \right)^T \qquad{(5)} \]

此时只需要求解(5),然后将其代入(4)即可得下式

\[ \frac{\partial f^{*}}{\partial \mathbf{x}} = \boldsymbol{\lambda}^{*T} \frac{\partial \mathbf{r}}{\partial \mathbf{x}} + \frac{\partial f}{\partial \mathbf{x}}. \qquad{(6)} \]

第二阶段

考虑\(x_c\)\(x_d\)之间的关系,将\(f^*\)转换为\(f^{**}\),由于\(x_c\)\(x_d\)表示,故我们讨论对$x_d的导数

\[ \frac{\partial f^{**}(\mathbf{x}_d)}{\partial \mathbf{x}_d} = \frac{\partial f^{**}(\mathbf{x}_c(\mathbf{x}_d), \mathbf{x}_d)}{\partial \mathbf{x}_d} = \frac{\partial f^{*}(\mathbf{x}_c(\mathbf{x}_d), \mathbf{x}_d)}{\partial \mathbf{x}_c} \frac{\partial \mathbf{x}_c}{\partial \mathbf{x}_d} + \frac{\partial f^{*}(\mathbf{x}_c(\mathbf{x}_d), \mathbf{x}_d)}{\partial \mathbf{x}_d} \qquad{(7)} \]

\(x_c\)\(x_d\)之间的关系通过\(R\)来表示,在平衡状态下\(R=0\),故可得到下式

\[ \begin{align*} \mathbf{0} &= \frac{\partial \mathbf{R}^{*}(\mathbf{x}_d)}{\partial \mathbf{x}_d} = \frac{\partial \mathbf{R}^{*}(\mathbf{x}_c(\mathbf{x}_d), \mathbf{x}_d)}{\partial \mathbf{x}_d} \\ &= \frac{\partial \mathbf{R}}{\partial \mathbf{x}_c} \frac{\partial \mathbf{x}_c}{\partial \mathbf{x}_d} + \frac{\partial \mathbf{R}}{\partial \mathbf{x}_d} \\ &\Rightarrow \frac{\partial \mathbf{x}_c}{\partial \mathbf{x}_d} = -\left[ \frac{\partial \mathbf{R}}{\partial \mathbf{x}_c} \right]^{-1} \frac{\partial \mathbf{R}}{\partial \mathbf{x}_d}, \end{align*} \qquad{(8)} \]

带入(7)可得

\[ \frac{\partial f^{**}}{\partial \mathbf{x}_d} = -\frac{\partial f^{*}}{\partial \mathbf{x}_c} \left[ \frac{\partial \mathbf{R}}{\partial \mathbf{x}_c} \right]^{-1} \frac{\partial \mathbf{R}}{\partial \mathbf{x}_d} + \frac{\partial f^{*}}{\partial \mathbf{x}_d}. \qquad(9) \]

同样通过伴随法求解,引入\(\lambda^{**}\)

\[ \left[ \frac{\partial \mathbf{R}}{\partial \mathbf{x}_c} \right]^T \boldsymbol{\lambda}^{**} = -\left( \frac{\partial f^{*}}{\partial \mathbf{x}_c} \right)^T, \]

则(9)可转变为下式

\[ \frac{\partial f^{**}}{\partial \mathbf{x}_d} = \boldsymbol{\lambda}^{**T} \frac{\partial \mathbf{R}}{\partial \mathbf{x}_d} + \frac{\partial f^{*}}{\partial \mathbf{x}_d}. \qquad(13) \]

Warning

特殊情况

如果辅助问题是通过​​最小化一个虚构能量\(I(x_c,x_d)\)来定义的,则在最优点处 \(\frac{\partial{I}}{\partial{x_c}}=0\)

\(R\)\(x_d\)的灵敏度转换为\(I\)\(x_d\)的灵敏度,如下

\[ \frac{\partial \mathcal{I}^{*}(\mathbf{x}_c^{*}(\mathbf{x}_d), \mathbf{x}_d)}{\partial \mathbf{x}_d} = \frac{\partial \mathcal{I}}{\partial \mathbf{x}_c} \frac{\partial \mathbf{x}_c^{*}}{\partial \mathbf{x}_d} + \frac{\partial \mathcal{I}}{\partial \mathbf{x}_d} = \frac{\partial \mathcal{I}}{\partial \mathbf{x}_d}, \]

此时求解即如下:

\[ \mathbf{x}_c^{*}(\mathbf{x}_d) = \arg \min_{\mathbf{x}_c} \mathcal{I}(\mathbf{x}_c, \mathbf{x}_d).\]

论文地址

论文地址