有限变形晶体塑性快速傅里叶变化CPFFT实现 原创

原始文献：《An elasto-viscoplastic formulation based on fast Fourier transforms for the prediction of micromechanical fields in polycrystalline materials》

DOI：10.1016/j.ijplas.2011.12.005

  在计算微观力学领域，如何高效预测多晶体内部的异质应力场量一直是核心难题。自 1990 年代 Moulinec 和 Suquet 提出基于 FFT 的线性解析法以来，谱方法凭借其无需网格划分、直接处理微观图像的优势，迅速成为挑战传统有限元法（FEM）的利器。

  然而，早期的 FFT 框架大多局限于刚塑性或线性弹性。2012 年，Ricardo A. Lebensohn 等人在 International Journal of Plasticity 上发表了具有里程碑意义的工作。该研究提出了一套严谨的弹性-黏塑性（EVP-FFT）公式，能够同时处理晶体的弹性各向异性与非线性滑移演化，为预测多晶材料在复杂载荷下的局部力学响应奠定了理论基础。

Lebensohn 等人的文章重点解决了以下几个力学与数值上的关键问题：

增广拉格朗日迭代 (Augmented Lagrangian)

  针对 EVP 本构中极强的非线性，文章引入了增广拉格朗日迭代程序。这种方法通过在傅里叶空间平衡相容性（Compatibility）与在实空间平衡本构关系，极大地提高了求解高对比度异质材料时的收敛稳健性。

有限应变运动学 (Finite Strain Kinematics)

  在有限变形框架下，总变形梯度被分解为弹性和塑性两部分。文章强调了在参考构型下求解第一类 Piola-Kirchhoff 应力平衡的重要性，这确保了在大旋转、大应变工况下计算的物理准确性。

应力共轭与本构更新

  为了保证能量守恒，文章在晶体本地坐标系下采用 Mandel 应力作为滑移驱动力，并配合隐式时间积分更新塑性变形梯度。

文章的模拟效果如下：

需要注意的是当前的这代积分方案和damask的快速傅里叶变化方案计算效果基本保持一致，整体也是使用fortran语言编写，并使用vtk格式用于输出，使用paraview可视化。

使用类似的思想，我们可以根据文章的公式实现对应的CPFFT的计算方案。这里展示使用matlab实现对应的CPFFT方案，matlab的显著优势可以很容易和相场和再结晶去结合，因此后续非常容易扩展。

使用FFT作为边值问题的求解器，使用固定点迭代完成内部的晶体塑性迭代。使用经典的位错密度模型计算硬化和热激活流动方程计算滑移系的剪切变形。

初始RVE模型使用neper建模，建立一个包含100个晶粒的多晶模型：

matlab导入几何模型网格：

并沿着X方向进行1.0%的拉伸变形，所有量纲使用m-s-pa。

拉伸变形结束后的累计剪切滑移结果：

拉伸变形结束后的统计储存位错密度分布结果：

拉伸变形结束后的几何必须位错密度分布结果：