2.6 弹-黏塑本构模型优选及参数识别分析方法
2.6.1 引言
目前的反分析研究大多是在材料本构模型假定已知的情况下，进行力学参数的反演。如果本构模型不能反映实际情况，参数无论如何精确地确定都难以合理反映实际的力学响应，从而不能起到有效的力学响应预测。由此可见，水利工程反分析中本构模型的优选是一个很重要的课题。
模型优选与参数识别是既有区别又相联系。传统的模型优选是指从具有某种属性的模型集合中优选出相对最佳、最准确的描述系统响应性态的模型。而一个模型要为人们所认识、接受和应用，需将模型以参数的形式具体表述，因而，参数识别又是模型优选的基础。目前，反分析中本构模型的优选主要由两大类组成：一类可称为模型的类内优选，它是指在一类模型中，优选出最佳的模型形式，例如，在弹塑性模型中优选不同的屈服函数；在黏弹性模型中优选不同的黏性方程等。其主要特点为，模型的种类确定，只优选其中的模型结构参数和模型参数。另一类称为模型的类间优选，该法不但优选模型的种类，而且同时优选其中的模型参数。对于实际工程问题，一般根据工程经验确定模型的种类，然后在确定的种类下进行模型结构参数和模型参数的优选（即进行类内优选）。本节介绍弹-黏塑性模型优选和参数识别分析。
2.6.2 弹-黏塑性模型优选基本原理
2.6.2.1 弹-黏塑性模型优选的一般通式
就具体的弹-黏塑性本构模型，在2.3.1节已经进行了概括，这里不再赘述。尽管在黏塑性本构模型方面还需要作进一步的研究，但本节无意提出新的模型，只是基于实测数据，从现有的模型中挑选出一个最佳模型，或从一个通式模型中获得一个最佳模型。本节将适合于连续、均质、各向同性介质的波兹纳本构方程的一般形式作为优选模型的通式，波兹纳方程的一般形式见式（2.4.6b）。该表达式右侧第一项为弹性应变率，第二项为黏塑性应变率。这一公式在形式上的结构式已确定，只要优选参数就可确定出具体模型。在这一模型中主要不同的表现是在黏塑性应变率部分，对未知函数Φ、F和Q的不同选择将表现为不同性质的类型。由于用一组实测数据同时要优选得到3个未知函数和模型的参数向量是比较困难的，这里只就关联流动法则即F≡Q的情形进行优选，这样待定的函数只有Φ和F。其中，F即为屈服准则函数，其通式见式（2.3.6）；屈服函数的函数Φ的常用形式见式（2.4.8）和式（2.4.9）。
2.6.2.2 弹-黏塑性模型优选步骤
由计算经验可知，同时优选的参数越多，容易出现反分析参数的不唯一问题。如果能对待反分析参数进行解耦反演，无疑能得到较合理的反分析值。由2.4.2节的分析可知，弹-黏塑性问题的有限元法也可用于弹塑性有限元分析。进行计算时，施加一个常荷载并对时间进行积分，直至所有的应变率均为零，得到的便是弹塑性解。
本节依据这个规律，提出弹-黏塑性模型优选两步法。
（1）首先利用弹-黏塑性收敛位移及脆性材料应力屈服函数的统一表达式来确定屈服函数的系数值，即优选屈服准则函数F的类型。
（2）然后由时变过程值利用波兹纳本构方程中的常用Φ函数表达式来确定脆性材料的常数N和流动系数γ，即优选屈服函数的函数Φ的类型。
由此可见，脆性材料应力屈服函数F的优选可转化为对屈服函数的统一表达式［式（2.3.6）］中各参数的优选问题。屈服函数的函数Φ的优选则转化为常数N和流动参数γ的优选［式（2.4.9）］。
2.6.2.3 参数估计
参数估计包括模型结构参数的确定和模型参数估计。对弹-黏塑性模型优选而言，参数估计就是对式（2.3.4）中弹性模量E和式（2.3.6）中强度参数α、K、β以及式（2.4.6）的流动系数γ和式（2.4.8）、式（2.4.9）常数N（或α）的估计。按前述分析，强度参数α、K、β及常数N（或α）为模型结构参数。为了求出这些参数的最佳估计，要基于一定的原则建立目标函数。为分析问题方便，这里采用在数据处理中广泛使用的最小二乘法来建立。
（1）弹性模量和屈服函数统一表达式的系数反演的目标函数为
式中：e（i）表示计算值与实测值的误差；Num表示测点数目。
（2）流动系数和常数N（或α）反演的目标函数为
式中：e（i，j）为tj时刻计算值与实测值的误差；n1为测点数，n2为计算中采用的时间段数。
为tj时刻第i个测点的位移实测值；
为tj时刻第i个测点的位移计算值。
适应度函数取为
本节采用遗传模拟退火算法进行弹-黏塑性本构模型优选和参数估计。
2.6.2.4 遗传模拟退火算法
由于遗传算法存在全局搜索能力强而局部搜索能力不足，且容易出现早熟收敛现象等缺点，而模拟退火算法虽然把握搜索过程的总体能力较差、运用效率低，但具有局部搜索能力强等优点。为此，本节将上述两种优化算法结合起来，采用遗传算法进行全局寻优，然后对当代群体中适应度较高的部分优良个体采用模拟退火算法来局部寻优。遗传模拟退火算法的具体步骤如下：
（1）编码：对需反演的参数进行实数编码。
（2）产生初始种群：依据参数取值范围随机生成N个符合初始条件的初始染色体，形成初始种群。
（3）适应度评价：按前述适应度函数计算个体的适应度。
（4）收敛判断：遗传算法收敛判断通常有如下两种方法，方法一是直接以演化代数进行控制，即给定最大演化代数参数，若演化过程达到此最大代数则终止演化；方法二是利用演化过程中最好染色体的变化情况来控制，如果最好染色体连续演化若干代（比如10代）保持不变（在给定精度范围内），则终止演化。方法一简单、直接，但最大演化代数依赖于所求问题的性质和所采用的遗传算法格式，使用中难以准确给定，给大了造成计算浪费，给小了又搜索不到最优解。方法二克服了上述缺点，且对于不太复杂的问题总能获得满意结果，但对于性质复杂的问题或者算法设计不当时，也会出现无休止迭代而收敛不到最优解的情况。鉴于此，采用上述双重演化终止条件，并且给定的最大演化代数相对足够大，以保证算法迭代总能搜索到最优解或近似最优解。若不收敛，则转向（5），若收敛，则转向（7）。
（5）遗传进化操作：采用前述的方法对父代种群进行选择、交叉和变异操作，产生子代种群。
（6）模拟退火算法的局部搜索：分别以当代群体中部分优良个体为初始点进行模拟退火算法的局部搜索，将各优良个体搜索到局部最优解代替当代群体中适应度相对较低的个体，形成新的子代种群；返回（3）循环进行。
（7）结果输出：输出最优个体的各个基因值和目标函数值。
遗传模拟退火算法的程序框图见图2.6.1。程序用Visual Fortran语言编写。
2.6.3 算例分析
以2.5.3节的承受内压作用的厚壁圆筒为例，该厚壁圆筒承受均布内压0.8MPa。有限元网格见图2.5.5。以下按2.6.2节的步骤进行模型优选。
图2.6.1 遗传模拟退火算法的程序框图
（1）优选目的。本算例优选的目的仅是为了说明2.6.2节提出的弹-黏塑性模型优选两步法的可行性以及验证研制程序的合理性。
（2）输入输出数据准备。假设材料服从广义米赛斯条件中内切圆锥的屈服形式［式（2.3.6）和表2.3.1］。为了进行模型优选研究，先进行弹塑性有限元计算，将得到的测点理论计算位移作为实际量测值，测点的位移值见表2.6.1。由于结构对称，测点x向和y向位移基本上满足对称。
表2.6.1 测点位移
（3）模型类属的确定。假设该模型为弹-黏塑性模型，即进行模型的类内优选。按2.6.2节的原理，本算例仅优选屈服准则函数和屈服函数的函数的类型，屈服准则函数的通式见式（2.3.6）；屈服函数的函数常用表达式见式（2.4.8）和式（2.4.9）。为分析问题方便，在进行屈服函数优选时，设定模型结构参数β=0，从而，仅需优选α和K两个模型结构参数（或称屈服函数统一表达式的系数，或简称屈服函数系数）；屈服函数的函数采用幂指数表达式［式（2.4.9）］，即需优选流动系数γ和模型结构参数N（或称常数N）。泊松比的参数敏感性较低，采用原值，不进行反分析。
（4）模型优选时数量级的同级化处理。为了避免模型优选时由于量纲差异很大而导致误差，本算例将模型结构参数和模型参数的数值转化为0～10之间，即弹性模量的单位取为104MPa，模型结构参数K的单位取为MPa，黏性系数的单位取为103MPa·s。设定各参数的取值范围，弹性模量E∈［2.5，3.5］×104MPa；模型结构参数0≤α≤0.5，0<K≤2.5MPa，N∈［1.0，4.0］；黏性系数η=1/γ∈［0.1，5.0］×103MPa·s。
（5）参数估计。设定优化反分析程序各控制参数的取值为：种群规模40，进化代数为100，进化的终止误差阀值为1×10-5。模拟退火算法中部分优良个体数取为10，接受次数取为30，步长上限为10，Tf=1×10-4，循环总数为10。采用遗传模拟退火算法进行参数估计，结果见表2.6.2。
表2.6.2 弹性模量和屈服函数系数优选值
确定弹性参数和屈服函数统一表达式的系数后，接下来利用黏塑性的过程值来优选黏性系数η和常数N。反分析的各控制参数的取值同前。弹-黏塑性计算采用隐式算法，累计计算时间90d（7.776×106s），时间步长为1d（86400s）。优选结果为η=1.1488×103MPa·s，N=1.1974。
（6）模型检验。将上述模型优选的参数值进行正分析可知，计算过程值与理论过程值拟合的精度高。这一方面是由于本节提出的模型优选两步法：参数合理解耦、逐步反演，反分析效果好；另一方面也可见，采用遗传模拟退火算法进行模型优选和参数估计有一定的优势，能一定程度减弱反分析的不唯一性。