ABAQUS蠕变分析全流程：解耦加载、参数拟合与CETOL设置

简介：这份文档围绕 ABAQUS 中材料蠕变问题的有限元计算流程展开，面向已具备一定非线性分析基础、需要处理塑料及金属合金长期受力变形问题的工程人员与研究生。内容先说明将分析拆为静态加载与蠕变过程两阶段的原因，再梳理 Power-law、Hyperbolic-sine law 等蠕变模型的参数获取与拟合思路，并提示蠕变分析步中时间增量、蠕变应变容差及积分算法的设置要点，有助于理解从材料曲线到仿真参数的完整链路。压缩包仅含 1 个 doc 文件，约 198KB，篇幅精炼，适合作为流程速查与参数设置参考。目前已有 1163 人学习，说明其在蠕变分析入门与工程排错场景中具有一定参考价值。读者可据此建立两阶段计算框架，弄清模型参数、容差范围与分析步设置之间的取舍，减少因增量步过小或参数不当造成的效率与精度问题。

1. 蠕变计算的解耦：为什么把加载和蠕变拆成两步算

做过有限元蠕变分析的人多半有个直觉：真正烧机时的是那个总时间动辄十万秒的蠕变分析步。可翻车的地方往往不在这里。一套完整的 ABAQUS 蠕变问题计算流程里，最容易被糊弄过去的反而是开头那个"与时间无关"的静态加载步——它先把应力场定死，后面的时间积分全在这个应力场上滚雪球。教科书喜欢把蠕变讲成"材料在恒载下慢慢变形"，但落到求解器里，蠕变是一组带时间项的率方程，它必须拿到一个已收敛、干净的初始应力场才能启动。所以真正值得搞懂的，不是那个十万秒怎么等，而是这个应力场凭什么可信、参数从哪来、容差凭什么定。

2. 静态加载步：蠕变分析为什么非要一个"与时间无关"的起点

2.1 解耦的代价与收益：蠕变应变和塑性应变为什么要分开算

蠕变本构是率相关的，塑性本构在率无关框架下是率无关的，两者同时激活时，本构积分的返回映射要在一个 Newton 迭代里同时处理屈服面、塑性流动和蠕变率，切线刚度矩阵的推导复杂度直接翻倍，条件数也更容易变差。工程上常见的做法是把这段耦合拆开：第一步只解弹塑性加载，蠕变不参与；第二步在固定荷载下让蠕变独自演化，塑性因为屈服面不再扩张而基本冻结。这样每一步的本构积分都只面对一种非弹性机制，收敛路径清晰，出错时也容易定位是哪一段的问题。

代价是这不是数学上等价的分解，而是一次有损近似。被丢掉的项是加载阶段本身积累的那部分蠕变应变。判断能不能这么拆，有个朴素的判据：加载持续时间要远小于蠕变的特征时间。压力容器打压几秒钟、结构自重施加几分钟，和几万小时的蠕变寿命比，前者积累的蠕变量通常在千分之一量级，忽略掉可以接受。反过来说，如果工况本身就是缓慢升压、持续几个小时才到工作载荷，那加载段和蠕变段的时间尺度已经同量级，强行拆开就会系统性低估总变形，这时候老老实实做耦合分析更靠谱。

还有一个容易被忽略的收益：拆开之后，静态步可以用较大的增量步快速走完几何非线性最剧烈的阶段，蠕变步则专心对付时间积分。两种物理过程的数值难度完全不同，硬塞进一个分析步里，增量步会被最苛刻的那个机制拖住。

2.2 几何非线性与增量步：静态步里两个开关的真实含义

静态步的设置核心就两件事：几何非线性开不开，增量步怎么长。前者决定应力用哪种度量，后者决定这一段能不能收敛。

蠕变本构里的等效应力 q 在 ABAQUS 里取的是真应力（Cauchy 应力）的偏量部分。小变形、小转动时工程应力和真应力差别可以忽略；一旦结构出现明显的大变形或大转动，比如薄壁件蠕变后的持续鼓胀，工程应力和真应力能差出几十个百分点，蠕变率又是应力的 n 次方，n 通常取 3 到 8，这点应力误差会被指数级放大。所以大变形工况下静态步必须打开 NLGEOM，让求解器在大变形框架下更新构形和应力。开了 NLGEOM 之后收敛会变难，增量步需要切得更细，这是必须付的账。

