基于等几何分析的形状优化

2026/05/20 14:00:00·2026/05/22 13:50:00

# GAMES302 第 9 讲：基于等几何分析的形状优化

本讲主讲人：徐岗（杭州电子科技大学）。主题从结构优化的一般模型出发，解释为什么等几何分析（Isogeometric Analysis, IGA）天然适合形状优化，并重点展开目标函数、MMA 优化、控制点灵敏度、线弹性/三维/超弹性场景与典型算例。

一、本讲问题：把“几何形状”变成可优化变量

形状优化（shape optimization）研究的是：在边界条件、载荷、材料参数和制造约束给定时，如何改变结构边界或外形，使结构性能指标达到最优。它和尺寸优化、拓扑优化同属结构优化，但三者改变对象不同：

尺寸优化：通常改变厚度、半径、截面尺寸等少量标量参数。
形状优化：改变边界曲线、曲面、控制点位置或权重，拓扑连通关系一般保持不变。
拓扑优化：允许材料分布、孔洞数量、连通关系发生改变，设计自由度更高，L10 会进入这一主题。

GAMES302 L09 的核心是：在 IGA 中，CAD 几何已经由 NURBS/B-spline 控制点、权重和节点向量描述；IGA 又直接用这些样条基函数离散 PDE。因此，形状变量可以直接取为控制点坐标或权重，优化后仍然得到一份 CAD 友好的样条模型。

形状优化的数学表述：目标函数、约束与状态方程共同构成 PDE 约束优化问题。

结构优化设计以力学原理和数学规划算法为基础。它通过最优化方法改变工程结构的尺寸、形状和拓扑构型，在给定设计空间和约束条件下实现最优性能设计。课程 PPT 强调，优化结果通常还需要后处理，例如局部光顺、局部补强、制造约束检查；但优化本身能给出人类设计师不容易直接想到的创新初稿，已经广泛用于航空航天、汽车、机械、土木和船舶工程。

一次完整形状优化通常包含三个环节：

设计模型：用样条曲线、曲面、实体或网格描述几何外形，并规定哪些几何自由度可以移动。
分析模型：对当前形状求解物理响应，例如线弹性位移场、应力场、模态频率或超弹性变形。
优化算法：根据目标函数和约束的灵敏度决定设计变量更新方向，迭代直到收敛。

传统 FEA 形状优化常见流程是 CAD 建模 → 网格生成 → 有限元求解 → 灵敏度分析 → 更新 CAD/网格。每次几何更新后可能都要重新网格化，且 CAD 几何与分析网格之间存在误差。IGA 的出发点正好对准这个痛点：用 CAD 中的 NURBS、T-splines 或相关样条基函数直接作为分析基函数。经典 IGA 论文 Cottrell、Hughes、Bazilevs 等指出，NURBS 可以精确表达许多 CAD 常见几何，并且同一基函数可用于几何映射和场变量近似，这也是本讲形状优化的底层逻辑。

二、前置知识回顾：从 L06-L08 到 L09

这一讲默认已经理解 L06-L08 中 IGA 求解 PDE 的基本链条。形状优化并不是新开一个孤立问题，而是在“PDE 求解器”外面套一层“设计变量迭代”。因此先把前置知识压缩成四条主线。

1. PDE 求解：状态方程决定结构响应

以线弹性为例，未知量是位移场 $$u$$ 。连续问题由三类方程构成：

平衡方程：内部应力与外力平衡，例如 $-\nabla \cdot \sigma = f$ 。
几何方程：应变由位移梯度给出，例如小变形下 $\varepsilon(u)=\frac12(\nabla u+\nabla u^T)$ 。
本构方程：应力与应变由材料矩阵关联，例如 $\sigma=D\varepsilon$ 。

边界条件分为 Dirichlet 约束和 Neumann 载荷：前者规定某些边界上的位移，后者规定外力或表面力。离散后常写为

$$K(s)u(s)=F(s),$$

其中 $$s$$ 是设计变量， $$K$$ 是刚度矩阵， $$F$$ 是载荷向量， $$u$$ 是位移自由度。形状变化会改变物理域 $\Omega(s)$ 、雅可比矩阵、积分权重、边界法向，有时还会改变载荷方向或载荷作用面积，所以 $$K$$ 和 $$F$$ 都可能依赖设计变量。

2. 弱形式：优化时真正微分的是积分表达

