心脏瓣膜流固耦合分析
理论与物理
心脏瓣膜模拟的背景
心脏瓣膜的流固耦合分析，在什么场景下会需要用到呢？
用于人工心脏瓣膜（机械瓣膜·生物瓣膜）的设计与评估。主要目的是评估瓣叶的开闭行为、壁面剪切应力（WSS）、血栓形成风险、溶血风险。
通过瓣膜的血流雷诺数在数千量级，处于过渡区，且瓣叶会发生大变形，因此流固耦合（FSI）是不可或缺的。
控制方程
血流的力学和牛顿流体不一样吗？
在大血管内，通常将血液近似为牛顿流体（$\mu \approx 3.5$ mPa·s），但在低剪切率区域，则需要Carreau-Yasuda模型等非牛顿模型。
瓣叶的结构分析使用超弹性模型。生物瓣膜采用Mooney-Rivlin模型，机械瓣膜的瓣叶则适用刚体运动模型。
其中 $W$ 是应变能密度函数，$I_1, I_2$ 是柯西-格林变形张量的不变量。
界面条件和通常的FSI一样吗？
基本上相同，但由于瓣叶是薄壳结构，需要额外处理接触判定（瓣膜关闭时瓣叶之间会接触）。使用浸入边界法或浸入有限元法，可以在不重新生成网格的情况下处理大变形。
心脏瓣膜是“完全被动系统”——仅靠血压差每天开闭10万次的精妙设计
对于设计人造物的工程师来说，心脏瓣膜的工作原理简单得令人惊讶。没有马达，也没有电信号，仅靠心室与主动脉之间的压力差，瓣膜就打开；当血液开始逆流时，依靠血液的动量，瓣膜就关闭——这是一个完全被动的流体力驱动系统。支配这种开闭的理论是FSI教科书式的例题，瓣膜的厚度（约0.5mm）和瓣尖的弹性模量直接关系到开放速度和关闭时的逆流量。特别是FSI分析表明，关闭时会产生“水击压力（水锤）”，在瓣尖瞬间产生0.1~0.3MPa的应力集中。这种反复应力会导致瓣膜钙化和疲劳破坏，因此理解该理论直接关系到人工瓣膜的寿命设计。
各项的物理意义
- 结构-热耦合项：温度变化引起的热膨胀诱发结构变形，变形又影响温度场。$\sigma = D(\varepsilon - \alpha \Delta T)$。【日常例子】夏天铁轨伸长导致缝隙变窄——温度升高→热膨胀→产生应力的典型例子。电子电路板在焊接后翘曲，也是不同材料热膨胀率差异造成的。发动机的气缸体因高温部和低温部的温差产生热应力，最坏情况下会导致裂纹。
- 流固耦合（FSI）项：流体压力·剪切力使结构变形，结构变形又改变流体区域，是双向相互作用。【日常例子】强风使吊桥缆索振动（涡激振动）——风力使结构摇晃，摇晃的结构改变气流，进而放大振动。心脏血流与血管壁的弹性变形、飞机机翼的颤振（气动弹性不稳定性）也是典型的FSI问题。有时单向耦合即可，但变形较大时双向耦合是必需的。
- 电磁-热耦合项：焦耳发热 $Q = J^2/\sigma$ 引起温度升高，温度变化又改变电阻，形成反馈回路。【日常例子】电炉的镍铬合金线通电后发热（焦耳热）变红——温度升高电阻改变，电流分布也变化。IH电磁炉的涡流发热、输电线路因温度升高导致垂度增加也是这种耦合的例子。
- 数据传递项：通过插值解决不同物理场间网格不匹配的问题。【日常例子】天气预报中，将“气温数据”和“风速数据”结合起来计算体感温度时，如果各自的观测地点不同就需要插值——在CAE的耦合分析中，结构网格和CFD网格通常也不一致，因此界面处的数据传递（插值）精度直接关系到结果的可信度。
假设条件与适用范围
- 弱耦合假设（单向耦合）：当一方物理场影响另一方，而反向影响可忽略时有效
- 需要强耦合的情况：FSI中的大变形、电磁-热耦合中温度依赖性较强时
- 时间尺度的分离：当各物理场的特征时间差异很大时，可通过子循环提高效率
- 界面条件的协调性：需确认耦合界面处的能量·动量守恒在数值上得到满足
- 不适用的案例：当三个以上物理场同时强耦合时，有时需要整体式方法
量纲分析与单位制
| 变量 | SI单位 | 注意事项·换算备忘 |
|---|---|---|
| 热膨胀系数 $\alpha$ | 1/K | 钢：约12×10⁻⁶、铝：约23×10⁻⁶ |
| 耦合界面力 | N/m²（压力）或N（集中力） | 确认流体侧与结构侧的力平衡 |
| 数据传递误差 | 无量纲（%） | 插值精度取决于网格密度比。建议控制在5%以下 |
数值解法与实现
离散化方法的选择
心脏瓣膜FSI分析中使用的数值方法有哪些种类？
主要有三种方法。
| 方法 | 流体 | 结构 | 特点 |
|---|---|---|---|
| ALE-FEM | FVM/FEM（动网格） | FEM | 界面精度高。大变形时需要重划网格 |
| IB法 | FDM/FVM（固定网格） | 纤维模型 | 无需重划网格。界面会模糊 |
| IFEM | FEM（固定网格） | FEM（浸入式） | 结构可使用FEM。实现稍复杂 |
像心脏瓣膜这样反复开闭的情况，用ALE法的话每个周期都要重划网格，很麻烦吧？
因此近年来，IB法或重叠网格法正逐渐成为主流。Griffith等人（IBAMR）的开源IB代码在心脏瓣膜研究中被广泛使用。
时间积分
心脏搏动周期要用多大的时间步长来求解呢？
心跳周期约为0.8秒。瓣膜开闭发生在数十毫秒内，因此需要 $\Delta t = 0.1$~$0.5$ ms 左右。一个周期需要1,600~8,000步。
为了排除初始瞬态影响，至少运行3~5个周期，统计量在稳定后的周期中获取。入口边界条件会设定MRI测量得到的流速波形或压力波形。
壁面剪切应力的评估
您提到WSS评估很重要，具体使用哪些指标呢？
时间平均WSS（TAWSS）和振荡剪切指数（OSI）是代表性指标。
低TAWSS（< 0.4 Pa）且高OSI（> 0.3）的区域被认为血栓风险较高。FDA（美国食品药品监督管理局）的指南中也推荐进行这些评估。
瓣尖接触的数值表达——如何将“瓣膜关闭”写入代码
心脏瓣膜FSI分析中，数值上最困难的就是“瓣尖之间的接触”处理。当三片瓣尖在中央完全闭合时，瓣尖间的间隙趋近于零。从流体分析角度看，“间隙变为零的瞬间压力趋于无穷大”存在奇异性，通常的CFD网格无法收敛。为了避免这个问题，浸入边界法（IBM）将瓣尖表达为“虚拟的体积力源”，在不改变网格的情况下近似模拟关闭。另一方面，CEL法中，当欧拉网格与瓣尖重叠时，通过改变流体的体积分数机制来表现关闭。两者都是通过不直接处理“零间隙”来规避奇异性的技巧——这是数值方法选择决定分析成败的典型案例。
整体式方法
将所有物理场作为一个联立方程组同时求解。对于强耦合问题稳定，但实现复杂且内存消耗大。
分区法（分离迭代法）
独立求解各物理场，在界面交换数据。易于实现，可复用现有求解器。适用于弱耦合。
界面数据传递
最近邻法（最简单但精度低）、投影法（具有守恒性）、RBF插值（对网格不一致鲁棒性强）。需要在守恒性和精度之间取得平衡。
子迭代
在每个耦合步内进行充分的迭代，确保界面条件的协调性。残差基准基于各物理场的典型值进行缩放。
Aitken松弛
自动调整耦合迭代的松弛系数。防止因过度松弛导致发散，是一种能加速收敛的自适应方法。
稳定性条件
注意附加质量效应（流固耦合中结构密度≈流体密度时）。不稳定时可应用Robin型界面条件或IQN-ILS法。
Aitken松弛的比喻
Aitken松弛类似于“平衡跷跷板”。一方推得太用力，另一方就会弹起，其反作用力又导致推得过猛——为了抑制这种振荡，自动调整推力大小的就是Aitken松弛。当耦合迭代振荡不收敛时，根据上次的修正量自动调整下次修正量的自适应方法。
实践指南
模型构建的实践步骤
请告诉我具体的模型构建步骤。
1. 基于CT/MRI图像的3D几何重建（使用Mimics, 3D Slicer等软件）
2. 瓣叶形状的CAD建模（患者特异性或人工瓣膜的CAD数据）
3. 网格生成：流体域使用多面体网格，结构使用壳单元
4. 边界条件设定：入口给定流量波形，出口使用Windkessel压力模型
5. 材料特性定义：瓣叶、主动脉壁、血液
Windkessel模型是什么？
是将血管系统下游近似为RC电路的模型。三元件Windkessel是标准模型，
$R_p$是近端阻力，$R_d$是远端阻力，$C$是顺应性。这样可以在主动脉瓣下游施加生理学上合理的压力边界条件。
网格质量的基准
需要多大程度的网格密度呢？
瓣叶附近至少需要配置3层以上的棱柱层。文献中通常总单元数在300万~1,000万左右。
| 区域 | 单元尺寸 | 备注 |
|---|---|---|
| 瓣叶表面 | 0.2~0.5 mm | 应力评估所需的精度 |
| 瓣口射流区域 | 0.3~0.8 mm | 解析速度梯度 |
| 主动脉窦 | 0.5~1.0 mm | 解析涡流 |
| 远场区域 | 1.0~3.0 mm | 确保计算效率 |
这种规模的计算需要HPC（高性能计算）吧。
典型的计算需要64~256核，耗时数天到一周左右。使用GPU兼容求解器（如Ansys Fluent Native GPU Solver等）有望加速计算。
なった
詳しく
報告