另一个容易被误解的是静态步里的"时间"。ABAQUS/Standard 的静态步中，时间是一个伪时间，只用来对荷载和边界条件做插值，不进入任何本构方程。所谓"与时间无关的加载过程"，准确意思是这一段物理时间不参与本构演化，而不是说求解器不看这个数。它决定了荷载从零加到满值分几步走，仅此而已。

** 静态加载步：只解弹塑性，蠕变尚未激活 *STEP, NAME=LOAD, NLGEOM=ON, INC=200 *STATIC 0.1, 1.0, 1E-05, 0.1 ** 边界条件与荷载 *BOUNDARY FIX, 1, 3, 0.0 *CLOAD TIP, 2, -1000.0 *OUTPUT, FIELD, FREQUENCY=1 *NODE OUTPUT U, RF *EL OUTPUT S, E, PE *END STEP

*STATIC数据行的四个数依次是初始增量、该步总时间、最小增量、最大增量。初始增量给 0.1 意味着第一步就走完十分之一，配合默认的 ramp 幅值，荷载在这一步内线性加到 10%。最小增量 1E-05 是收敛失败时允许切到的最小值，再往下切求解器就报错退出。最大增量 0.1 是上限，防止荷载跳变太快导致 Newton 迭代从远离解的初值出发。这套参数对大多数弹塑性预加载够用；如果结构里有接触或者屈服面扩张剧烈，把初始增量降到 0.01 会更稳。

*EL OUTPUT里特意输出了PE，等效塑性应变。它是状态变量，会原封不动传给下一个分析步，蠕变步里判断塑性是否还在发展就靠它。

2.3 应力场的交接：连续分析步和拆成两个作业差在哪

最省事的做法是把静态步和蠕变步写进同一个 inp 文件，两个*STEP首尾相接。ABAQUS 会在第一个步收敛后，把节点位移、单元应力、所有状态变量（等效塑性应变、背应力、*DEPVAR 定义的用户变量）作为第二个步的初值，不需要任何额外操作。这是首选路径，因为不存在数据丢失的环节。

如果因为计算资源或者调试需要把两段拆成两个作业，交接方式就变成了 restart。第一个作业在静态步末尾用*RESTART, WRITE, FREQUENCY=1写出重启动文件，第二个作业读进来继续。这里有个血泪经验：restart 文件里的状态变量完整性依赖于材料定义一致，如果两次作业之间改了*DEPVAR的数目，或者加删了材料子程序，读进来的状态变量会对不上号，求解器不报错，但蠕变率会从一个莫名其妙的历史值开始算，结果看着像对的，其实是黑匣子。

至于*IMPORT，它主要用于跨求解器传递（比如 Standard 到 Explicit），在 Standard 到 Standard 的纯蠕变接力里用得少，配置反而更繁琐。同一份 inp 写两个 step，是最不容易出错的选择。

3. Power-law 蠕变模型：为什么微分形式必须积分成时间显式形式才能拟合

3.1 微分形式与积分形式：一个 t^m 引出的形式选择

ABAQUS 内置的蠕变模型里，Power-law 是最常用的一支，它分 time hardening form 和 strain hardening form 两种写法。time hardening 的微分形式是：

ε̇_cr = A · q^n · t^m

其中 ε̇_cr 是等效蠕变应变率，q 是等效偏应力，t 是总时间，A、n、m 是三个待定材料常数。这个式子是求解器真正在用的形式——每个增量步里，本构积分拿当前应力和当前时间算出蠕变率，再乘时间增量累加到蠕变应变上。

问题出在参数标定这一环。实验给的是恒应力下的蠕变曲线族：几条不同应力水平的 ε-t 曲线，横轴时间、纵轴蠕变应变。要对微分形式做最小二乘拟合，得先对时间求导，而实验曲线上的数值导数噪声极大，尤其在早期段，微分一次误差能放大一个数量级。所以标准做法是先解析积分，把微分方程变成时间显式表达式：