IGA 和 FEM 一样不会直接求强形式，而是把 PDE 转成弱形式。例如线弹性问题可写为：求 $u \in V$ ，使得对任意测试函数 $v \in V_0$ ，有

a(u,v;\Omega)=\ell(v;\Omega),

其中

a(u,v;\Omega)=\int_\Omega \varepsilon(v)^T D\varepsilon(u)\,d\Omega, \qquad \ell(v;\Omega)=\int_\Omega f\cdot v\,d\Omega+\int_{\Gamma_t} t\cdot v\,d\Gamma.

这一步很关键：形状优化里的灵敏度，本质上就是对这些随几何变化的积分做链式法则。物理域改变时，积分区域、雅可比行列式 $$|J|$$ 、形函数对物理坐标的导数、边界积分测度都会变化。

3. 控制点：几何映射也是设计变量映射

NURBS 曲线、曲面或实体通过控制点 $$P_A$$ 、权重 $$w_A$$ 和基函数 $R_A(\xi)$ 给出物理坐标：

x(\xi)=\sum_A R_A(\xi)P_A.

如果选择某个控制点坐标 $P_{A,x}$ 作为设计变量，那么物理点对该变量的导数很直接：

\frac{\partial x(\xi)}{\partial P_{A,x}}=R_A(\xi)e_x.

这就是 IGA 形状优化的优势来源：移动一个控制点不会只移动一个网格节点，而是通过局部支撑的样条基函数平滑影响一片几何区域。NURBS 权重也可以作为设计变量，尤其对圆弧、圆孔、二次曲线这类有理几何很重要，但权重优化更容易带来尺度和非线性问题，实践中经常先固定权重，只优化控制点坐标。

4. 灵敏度分析：优化算法需要梯度，而不只是函数值

梯度型优化算法，例如移动渐近线法（MMA）、序列二次规划（SQP）或内点法，需要目标函数和约束函数对设计变量的导数。若只有目标函数值，就只能使用有限差分或启发式搜索，效率会很低。L09 因此把“结构灵敏度分析”作为核心：给定 $$J(s,u(s))$$ ，要求

\frac{dJ}{ds_i}=\frac{\partial J}{\partial s_i}+\frac{\partial J}{\partial u}\frac{\partial u}{\partial s_i}.

难点在于 $\partial u/\partial s_i$ 。它不是直接给定的，而要由状态方程 $$K(s)u(s)=F(s)$$ 的微分推出。这正是后面直接微分法和伴随法的入口。

三、为什么 IGA 适合形状优化

把 IGA 用于形状优化，不只是“换一种有限元基函数”。它改变了设计模型与分析模型之间的关系。

基于 NURBS 控制点的形状参数化：控制点的移动会通过样条基函数平滑地传播到几何边界。

1. 几何精确性：优化过程中尽量不丢 CAD 几何

CAD 中圆、圆弧、圆柱、圆锥曲面等常用几何可以由 NURBS 精确表达。传统网格化会把曲边近似成分片直线或分片多项式曲面，几何误差会进入分析结果；而 IGA 直接在样条几何上积分和求解，更适合需要边界高精度的形状优化，例如带孔板、叶片、扳手开口、流固边界等。

2. 设计变量少而语义明确

若用有限元网格节点作为形状变量，变量数量可能很大，而且移动节点容易造成网格畸变。IGA 中粗控制网格可以表达整体形状趋势，少量控制点即可改变边界曲率、厚度和过渡区域。课程中悬臂梁、1/4 带孔平板、开口扳手等例子都采用控制点坐标或权重作为设计变量。

3. 设计模型和分析模型可以分离但可传递

PPT 后半部分提到：形状优化时可以让设计模型使用粗网格控制形状变化，让分析模型使用细网格计算更准确的位移和灵敏度。NURBS 节点插入给出了粗细模型之间的双射关系，因此分析模型上的灵敏度可以传回粗控制点，粗控制点更新后又能映射到细模型。这避免了“设计变量过多”和“分析精度不足”之间的冲突。

4. 参数域保持稳定

IGA 形状优化中，设计变量改变物理域 $\Omega$ ，但参数域 $\hat\Omega$ 通常保持不变。几何映射 $x(\xi;s)$ 变化，积分仍在标准参数单元上进行：

\int_{\Omega(s)} g(x)\,d\Omega=\int_{\hat\Omega} g(x(\xi;s))|J(\xi;s)|\,d\xi.

因此求导可以集中在 $x(\xi;s)$ 、 $J(\xi;s)$ 、 $B(\xi;s)$ 等量上，而不必每次重建网格拓扑。这也是课程里说“几何域会因设计变量变化而变化，但参数域不会受到影响，从而使单元保持良好质量”的意思。

四、形状优化的数学模型

一般形状优化可写成 PDE 约束优化问题：

\begin{aligned} \min_{s}\quad & J(\Omega(s),u(s),s) \\ \text{s.t.}\quad & K(s)u(s)=F(s), \\ & g_j(\Omega(s),u(s),s)\le 0,\quad j=1,\dots,m,\\ & s_i^{\min}\le s_i\le s_i^{\max}. \end{aligned}

这里 $$s$$ 可以是外边界控制点坐标、控制点权重、FFD 控制点、厚度变量或某些 CAD 参数。 $$J$$ 是目标函数， $$g_j$$ 是约束函数， $$K u=F$$ 是离散后的状态方程。

课程 PPT 中总结了四类最常见目标函数：柔度、节点位移/关键点位移、最大应力、体积或面积。其中体积和最大应力也常作为不等式约束，而不是目标。比如工程上常见写法是：在体积不超过初始体积 35% 的条件下，使加载点位移最小；或在最大应力不超过许用值的条件下，使体积最小。

优化算法方面，本讲重点介绍 MMA（Moving Asymptotes Method）。MMA 的想法是把复杂非线性优化问题在每次迭代中近似成一组更容易求解的凸子问题，并通过上下移动的渐近线控制每个设计变量的更新范围。它适合大规模结构优化，因为它只需要函数值和一阶梯度，能处理很多不等式约束。

五、四类目标函数：柔度、重量/体积、应力、频率/模态

1. 柔度/刚度目标：最常用的结构性能指标

柔度（compliance）可理解为外力做功或结构在载荷下的“软弱程度”。离散线弹性中，若 $$K u=F$$ ，柔度常写为

$$C=F^T u=u^T K u.$$

柔度越小，在同样外力下位移越小，结构越“硬”。因此柔度最小化常被称为刚度最大化。典型问题是：

\min_s C(s)=F(s)^T u(s)\quad \text{s.t.}\quad K(s)u(s)=F(s),\quad V(s)\le V^*.

如果外力不随形状变化，柔度导数有一个很简洁的形式。由 $$C=F^Tu$$ 和 $$Ku=F$$ 可得，在自伴随线弹性问题中，柔度对设计变量的总导数常可写成

\frac{dC}{ds_i}=2\frac{\partial F^T}{\partial s_i}u-u^T\frac{\partial K}{\partial s_i}u,

当 $$F$$ 与 $$s_i$$ 无关时进一步化为 $-u^T K_{,i}u$ 。这就是为什么柔度优化常被用作教学中的第一类例子：物理意义清晰，伴随变量也容易处理。

2. 重量/体积目标：成本、质量与制造约束

体积或面积既可以作为约束，也可以作为目标。二维问题中面积为

A(s)=\int_{\Omega(s)}1\,d\Omega=\int_{\hat\Omega}|J(\xi;s)|\,d\xi.

三维问题中体积为

V(s)=\int_{\Omega(s)}1\,d\Omega=\int_{\hat\Omega}|J(\xi;s)|\,d\xi.

如果材料密度 $\rho$ 为常数，重量就是 $W=\rho gV$ ，质量就是 $M=\rho V$ 。工程上常见两种模型：

性能优先：最小柔度或最小位移，约束 $V\le V^*$ ，控制材料用量。
轻量化优先：最小体积或重量，约束最大应力、最大位移或最低固有频率。

体积灵敏度来自雅可比行列式的导数。若 $$J$$ 是几何映射的雅可比矩阵，则

\frac{\partial |J|}{\partial s_i}=|J|\operatorname{tr}(J^{-1}J_{,i}),

所以 IGA 中只要能写出 $x_{,i}$ ，就能进一步得到 $J_{,i}$ 和体积导数。

3. 应力约束：避免局部破坏

柔度最小化降低的是整体变形能，但不保证每个局部点的应力都安全。带孔板、圆角、开口扳手、叶片根部等位置很容易出现应力集中。应力优化一般有两种写法：

\min_s \sigma_{\max}(s) \quad \text{或} \quad \sigma_{\max}(s)\le \sigma_{\text{allow}}.