ε_cr = A/(m+1) · q^n · t^(m+1)

积分的前提是应力 q 在恒载试验中保持不变，这正好对应曲线族的测试条件。积分之后，左边是实验直接测到的应变量，右边是三个参数的显函数，最小二乘可以直接用，不需要任何数值微分。这一步看着只是数学变形，实际上是整个标定流程能站得住的前提。

m 的符号还携带物理信息。m 小于零对应蠕变率随时间衰减的一次蠕变阶段，m 等于零对应稳态蠕变，m 大于零对应加速段。单一 Power-law 无法描述完整的三阶段曲线，特别是三次蠕变，它只能覆盖一次段和稳态段。如果实验曲线尾巴明显上翘，那就是损伤或者颈缩在起作用，得换更复杂的模型或者引入损伤变量，硬用 Power-law 拟合会把整个参数带偏。

3.2 time hardening 与 strain hardening：恒载下等价，变载下分家

两种形式的差别，本质上是"用什么量来代表历史"。

time hardening 认为蠕变率只取决于当前应力和当前时间，应力一变，之前那段历史就被"忘记"了——新的时间点从零开始按新应力算率。strain hardening 则用当前的累积蠕变应变代替时间，把历史压缩进这个标量里，应力变化后蠕变率会沿着与已积累应变量对应的位置继续走。ABAQUS 里的 strain hardening 形式是：

ε̇_cr = { A · q^n · [(m+1) · ε_cr^m] }^(1/(m+1))

把它和 time hardening 的积分式对照就能看出，这是把积分式反解回率形式、再把 t 替换成 ε_cr 得到的。所以在恒应力条件下，两种形式给出完全相同的曲线，参数可以互相换算。分歧只在应力路径变化时出现：卸载、再加载、应力松弛、多轴加载下主应力方向旋转，这些场景里 time hardening 会低估或者高估历史效应。

原文提到 time hardening 适用于"蠕变过程应力变化范围不太大"的简单工况，这个判断是对的。实际工程里，恒温恒载的长期变形评估占大多数，用 time hardening 完全够，形式简单、参数少、拟合稳定。只有在载荷谱复杂、需要模拟加载历史和应力松弛时，才值得上 strain hardening，代价是拟合难度上升，因为它的参数不再能从业已积分的曲线族里直接线性求解。

还有一支 hyperbolic-sine 形式，蠕变率里带[sinh(B·q)]^n。它的价值在应力跨度大的场合：幂律在应力外推到高区间时会指数爆炸，双曲正弦的增长率更温和，能同时覆盖低应力蠕变和高应力蠕变两段。

3.3 从曲线族反解 A、n、m：拟合里最容易翻车的一步

拿到三条不同应力水平的蠕变曲线，怎么把 A、n、m 揪出来，是这份流程里技术含量最高的部分。

import numpy as np from scipy.optimize import least_squares # 实验数据：每个应力水平一组 (time, creep_strain) # 时间单位：小时；应力单位：MPa data = [ (10.0, np.array([0.0, 100.0, 500.0, 1000.0, 5000.0]), np.array([0.0, 0.0012, 0.0028, 0.0041, 0.0089])), (20.0, np.array([0.0, 100.0, 500.0, 1000.0, 5000.0]), np.array([0.0, 0.0021, 0.0049, 0.0072, 0.0156])), (30.0, np.array([0.0, 100.0, 500.0, 1000.0, 5000.0]), np.array([0.0, 0.0034, 0.0080, 0.0118, 0.0255])), ] def creep_integral(t, q, A, n, m): # Power-law time hardening 的积分形式 # t 必须为正，拟合时把 t=0 的点排除或加极小量 return A / (m + 1.0) * np.power(q, n) * np.power(t, m + 1.0) def residual(p): A, n, m = p res = [] for q, t, eps in data: mask = t > 0 # 跳过起点的零时刻 pred = creep_integral(t[mask], q, A, n, m) # 相对误差加权：避免数值大的长时间点主导整个残差 res.append((pred - eps[mask]) / np.maximum(eps[mask], 1e-6)) return np.concatenate(res) # 初值很关键：n 一般在 3~8，m 在 -1~0 之间 sol = least_squares(residual, x0=[1e-6, 5.0, -0.8], bounds=([0.0, 0.1, -0.99], [1.0, 20.0, 3.0])) A, n, m = sol.x print("A=%.6e n=%.4f m=%.4f" % (A, n, m))