直接使用最大值函数会遇到不可导问题，因为最大应力点可能在迭代中跳变。课程 PPT 给出常用的 p-norm 近似：

\|\sigma\|_p=\left(\sum_{q=1}^{n_{pts}} |\sigma_q|^p\right)^{1/p},\qquad p\to\infty.

当 $$p$$ 足够大时，最大应力点的贡献会远大于其他点。PPT 中提到实践测试里 $$p$$ 取 40 以上基本可满足需要。实际工程还常用 KS 函数（Kreisselmeier-Steinhauser function）平滑最大值，因为它数值稳定性更好。

4. 频率/模态目标：让结构避开共振

课程页原始短笔记没有展开模态目标，但完整形状优化中这类目标很重要。结构自由振动满足广义特征值问题：

K(s)\phi_k(s)=\lambda_k(s)M(s)\phi_k(s),\qquad \omega_k=\sqrt{\lambda_k}.

其中 $$M$$ 是质量矩阵， $\phi_k$ 是第 $$k$$ 阶模态。常见目标包括：

最大化第一阶固有频率 $\omega_1$ ，避免低频共振。
约束某些频率不落入工作频带。
控制模态间隔，避免频率交叉导致振型不稳定。

若特征值为单重且模态按 $\phi_k^TM\phi_k=1$ 归一化，则特征值灵敏度可写成

\frac{d\lambda_k}{ds_i}=\phi_k^T\left(K_{,i}-\lambda_k M_{,i}\right)\phi_k.

形状变化会同时改变刚度矩阵和质量矩阵，因此频率优化比柔度优化更敏感，也更容易遇到模态交换问题。复习时可把它和 L09 的柔度/体积/应力放在同一框架下：所有目标都需要 $K_{,i}$ 、 $F_{,i}$ 、 $M_{,i}$ 或几何积分导数。

六、灵敏度分析：从控制点导数到伴随法

结构灵敏度分析求的是目标函数或约束函数对设计变量的导数。PPT 强调，它在结构优化、可靠性评估、参数识别中都是先决条件，优化效率很大程度取决于灵敏度计算的效率和精度。

形状灵敏度推导：核心是对

$K(s)u(s)=F(s)$

和目标函数同时求导。

1. 总导数分解

设目标函数为 $$f(s,u(s))$$ ，则对第 $$i$$ 个设计变量有

\frac{df}{ds_i}=\frac{\partial f}{\partial s_i}+\frac{\partial f}{\partial u}\frac{\partial u}{\partial s_i}.

很多结构目标函数没有显式依赖 $$s_i$$ ，但通过几何积分、载荷边界、面积体积仍可能隐式依赖形状。PPT 中说“ $$f$$ 对 $$x_i$$ 微分的值通常为 0”指的是某些目标函数形式下显式项可以消失，但在完整推导里不要机械忽略几何项。

2. 直接微分法

对平衡方程 $$Ku=F$$ 求导：

K_{,i}u+Ku_{,i}=F_{,i}.

于是

u_{,i}=K^{-1}(F_{,i}-K_{,i}u).

代回目标函数导数：

\frac{df}{ds_i}=f_{,i}+f_u K^{-1}(F_{,i}-K_{,i}u).

这就是直接微分法。它的优点是概念直观；缺点是如果设计变量很多，每个 $$s_i$$ 都要解一次类似 $K u_{,i}=F_{,i}-K_{,i}u$ 的线性系统，成本非常高。IGA 形状优化中控制点可能有几十到几千个自由度，直接法往往不划算。

3. 伴随法

伴随法把 $K^{-1}$ 从每个设计变量的计算中挪出来。定义伴随变量 $\lambda$ 满足

K^T\lambda=f_u^T.

则

f_u K^{-1}(F_{,i}-K_{,i}u)=\lambda^T(F_{,i}-K_{,i}u).

最终

\frac{df}{ds_i}=f_{,i}+\lambda^T F_{,i}-\lambda^T K_{,i}u.

关键优势是：对一个目标函数，只需解一次伴随方程，就能计算所有设计变量的梯度。NGSolve 的 PDE-constrained shape optimization 教程也以类似方式构建状态方程、伴随方程和形状导数；它强调线性问题中伴随算子的左端常是状态算子的转置。这个观点和 L09 的矩阵推导是一致的。

4. IGA 中 $K_{,i}$ 与 $F_{,i}$ 从哪里来

线弹性单元刚度矩阵可写成

K_e=\int_{\hat\Omega_e} B(\xi;s)^T D B(\xi;s)|J(\xi;s)|\,d\xi.

对设计变量求导：

K_{e,i}=\int_{\hat\Omega_e}\left(B_{,i}^TDB+B^TDB_{,i}\right)|J|+B^TDB|J|_{,i}\,d\xi.

载荷向量若为体力或边界力，也有类似形式：

F_e=\int_{\hat\Omega_e}N^T f|J|\,d\xi+\int_{\hat\Gamma_e}N^T t|J_\Gamma|\,d\eta.

因此 $F_{,i}$ 取决于形状变化是否改变体力方向、边界力方向、边界长度/面积以及法向。课程后面提到“设计相关载荷”时，这一项尤其不能忽略。

5. 从控制点到几何导数

IGA 中 $$B$$ 是基函数对物理坐标的导数组合，而物理导数来自参数导数与雅可比矩阵的链式关系。设 $\hat G$ 存储基函数对参数坐标的导数， $$G$$ 存储基函数对物理坐标的导数，则通常有

G=J^{-1}\hat G,\qquad J=\hat G X.

其中 $$X$$ 是控制点坐标矩阵。对设计变量求导：

G_{,i}=-(J^{-1}J_{,i})J^{-1}\hat G,

因为 $\hat G$ 只依赖参数域和基函数，不依赖物理控制点。 $B_{,i}$ 再由 $G_{,i}$ 线性组合得到。雅可比行列式导数为

|J|_{,i}=|J|\operatorname{tr}(J^{-1}J_{,i}).

最后，若设计变量正是第 $$A$$ 个控制点在第 $$d$$ 个方向的坐标，则

X_{,i}=E_{A,d},\qquad x_{,i}(\xi)=R_A(\xi)e_d.

这条链路把抽象的形状灵敏度落到了控制点坐标上：控制点移动 → 物理坐标变化 → 雅可比变化 → $$B$$ 与 $$|J|$$ 变化 → $$K$$ 、 $$F$$ 、目标函数变化。

6. 有限差分、半解析法与复变函数法的位置

PPT 列出了结构灵敏度的多种方法：解析法、有限差分法、半解析法、复变函数法、伴随变量法。它们的取舍如下：

解析法：精度高，适合推导清楚的线弹性和标准单元，但开发成本高。
有限差分法：实现简单，对黑箱求解器友好，但步长选择困难，变量多时成本高。
半解析法：部分导数用解析表达，复杂项用差分折中。
复变函数法：在函数解析性满足条件时可得到很高精度的一阶导，但工程实现有门槛。
伴随变量法：目标/约束数量少、设计变量多时非常高效，是大规模形状优化的主力。

超弹性材料中刚度矩阵与位移相关，材料非线性更强。PPT 说明此时上述线弹性位移导数过程不能直接照搬，课程示例采用中心差分近似位移导数，再代入不同目标函数的导数表达。

七、标准例题与课程算例

例题 1：悬臂梁柔度最小化的灵敏度方向

题目：二维悬臂梁左端固支，右端受竖向集中力。设计变量为上、下边界若干 NURBS 控制点的 $$y$$ 坐标。目标是在体积约束 $V\le V^*$ 下最小化柔度 $$C=F^Tu$$ 。说明一次梯度迭代中需要计算哪些量。

解题框架：

用当前控制点生成几何映射 $x(\xi;s)$ ，装配 $$K(s)$$ 与 $$F(s)$$ 。
求解状态方程 $$K(s)u(s)=F(s)$$ ，得到位移场。
计算柔度 $$C=F^Tu$$ 与体积 $V=\int_{\hat\Omega}|J|d\xi$ 。
对每个控制点变量 $s_i=P_{A,y}$ ，计算 $x_{,i}=R_Ae_y$ ，再得到 $J_{,i}$ 、 $B_{,i}$ 、 $K_{,i}$ 和 $V_{,i}$ 。
若外力不随形状改变，柔度灵敏度可用 $C_{,i}=-u^TK_{,i}u$ ；若载荷边界也变化，则加入 $2F_{,i}^Tu$ 。
把 $C_{,i}$ 和 $V_{,i}$ 交给 MMA，MMA 根据移动渐近线和约束给出下一步控制点更新。

要点：悬臂梁例子看起来只是“边界往上或往下移动”，实际每个控制点移动都会改变一片区域的刚度贡献。柔度梯度为负的方向表示增加该变量会降低柔度，MMA 会在体积约束和上下界内选择可行更新。