这段代码有几个刻意为之的地方。第一，残差用相对误差而不是绝对误差，因为曲线族的时间跨度常常跨三四个数量级，最长时间点的应变量是早期点的十几倍，用绝对误差拟合等价于只关心曲尾，短时段的精度会被牺牲掉。第二，三个参数一起拟合，跨所有应力水平共享，这比逐条曲线单独拟合再取平均要稳妥得多，后者会把实验噪声放大成参数离散。第三，n和m都设了边界，m下界卡在 -0.99 是为了避免m+1趋零导致积分式除零，同时这个区间覆盖了绝大多数金属和塑料的实测值。

提示：拟合出来的 A 不是无量纲数，它的量纲由 n 和 m 决定。应力换单位、时间换单位，A 的数值必须重算，不能按比例缩放。

拿到 A、n、m 之后，在 CAE 的材料模块里添加 Mechanical → Creep，选 Power-law，输入这三个数，*CREEP, LAW=TIME就会写进 inp。至此材料侧的准备结束，剩下的全是分析步的事。

4. 蠕变分析步：粘性分析步凭什么敢忽略惯性

4.1 *VISCO 的时间轴：总时间、初始增量和自动增量

蠕变的两大特征是时间相关、惯性效应不明显。前者要求本构方程里出现真实时间，后者意味着加速度项可以整个扔掉。ABAQUS 针对这一类过程提供了专门的粘性分析步*VISCO，它做的是准静态求解，不组装质量矩阵，时间只用来推进本构，不用来做动力学积分。

** 蠕变分析步：真实物理时间在这里登场 *STEP, NAME=CREEP, NLGEOM=ON, INC=10000 *VISCO, CETOL=1E-5 1E-3, 1.0E6, 1E-8, 1.0E4 *BOUNDARY FIX, 1, 3, 0.0 *END STEP

*VISCO数据行同样是四个数：初始时间增量、该步总时间、最小时间增量、最大时间增量。这里的总时间 1.0E6 是真实的物理时间，单位必须和材料参数里时间的单位严格一致——如果 A 是按小时拟合的，这里填的就是小时。这是单位制错误最常见的爆发点，很多"结果差六个数量级"的事故都出在这一行。

初始增量给 1E-3 是相对总时间的一个极小值，目的是让第一步落在瞬态区间里，把一次蠕变阶段的率变化捕捉住。如果第一步就跨到 1E3，隐式积分会直接跳到接近稳态的蠕变应变上，早期那段高蠕变率被整个抹掉，后面的曲线整体下移。

INC=10000是增量步数上限，蠕变分析动辄几千步，默认的 100 步根本不够，跑到一半报"too many increments"是常态。

4.2 CETOL 容差：一个决定增量步大小的隐式旋钮

CETOL是蠕变分析里最关键的数值参数，它控制的是每个增量步内蠕变应变增量的局部积分误差。求解器在试探步长时会估算这个误差，超过 CETOL 就把增量步砍半重来，远低于它就尝试放大步长。

容差和增量步大小是反比关系，而且不是线性的。CETOL 从 1E-4 收紧到 1E-6，增量步数可能翻五到十倍，计算时间同比例上涨。原文给出的 1E-4 到 1E-6 区间是经验之谈，默认值在 1E-5 量级。实际取值要看你要什么：