例题 2：控制点移动如何影响边界形状

题目：NURBS 曲线 $x(\xi)=\sum_A R_A(\xi)P_A$ 描述结构上边界。若只移动第 $$B$$ 个控制点 $$P_B$$ 的 $$y$$ 坐标 $\Delta s$ ，边界上参数点 $\xi$ 的位移是多少？为什么这比移动有限元边界节点更平滑？

解：新的边界为

x^{new}(\xi)=\sum_{A\ne B}R_A(\xi)P_A+R_B(\xi)(P_B+\Delta s\,e_y).

所以物理点变化为

\Delta x(\xi)=R_B(\xi)\Delta s\,e_y.

NURBS/B-spline 基函数具有局部支撑且通常具备较高连续性，因此影响范围局部、过渡平滑。有限元节点直接移动容易造成相邻单元拉伸、翻转或边界折线不光滑；控制点移动则通过基函数自动形成平滑形变。这也是 IGA 形状优化适合 CAD 边界设计的根本原因。

课程中的典型数值例子

PPT 展示了多个算例：

1/4 带孔平板：在材料体积约束下使柔度最小。初始模型为少量控制点的双二次 NURBS 曲面，设计变量包括若干控制点的 $$x/y$$ 坐标和权重。这个例子常用来展示孔边应力集中与圆弧几何保持。
悬臂梁：设计目标是在体积约束下使荷载点竖向位移最小。PPT 中给出高度范围、控制网格和收敛准则，适合理解控制点移动与边界形变。
开口扳手：在给定材料体积下，使两种加载情况下的加载点垂直位移最小。开口尺寸固定，手柄区域由双二次控制网络表示，体积限制为初始体积的 35%。这是多工况目标的代表。
三维单片/多片模型：包括三维悬臂梁、多片 1/4 带孔板、多片开口扳手，展示 IGA 多片几何上的优化能力。
涡轮叶身结构：以柔度最小化为目标，展示复杂 CAD 曲面和工程边界条件下的形状优化。
超弹性圆角板、体积梁、开口扳手：说明非线性材料中灵敏度计算更复杂，可能需要差分近似。

悬臂梁形状优化迭代过程：形状更新、IGA 求解、灵敏度计算和线搜索反复交替。

八、等几何形状优化算法框架

结合 PPT 中的算法框架，可以把 IGA 形状优化写成如下流程：

读入样条几何模型，包括基函数信息、节点向量、控制点和权重。
确定初始设计模型与分析模型，设置可动控制点、固定边界、对称边界和设计变量上下界。
进行 IGA 求解，得到当前位移解 $U^{(k)}$ 。
计算目标函数、约束函数和形状梯度；必要时求解伴随方程或形状位移方程。
设置初始下降步长或交给 MMA/SQP 产生候选设计变量 $s^{(k+1)}$ 。
更新控制点坐标，例如 $X^{(k+1)}=X^{(k)}+t^{(k)}V^{(k)}$ 。
重新进行 IGA 求解，计算新目标函数值。
若目标下降且约束满足，接受更新；否则缩小步长或调整移动渐近线。
检查收敛条件，例如 $|J^{(k+1)}-J^{(k)}|/|J^{(k)}|<\varepsilon$ 或设计变量变化足够小。

这里有两个容易忽视的工程细节：

固定几何特征：开口扳手例子中螺栓形状和尺寸固定，不应被优化变量破坏。实际问题中孔径、接口、装配边界、载荷作用点经常需要固定。
可制造性与网格质量：控制点虽然给出平滑边界，但仍可能导致局部自交、厚度过小、曲率过大。优化约束应加入最小厚度、曲率、对称性、加工方向等条件。

MMA 的作用可以理解为“梯度信息的稳定使用器”。直接沿负梯度下降可能违反约束或震荡；MMA 用上下渐近线限制变量变化，并为每次迭代构造近似子问题。课程中给出的 $\alpha_i$ 、 $\beta_i$ 、上下渐近线更新规则，就是为了让变量在连续几次同向变化时放宽步长，在方向反复变化时收紧步长。

九、后续衔接：L10 拓扑优化

L09 的形状优化和 L10 的拓扑优化共享同一个优化框架：状态方程、目标函数、约束函数、灵敏度分析、梯度型更新。但它们的设计变量和几何自由度不同。

形状优化：边界移动，拓扑不变。适合已有 CAD 初稿的外形精修，例如圆角、孔边、叶片轮廓、扳手手柄。
拓扑优化：材料密度或相场变量在设计域内变化，可生成新孔洞和新连通关系。适合概念设计阶段，例如找传力路径和材料布局。

L10 中常见的 SIMP 方法会把每个单元或积分点的材料密度 $\rho$ 作为设计变量，刚度写成 $E(\rho)=\rho^pE_0$ 。这和 L09 的控制点变量非常不同：L09 优化的是几何映射 $x(\xi;s)$ ，L10 优化的是材料分布 $\rho(x)$ 。但柔度目标、体积约束和伴随灵敏度形式会再次出现，因此 L09 是理解 L10 的必要铺垫。

一个实际工程流程常常是：先用拓扑优化得到大致传力路径，再用 CAD 重建光滑边界，最后用 IGA 形状优化精修边界和局部应力。这也是为什么本讲强调 CAD 控制点和分析模型统一表示，因为它正好位于“概念拓扑”与“可制造 CAD 几何”之间。

十、复习速查表

概念	公式/关键词	考试与复习要点
状态方程	$$K(s)u(s)=F(s)$$	形状变量通过 $$K$$ 、 $$F$$ 和几何积分影响位移。
弱形式	$a(u,v)=\ell(v)$	形状导数本质上是对积分域、雅可比、物理导数求导。
NURBS 几何	$x(\xi)=\sum_A R_A(\xi)P_A$	控制点坐标和权重可作为设计变量。
控制点灵敏度	$x_{,i}=R_Ae_d$	第 $$A$$ 个控制点第 $$d$$ 维移动对物理点的影响由基函数决定。
柔度	$$C=F^Tu=u^TKu$$	最小柔度等价于最大刚度，常配体积约束。
体积/面积	$V=\int_{\hat\Omega}\|J\|d\xi$	导数依赖 $\|J\|_{,i}=\|J\|\operatorname{tr}(J^{-1}J_{,i})$ 。
应力目标	$\\|\sigma\\|_p=(\sum \|\sigma_q\|^p)^{1/p}$	用 p-norm 或 KS 函数平滑最大应力。
模态目标	$K\phi=\lambda M\phi$	频率优化需考虑 $K_{,i}$ 与 $M_{,i}$ ，注意模态交换。
直接微分法	$Ku_{,i}=F_{,i}-K_{,i}u$	每个设计变量一次线性求解，变量多时成本高。
伴随法	$K^T\lambda=f_u^T$	一个目标函数通常只需一次伴随求解，适合多设计变量。
MMA	Moving Asymptotes	用移动渐近线构造近似子问题，稳定处理约束优化。
L10 衔接	SIMP / 拓扑优化	形状优化改边界，拓扑优化改材料分布和连通关系。

十一、参考来源

本地 PPT / 源笔记

GAMES302《等几何分析》Lecture 09：基于等几何分析的形状优化，徐岗，杭州电子科技大学。本次更新已读取本地 PDF：/Users/zhengxinyu/MyDoc/高级图像学/GAMES302-slides/GAMES302-09.pdf。
本地课程源笔记：/Users/zhengxinyu/MyDoc/高级图像学/课程笔记/games302-09-iga-shape-optimization.md。原笔记提供了形状优化概述、四类目标函数和 IGA 优势的基础文本。

网页参考

Cottrell, J. A., Hughes, T. J. R., Bazilevs, Y. 等关于 Isogeometric Analysis 的经典资料与论文条目：强调 CAD、有限元、NURBS、精确几何和网格细化的统一，是 IGA 适合形状优化的理论背景。
Wall, W. A., Frenzel, M. A., Cyron, C. Isogeometric structural shape optimization / NURBS based isogeometric shape optimization 相关论文条目：讨论 NURBS 控制点位置与权重的解析灵敏度，支撑“控制点作为设计变量”的推导。
NGSolve Documentation, “PDE-Constrained Shape Optimization (fully automated)”：给出状态方程、伴随方程和形状导数的教程式实现，作为伴随法和 PDE 约束优化框架的网页参考。
ASME Digital Collection, “Isogeometric Shape Optimization for Design Dependent Loads”：强调 NURBS 设计模型可以用较少控制点表达复杂几何，并讨论设计相关载荷下 IGA 形状优化需要处理的载荷灵敏度。