容差收紧到一定程度会遇到硬墙：当所需的时间增量小于数据行第四个数（最小增量）时，求解器不再切分，直接判定不收敛。这时候要么放宽 CETOL，要么把最小增量再调小两个数量级，后者会让计算时间爆炸。所以合理的流程是先用 1E-4 试算一遍，看.sta文件里增量步的分布——如果大部分步长都被容差压得很小，说明模型本身的蠕变率变化剧烈，这时候收紧 CETOL 收益不大，反而是网格或者材料参数有问题。

4.3 积分算法与收敛：隐式蠕变积分在求解器内部做了什么

*VISCO用的是隐式时间积分。每个增量步内，求解器对蠕变率方程做隐式离散，把末时刻的蠕变应变写成应力和时间的隐函数，然后和平衡方程一起做 Newton 迭代。隐式的最大好处是无条件稳定：时间步长再大也不会像显式那样发散，只是精度会掉。这就是为什么蠕变分析可以放心用自动增量，而不用像显式动力学那样盯着稳定时间步。

代价是每个增量步都要更新一致切线刚度。蠕变柔度矩阵和弹性刚度叠加后，整体刚度矩阵的条件数会变差，特别是当蠕变应变远大于弹性应变时，矩阵接近奇异，Newton 迭代的收敛速度会明显下降。这时候经常看到.msg文件里出现"attempts exceeded"，然后增量步被自动砍半。

有一个细节值得注意：蠕变步里的塑性和蠕变是可以同时激活的。但因为静态步结束时荷载已经加满，蠕变步里没有新的加载，屈服面不再扩张，塑性应变基本保持不变，实际参与演化的只有蠕变。这正是第 2 章解耦策略能成立的根本原因——不是求解器不支持耦合，而是这条路径上的荷载条件让耦合自动失效了。

5. 蠕变分析避坑：五处最容易让结果"看起来对"的设置

5.1 结果差六个数量级：蠕变参数的量纲随 n 和 m 一起变

现象：模型能算完，收敛曲线也漂亮，位移云图的分布形态和预期一致，但数值整体偏大或偏小几个数量级，检查边界条件、网格、荷载都查不出问题。

原因：A 的量纲不是固定的。从ε̇ = A·q^n·t^m反推，A 的量纲是 应变/(应力^n · 时间^m)。实验数据用 MPa 和小时拟合出来的 A，直接搬到用 Pa 和秒建模的模型里，由于 n 通常在 5 附近、m 在 -0.8 附近，换算系数是 (1e6)^5 × (3600)^(-0.8) 这种量级的东西，差出十几个数量级都不奇怪。

解决：建模前先把单位制钉死，写在一张纸上贴在显示器边。拟合脚本里显式注释清楚输入应力和时间的单位，A 的计算结果标注完整量纲。换单位时不要手工缩放 A，直接重跑拟合脚本，让计算机去算那个幂次。

5.2 增量步被砍到 1E-9：CETOL 和最大增量在打架

现象：.sta文件里时间推进极慢，跑了几个小时物理时间才走了总时间的千分之一，.msg里全是 attempt 和 cutback 记录，最后报"time increment has become too small"。

原因：两个方向的力量在互相拉扯。CETOL 判据要求步长足够小才能满足精度，而数据行里的最大增量或者*CONTROLS里的上限又在压制增长，两者的平衡点落到了最小增量以下。另一种情况是初始增量给得太大，第一步就超差，求解器疯狂砍半，砍到极限也没找到能收敛的步长。

解决：先把 CETOL 放宽到 1E-4 看能不能跑通，确认模型本身没有别的问题；再逐步收紧，同时观察增量步分布。如果收紧到 1E-5 就开始卡，说明这个工况的蠕变率变化就是很剧烈，这时候应该检查材料参数是不是拟合得过于陡峭，而不是继续跟容差较劲。

5.3 位移曲线是直线：默认 ramp 幅值把蠕变算成了慢加载

现象：蠕变步的位移-时间曲线是一条漂亮的斜直线，完全没有一次蠕变的减速段，也没有向稳态过渡的弯曲。

原因：ABAQUS 的默认幅值行为是 ramp，荷载在分析步总时长内从初始值线性加到目标值。如果蠕变步里又写了一条*CLOAD，哪怕数值和上一步完全一样，求解器也会认为这是一个从零开始的新加载，用整个 1E6 小时把荷载慢慢加上去。结果算出来的根本不是蠕变，是一个超慢速的静力加载过程。

解决：蠕变步里通常不应该再出现荷载定义，边界条件和荷载会自动从上一个步继承。如果确实需要在蠕变过程中改变荷载，显式写*AMPLITUDE, NAME=INST, DEFINITION=STEP或者在*STEP里加AMPLITUDE=STEP，让荷载瞬时施加，不要占用蠕变时间轴。

5.4 卸载后变形不回去：time hardening 丢了应力历史

现象：模拟一个加载-保载-卸载的循环，保载段的蠕变变形合理，但卸载后的残余变形比实测大出一截，回弹量偏小。

原因：time hardening 形式下，蠕变率只认当前应力和当前时间。卸载瞬间应力降到零，蠕变率也应该归零，但因为历史信息没有被携带，再加载时它会按老的时间轴继续算，等效于把卸载期间的时间白白送给了蠕变。strain hardening 用累积蠕变应变代替时间，这种情况下的表现要好得多。

解决：载荷谱里有明显卸载、反向加载或者应力松弛的工况，直接换LAW=STRAIN。代价是参数要重新拟合，因为 strain hardening 的积分式不能直接从业已积分的曲线族里解出来，通常需要迭代。

5.5 大变形下蠕变率虚高：两个 step 的 NLGEOM 不一致

现象：薄壁结构在长时间保载后鼓胀明显，实测的鼓胀速率比仿真慢，而且仿真后期的蠕变率越来越离谱。

原因：静态步开了 NLGEOM，蠕变步忘了开，或者反过来。蠕变本构里的等效应力用的是真应力，如果蠕变步运行在小变形框架下，求解器会拿工程应力去算蠕变率。鼓胀 10% 时真应力和工程应力差 10%，n 取 5 就是 60% 的率误差，而且这个误差会随变形累积自我放大。

解决：两个 step 的 NLGEOM 必须一致，大变形工况下都开。这不是可选项，是物理正确性的前提。

6. 取舍：解耦近似、容差和算力之间的那条线

这套流程的本质，是一连串有意识的近似叠加出来的结果。静动态解耦丢掉了加载阶段的蠕变应变，Power-law 丢掉了三阶段曲线的完整描述，隐式积分丢掉了瞬态精度换稳定性，CETOL 又在精度和增量步数之间划了一道线。每一处近似单独看都合理，叠在一起之后，最终结果的误差量级需要心里有数。

我一般会这么判断：解耦带来的误差通常是千分之一量级，除非加载时间本身就很长；模型形式带来的误差在百分之一到百分之十，取决于应力变化范围和是否进入三次蠕变；数值容差带来的误差在千分之一到百分之一，CETOL 收紧一档大致能压掉一半。也就是说，真正决定结果可信度的不是容差设多小，而是材料模型选得对不对、参数拟合得准不准。花三天把 CETOL 从 1E-5 调到 1E-7，收益远不如花半天重新检查一遍曲线族的拟合残差分布。

至于 time hardening 和 strain hardening 之间的选择，我的经验是看应力路径的复杂度而不是看教科书推荐。恒温恒载的长期变形评估，time 形式又快又稳，没必要上 strain；一旦载荷谱里出现卸载、保载交替、多轴旋转，strain 形式的历史记忆就是必需项，多花的拟合时间会在结果可信度上还回来。

验证方法上，有个成本很低但很有效的手段：把同一个模型用两套不同的 CETOL 各跑一遍，对比关键节点的位移和应力。如果两条曲线的差异在 5% 以内，说明数值收敛性没问题，剩下的误差都归材料模型；如果差异超过 20%，那问题在数值侧，跟材料参数没关系，先解决收敛再谈精度。

从那以后我每次建蠕变模型，都先把曲线族的单位和拟合残差过一遍，再决定 CETOL 从 1E-4 起步还是直接上 1E-5，最后一定用两档容差对跑一遍才敢出结论。希望帮到你。

本文还有配套的精品资源，点击获取