CN115828762A - 受弯梁材料本构参数模型修正及判别方法 - Google Patents

Abstract

本发明公开了一种受弯梁材料本构参数模型修正及判别方法。步骤如下：选择合适的材料本构方程；选取材料本构方程中的待修正参数，并确定合理的待修正参数取值区间；基于平截面假定，由截面应力平衡列轴力平衡方程，求解截面受压区高度；判断轴力平衡方程中是否存在线性相关项；由力矩平衡列弯矩平衡方程，计算截面抗弯刚度，截面受压区高度为常数时，判别式EIB≠0，能进行材料本构参数模型修正；截面受压区高度不为常数时，判断式与参数个数相关；采用理论模型或者响应面模型通过优化算法对待修正参数进行识别，完成模型修正。本发明提出一种受弯梁材料本构参数模型修正及判别方法，为非线性模型修正方法在实际工程中的应用提供了理论基础。

Description

受弯梁材料本构参数模型修正及判别方法

技术领域

本发明属于结构模型修正领域，涉及梁结构本构模型参数修正技术，具体地说是涉及一种受弯梁材料本构参数模型修正及判别方法。

背景技术

模型修正是结构健康监测系统的重要组成部分，模型修正技术能有效减小模型与实际结构的误差，在土木工程领域中应用广泛，也是目前比较热门的研究方向。已有模型修正研究多以构件刚度参数为修正量，从材料本构模型的层面对本构参数进行修正的研究甚少。本发明开展混凝土梁桥静力响应数值模拟，建立本构关系参数与结构静力响应值(应变、刚度、挠度)的响应面，提出修正材料本构关系参数的方法，从本构层面进行有限元模型修正，识别本构参数，同时提出能否进行非线性本构模型参数修正的判据。

发明内容

针对受弯梁能否进行非线性模型修正，识别材料本构参数的问题，本发明提出一种受弯梁材料本构参数模型修正及判别方法。

本发明所述受弯梁材料本构参数模型修正及判别方法，步骤如下：

(1)根据工程经验，选择合适的材料本构方程；

(2)选取材料本构方程中的待修正参数，并确定合理的待修正参数取值区间，a)仅一个待修正参数时，可进行材料本构参数模型修正，转步骤(6)；b)两个及以上待修正参数时，继续下述步骤；

(3)基于平截面假定，假设截面曲率为

由截面应力平衡，列轴力平衡方程，求解截面受压区高度

表示关于自变量

的函数；

(4)以待修正参数为变量，判断轴力平衡方程中是否存在线性相关项，若存在，则不能进行材料本构参数模型修正，转步骤(7)；若不存在，继续下述步骤；

(5)由力矩平衡，列弯矩平衡方程，求得弯矩M，计算截面抗弯刚度

若对多个不同截面曲率

k＝1,2,3...，a)截面受压区高度

为常数时：若判别式EIB＝0，不能进行材料本构参数模型修正，转步骤(7)；判别式EIB≠0，能进行材料本构参数模型修正，继续下述步骤；

b)截面受压区高度

不为常数时：两个待修正参数，判别式EIA＞0.4％，三个待修正参数，判别式EIR＞0.7％，四个待修正参数，判别式EIR＞1.5％，能进行材料本构参数模型修正，继续下述步骤；否则，不能进行材料本构参数模型修正，转步骤(7)；

(6)采用理论模型或者响应面模型，通过优化算法对待修正参数进行识别，完成材料本构参数的修正，模型修正时目标响应值的个数不小于待修正参数的个数；

(7)模型修正结束；

具体的，步骤(5)中，判别式分别为：

其中，带下标的EI表示对不同的

值，通过理论模型计算得到的截面抗弯刚度，

EI为响应值，n表示待修正参数的数量，下标k表示顺序号，k＝1,2,3...,n-1；

EI_k：表示在选取的多组EI响应值中，按照曲率由小到大排序，曲率为

时对应计算得到的截面抗弯刚度EI_k；

EI_k+1：表示在选取的多EI组响应值中，按照曲率由小到大排序，曲率为

时对应计算得到的截面抗弯刚度EI_k+1，与EI_k相邻；

EI_max：表示在选取的多EI组响应值中，按照曲率由小到大排序，最大的曲率值对应计算得到的响应值；

EI_min：表示在选取的多组EI响应值中，按照曲率由小到大排序，最小的曲率值对应计算得到的响应值。

具体的，步骤(1)中材料本构方程可根据相应的设计规范选择，如混凝土、普通钢筋、预应力钢筋本构模型可参照《JTG 3362-2018公路钢筋混凝土及预应力混凝土桥涵设计规范》选用。

具体的，步骤(2)中材料本构方程中的待修正参数取值区间要包含待修正参数的实际值。

具体的，步骤(3)中，轴力平衡方程与梁截面形状、组成截面材料的种类数量、材料本构方程相关，以单一材料矩形截面梁为例：

∑N＝0；

其中，∑N为截面轴力之和，N表示轴力，B、H分别为矩形截面的宽度和高度，z为截面受压区高度，

分别为材料的受压应力、受拉应力，由材料本构方程计算，

为距离截面中性轴y位置的应变，

为截面曲率，y为距离截面中性轴的位置；

假设截面曲率

已知，可求得应变

通过材料本构方程，可求得受压应力

受拉应力

轴力平衡方程中只有截面受压区高度z未知，故可求解截面受压区高度

表示关于自变量

的函数。

具体的，步骤(5)中，以单一材料矩形截面梁为例：

M_i＝∑M；

其中，∑M为截面弯矩之和，M表示弯矩，M_i为与截面曲率

对应的总弯矩。

具体的，步骤(6)中，理论模型指采用数学方法建立的精确理论数学模型，响应面模型指由试验设计方法生成样本，拟合得到响应与待修正参数之间的响应面模型，是近似模型，与理论模型相比，有一定的误差，响应面模型可采用多项式模型、克里格模型、径向基函数、BP神经网络等模型。

具体的，步骤(6)中，模型修正的目标响应值取截面抗弯刚度EI、梁的挠度w或者梁的挠度倒数w^-1等响应。

具体的，步骤(6)中，模型修正可采用单目标或多目标优化算法。

本发明基于理论模型和响应面模型，以截面抗弯刚度EI、梁的挠度w或者梁的挠度倒数w^-1等响应为目标响应，进行非线性响应面模型修正，分析了多种材料本构方程，总结了受弯梁能否进行非线性模型修正的判别方法，为非线性模型修正方法在实际工程中的应用提供了理论基础。

附图说明

图1是本发明方法的流程框图。

图2是本发明方法的理论模型的参数修正步骤。

图3是本发明方法的矩形截面梁受弯分析。

图4是本发明方法的响应面模型的参数修正步骤。

图5是本发明方法的简支梁示意图(单位:mm)。

图6是本发明方法的梁截面示意图(单位:mm)。

图7是本发明方法的受弯状态下材料本构模型参数修正的判别方法示意图。

图8是本发明方法实施例一的结构设计布置示意图(单位:mm)。

图9是本发明方法实施例二的结构设计布置示意图(单位:mm)。

具体实施方式

下面结合附图及实施例对本发明做进一步说明，下面的描述涉及附图时，除非另有表示，不同附图的相同数字表示相同或相似的要素。

本发明所述受弯梁材料本构参数模型修正及判别方法，实现流程如图1，具体内容如下：

一、理论模型及响应面模型的材料本构参数修正方法

1)梁结构基于理论模型的材料本构参数修正方法(EI响应)

本节算例以UHPC本构(上升段)为例进行数值模拟分析，建立简支梁数值模型，以EI为响应采用基于理论模型的参数修正方法对本构参数进行修正。

数值算例采用简支梁(受弯结构)模型，通过数值算例对本构参数修正方法进行详细介绍，应用两种方法实现相关参数修正，并比较两种方法的修正效果。简支梁为不配筋的超高性能混凝土梁，矩形截面(B×H：200mm×400mm)，计算跨径为4100mm。

为分析方便算例采用本构模型的上升段(极限应变之前)作为分析本构模型，如下式：

式中

式中f_c、f_t分别为混凝土的抗压、抗拉强度值；ε₀为混凝土压应力刚达到f_c时混凝土的应变值，即峰值应变；

表示梁结构受弯时梁轴线变形后的曲率，y_i表示距中性轴的距离。

通过上式可知：受压本构主要由三个参数决定，分别是f_c、ε_c、ε_c0；受压本构主要也是由三个参数决定，分别是f_t、ε_t、ε_t0；在本算例中将f_c和f_t选取为待修正参数。

①方法介绍及修正步骤

本节方法建立响应与本构参数之间的理论模型，利用数学中的最优化问题去解决工程应用中的参数修正问题，同时需要运用MATLAB软件辅助计算来实现。

以简支梁桥结构作为分析对象，该方法首先根据结构的材料本构关系选取待修正参数，然后根据平截面假定及结构的受力状态推导出截面的轴力平衡及弯矩平衡方程，以EI为目标响应值建立理论模型，然后根据结构设计已知条件及平衡方程运用MATLAB软件编程或非线性有限元分析软件对结构进行响应计算，最后根据理论模型运用MATLAB最优化问题求解办法，求解出待修正参数的修正值并进行误差分析。具体修正步骤如图2所示。

②理论分析

根据平截面假定以及本构模型，可得到本构上升段时截面的应力、应变分布图，如图3。

根据平衡条件，设截面受压区高度为z，距截面中性轴距离为y的位置，由静力关系可以写出两个平衡方程：

将本构方程代入式(2)计算可得轴力平衡方程：

将本构方程代入式(3)计算可得弯矩平衡方程：

通过式(4)可知受压区高度满足以下函数关系：

根据材料力学可知有以下关系，式中EI表示截面抗弯刚度：

那么通过式(4)、(5)和(7)可知EI需满足以下函数关系：

由以上推导，根据最优化问题(无约束最优化问题)，当EI_0i为已知目标响应值时，EI_i是根据以上推导得到的函数计算值，在参数给定的参数区间内求解参数f_c和f_t的修正值，可以建立一个关于EI_i的数学理论模型，如式(9)，目标函数采用平方和最小，当目标函数最小时说明EI_i与EI_0i无限接近，而此时EI_i是对应的参数值f_c和f_t通过理论计算得到，则此时的f_c和f_t就是拟合得到的参数的修正值，此拟合步骤要应用MATLAB实现相关计算。

式中：EI_0i为已知目标响应值，EI_i为计算响应值。

同时可以观察最后输出目标函数的数值检验模型的精度，当目标函数f越接近于0说明建立的理论模型的吻合度越高。

③试验数值计算

本试验设计两组理论试验，A组试验两参数的取值接近，模拟拉压本构参数接近的材料；B组试验两参数的取值相差比较大，模拟拉压本构参数相差大的材料(即模拟UHPC材料的拉压强度性能)。先求解两组EI响应，然后根据响应值结合理论模型求解参数的修正值：

(1)A组试验：以

为已知量，峰值应变ε₀＝0.002，给f_c、f_t一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，30]，通过两平衡方程进行计算分别计算出响应值EI₁、EI₂和EI₃。

(2)B组试验：以

为已知量，峰值应变ε₀＝0.002，给f_c、f_t一个包含理论值的区间为f_c∈[120，270]、f_t∈[10，30]，通过两平衡方程进行计算分别计算出响应值EI₁、EI₂和EI₃。

A组中取一组理论值f_c＝28MPa、f_t＝21MPa，B组中取一组理论值f_c＝150MPa、f_t＝13MPa，应用MATLAB计算，截面受压区高度z的计算结果见表1，弯矩M的计算结果见表2，响应值EI_i的计算结果见表3。

表1两组试验在不同

值时计算得出的受压区高度z(单位：mm)

表2两组试验在不同

值时计算得出的弯矩M(单位：kN·m)

表3两组试验在不同

值时计算得出的EI_i目标响应值(单位：N·mm²)

④参数灵敏度分析

灵敏度分析采用局部灵敏度分析，对某个参数变化增幅1％，比较响应值的变化，分析判断其显著性。两组试验中各参数的局部灵敏度计算结果如下表4和表5，表中响应为EI_i(单位为：N·mm²)。

表4试验A组f_c＝28MPa，f_t＝21MPa参数的局部灵敏度分析

表5试验B组f_c＝150MPa，f_t＝13MPa参数的局部灵敏度分析

通过灵敏度分析，由表4、表5可知：在两组试验中，f_t的灵敏度均大于f_c的灵敏度，初步判断在应用该方法修正参数时，f_t参数比f_c参数对EI响应计算的影响更显著，而f_t的修正效果比f_c的修正效果较好。

⑤优化求解及效果

根据建立的理论模型，通过前面计算的响应值优化拟合计算得到f_c、f_t的修正解，比较修正解与理论解的误差。输入目标响应EI_i(三个或任意两个)，计算结果见表6。

表6两组试验参数修正的计算结果表

两组试验理论模型修正参数时输出的目标函数皆符合f＝0；计算结果显示，经过理论模型对参数修正，参数修正结果非常好，参数相对误差均为0，修正值为精确解。综合说明，对梁结构的基于理论模型的本构材料参数修正方法是有效的，能够实现对材料参数进行精确修正，得到精确解。

2)梁结构基于响应面模型的材料本构参数修正方法(EI响应)

本节同样以UHPC本构(上升段)为例，建立简支梁数值模型，以EI作为响应值，采用基于响应面模型的参数修正方法对本构参数进行修正，详细介绍该方法并阐述该方法的有效性及修正效果。

①方法介绍及修正步骤

响应面法在模型修正中应用得相当广泛，本节该修正方法就是建立响应与本构参数之间的响应面模型，实现解决数值算例及工程应用中的参数修正问题。该方法首先根据设计结构的材料本构关系选取待修正参数，同时确定待修正参数的取值区间，然后选择试验设计方法进行样本设计，接着根据已有条件对结构计算及试验样本计算，再然后确定目标函数、以EI或w为目标响应值建立响应面模型，同时对模型进行评估，最后根据响应面模型优化求解出待修正参数的修正值并进行误差分析。结构计算需要运用MATLAB软件或有限元分析软件，求解修正值过程需要运用MATLAB软件和Design-Expert 12软件辅助计算来实现。具体修正步骤如图4所示。

②试验设计

根据本构模型，将f_c和f_t选取为待修正参数，同样设计两组试验，采用CCD实验设计方法进行试验设计，由两水平两因素、采用一个中心点，生成9组样本序列(通过Design-Expert12软件实现)：

(1)A组试验样本：理论值取f_c＝28MPa、f_t＝21MPa，给f_c、f_t一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，30]；参数水平变量表见表7所示、CCD实验设计表见表8所示。

(2)B组试验样本，理论值取f_c＝150MPa、f_t＝13MPa，给f_c、f_t一个包含理论值的区间为f_c∈[120，270]、f_t∈[10，30]。参数水平变量表见表9所示、CCD实验设计表见表10所示。

表7A组参数水平变量表

表8A组CCD实验设计表

表9B组参数水平变量表

表10B组CCD实验设计表

③试验数值(EI响应)计算

根据所拟定简支梁模型的截面尺寸及跨径，以

为已知量，试验设计以EI作为响应值建立响应面模型做参数修正分析，对9组样本序列通过本构关系计算结构的各响应值，EI_i响应值的计算可用MATLAB软件或有限元分析软件进行计算，设计两组理论试验的9组样本序列EI_i响应计算结果如表11、表12。

表11A组试验样本

值不同时的EI_i响应值(N·mm²)

表12B组试验样本

值不同时的EI_i响应值(N·mm²)

④参数灵敏度分析

对试验样本的响应值做方差分析(应用Design-Expert 12软件)，通过软件得到F检测法各参数一次项、二次项和交叉项的F值和P值，进而判断各参数的灵敏度及显著性。

两组试验样本序列，EI_i作为响应值时，各参数项的F、P值如表13、表14所示。

通过方差分析，由表13、表14可知：在A组试验中，f_c与f_t的显著性相当，初步判断在应用该方法修正参数f_c、f_t时，两参数修正效果相当；在B组试验中，f_t的显著性大于f_c的显著性，初步判断在应用该方法修正参数f_c、f_t时，f_t的修正效果比f_c的修正效果较好。

表13A组试验样本EI_i为响应值各参数项的F、P值

表14B组试验样本EI_i为响应值各参数项的F、P值

⑤刚度响应(EI)与参数的响应面拟合

响应面模型本节所应用的就是二次多项式(含交叉项)模型来拟合，形式如式(10)，式中y为响应特征值，x_i(i＝1,2,…n)为待修正参数，α为待定系数，待定系数应用最小二乘法回归分析求得。

A组试验中f_c、f_t与EI_i的二次多项式响应面方程如下：

B组试验中f_c、f_t与EI_i的二次多项式响应面方程如下：

对拟合的响应面模型需要进行精度检测，精度检测用来评价模型与实际有限元模型的接近程度，进行精度检测时一般采用R-squared和Adjusted R-squared两个指标。R-squared(即R²)描述的变量之间相关程度的指标量；Adjusted R-squared(即Adj R²)是考虑了变量数目的影响根据样本数量和变量数量进行调整后的R²指标量。两个判定系数取值范围都为[0，1]，越接近1表明拟合的计算模型在试验样本空间内越能描述变量之间的相关程度，具体计算公式如下。

式中，y_R(i)：表示模型(响应面模型)的某次计算值；y(i)：表示理论(有限元分析)的某次计算值；y：表示理论(有限元分析)计算结果的平均值；NN：表示试验样本空间的数量；n：表示变量(待修正参数)的数量。

A、B两组试验样本空间R²和Adj R²的计算结果见下表15响应面精度检测表：

表15A、B两组试验样本空间的响应面精度检测表

由上表可知，两组试验的R²和Adj R²的值都接近于1，说明拟合的响应面模型在试验样本空间内能代替有限元模型，能较准确的描述响应值和参数之间的相关关系。

⑥刚度响应(EI)优化求解及效果

根据建立的响应面模型，通过计算理论组的EI目标响应值，输入目标响应EI_i(三个或任意两个)，应用Design-Expert 12软件优化拟合求解计算f_c、f_t的修正解，比较修正解与理论解的误差。

(1)A组试验样本：理论组取f_c＝28MPa、f_t＝21MPa，输入计算得到对应的EI响应值，通过优化拟合计算得到修正解，计算结果见表16。

表16A组试验：f_c＝28MPa，f_t＝21MPa参数修正的计算结果表

表17B组试验：f_c＝150MPa，f_t＝13MPa参数修正的计算结果表

(2)B组试验样本，理论组取f_c＝150MPa、f_t＝13MPa，输入计算得到对应的EI响应值，通过优化拟合计算得到修正解，计算结果见表17。

由表16、表17表计算结果显示，经过响应面模型修正法对参数拟合修正，两组试验参数修正结果都较好，相对误差控制在5％以下，其中A组相对误差控制在1％以下。综合说明，对梁结构的基于响应面模型的本构材料参数修正方法是有效的，能够实现对材料参数进行较精确修正，求出修正解。同时，在A组试验中，f_c与f_t两参数修正效果相当；在B组试验中，f_t的修正效果比f_c的修正效果较好，两组修正结果对比说明参数的修正效果与参数显著性有直接联系。

二、结构受弯状态下材料本构模型参数修正及判别方法

不同本构组算例采用简支梁(受弯状态)作为分析对象，尺寸拟定同上节，如图5所示，跨径为4100mm。截面都为矩形截面(200mm×400mm)，根据需要分为三种形式：

(1)结构中不配筋，截面如图6(1)所示；

(2)结构只在截面受拉侧配受力筋，采用钢筋类型为HRB335带肋钢筋、公称直径为d₁，截面配筋面积为As，保护层厚度为t，截面如图6(2)所示；

(3)结构在截面受拉、受压双侧配受力筋，采用钢筋类型为HRB335带肋钢筋，受拉侧配筋公称直径为d₁、截面配筋面积为As₁、保护层厚度为t₁，受压侧配筋公称直径为d₂、截面配筋面积为As₂、保护层厚度为t₂，截面如图6(3)所示。

1)本构模型数值算例理论分析及响应计算

本构组一

本构组一假设由单种材料组成，其中拉压本构模型完全相同，均为非线性本构，选取两个待修正参数；假设的材料本构模型形式如式(19)，结构为受弯状态则式中

式中

根据平截面假定及平衡条件，设截面受压区高度为z，由前面推导的静力关系式(2)、(3)可得到的两个平衡关系如下：

联立式(20)、(21)进行化简并求解z、EI：

经计算，可得到z的三个解如式(23)，但有效解只有一个即H/2，由此可判断z与

值的变化无关，而由式(22)可判断EI会随

值得变化而改变。

采用图6(1)中拟定的不配筋截面，梁结构不配筋，选取f_c、f_t为待修正参数，设计一组理论值：f_c＝28MPa、f_t＝21MPa，峰值应变ε₀＝0.002，对模型进行理论分析计算。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应的计算结果见表18。

表18试验组在不同

值时计算得出的z_i、M_i、EI_i

本构组二

本构组二假设由两种材料组成，两种材料拉压本构模型相同，均为线性本构，选取三个待修正参数；假设的材料本构模型形式如式(24)，结构为受弯状态则式中

式中

根据平截面假定及平衡条件，设截面受压区高度为z，由前面推导的静力关系式可得到的两个平衡关系如下：

对式(25)进行化简，式(26)已是最简形式：

联立式(26)、(27)进行化简并求解z、EI，经计算z为常数，EI线性相关：

采用图6(2)中拟定的拉侧配筋截面，在梁结构受拉侧配2根受力钢筋，带肋钢筋公称直径d＝22mm，截面配筋面积As＝760.265mm²，选取f_c、f_t、f_s为待修正参数，设计一个理论值：f_c＝30MPa、f_t＝20MPa、f_s＝25MPa，峰值应变ε₀＝0.002，对模型进行理论分析计算。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应的计算结果见表19。

表19试验组在不同

值时计算得出的z_i、M_i、EI_i

本构组三

本构组三假设由单种材料组成，其中拉压本构模型形式上相同，表达形式中的参数不同，均为非线性本构，选取两个待修正参数；假设的材料本构模型形式如式(29)，结构为受弯状态则式中

式中

根据平截面假定及平衡条件，设截面受压区高度为z，由静力关系推导可得到的两个平衡关系如下：

联立式(30)、(31)进行化简并求解z、EI：

采用图6(1)中拟定的不配筋截面，梁结构不配筋，选取f_c、f_t为待修正参数，设计一个理论值：f_c＝28MPa、f_t＝21MPa，峰值应变ε₀＝0.002，对模型进行理论分析计算。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应。

本构组四

本构组四假设由单种材料组成，其中拉压本构模型不同，受压为非线性本构、受拉为线性本构，选取两个待修正参数；假设的材料本构模型形式如下式，结构为受弯状态则式中

式中

根据平截面假定及平衡条件，设截面受压区高度为z，由静力关系推导可得到的两个平衡关系同式(4)、(5)，联立进行化简并求解z、EI：

采用图6(1)中拟定的不配筋截面，梁结构不配筋，选取f_c、f_t为待修正参数，设计一个理论值：f_c＝30MPa、f_t＝21MPa，峰值应变ε₀＝0.002，对模型进行理论分析计算。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应。

本构组五

本构组五假设由两种材料组成，第一种材料拉压本构模型不同，受压为非线性本构、受拉为线性本构，第二种材料以截面受拉侧配筋体现为线性本构，选取三个待修正参数；假设的材料本构模型形式如下式，结构为受弯状态则式中

式中

根据平截面假定及平衡条件，设截面受压区高度为z，由静力关系推导可得到的两个平衡关系如下：

联立两式进行化简并求解z、EI：

采用图6(2)中拟定的拉侧配筋截面，在梁结构受拉侧配2根受力钢筋，带肋钢筋公称直径d＝22mm，截面配筋面积As＝760.265mm²，选取f_c、f_t、f_s为待修正参数，设计一个理论值：f_c＝30MPa、f_t＝20MPa、f_s＝25MPa，峰值应变ε₀＝0.002，对模型进行理论分析计算。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应。

本构组六

本构组六假设由两种材料组成，第一种材料拉压本构模型不同、但均为非线性本构，第二种材料以截面拉压两侧配筋体现、均为线性本构、但表达形式中参数不同，选取四个待修正参数；假设的材料本构模型形式如下式，结构为受弯状态则式中

式中

根据平截面假定及平衡条件，设截面受压区高度为z，由静力关系推导可得到的两个平衡关系如下：

联立两式进行化简并求解z、EI：

采用图6(3)中拟定的拉压双侧配筋截面，在梁结构受拉、受压双侧配筋，在受拉侧配2根受力钢筋，带肋钢筋公称直径d₁＝22mm，截面配筋面积As₁＝760.265mm²；在受压侧配2根受力钢筋，带肋钢筋公称直径d₂＝10mm，截面配筋面积As₂＝157.08mm²；选取f_c、f_t、E_s1、E_s2为待修正参数，设计一个理论值：f_c＝26.64MPa、f_t＝2.259MPa、E_s1＝1.8×10⁵MPa、E_s2＝2×10⁵MPa，峰值应变ε₀＝0.002，对模型进行理论分析计算。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应。

2)本构模型数值算例不同工况下的参数修正结果

基于理论模型对受弯梁材料本构模型参数进行修正，采用多种判别式对本构模型参数修正进行分析探讨。结合非线性度的概念，采用1-范数对响应进行分析计算得到非线性度判别式，如式(43)所示，式中ψ表示非线性度，N(EI)表示非线性系统的EI响应值，L(EI)表示线性系统的EI响应值。经计算分析发现，用非线性度判别式对修正效果进行判别时，判别效果不明显，无法准确分类判断。同时，还结合范德蒙行列式等提出多种判别式，但判别效果均不明显，因此均未采用。

最终提出三个对EI响应进行分析计算的判别公式，计算方法如式(44)～(46)，式(44)中EIB指标适用于受压区高度为定值情况，式(45)中EIA指标适用于受压区高度不为定值且待修正参数为两个时的情况，式(46)中EIR指标适用于受压区高度不为定值且待修正参数≥3时的情况：

其中，带下标的EI表示对不同的

值，通过理论模型计算得到的截面抗弯刚度，

EI为响应值，n表示待修正参数的数量，下标k表示顺序号，k＝1,2,3...,n-1。

EI_k：表示在选取的多组EI响应值中，按照曲率由小到大排序，曲率为

时对应计算得到的截面抗弯刚度EI_k。

EI_k+1：表示在选取的多组EI响应值中，按照曲率由小到大排序，曲率为

时对应计算得到的截面抗弯刚度EI_k+1，与EI_k相邻。

EI_max：表示在选取的多组EI响应值中，按照曲率由小到大排序，最大的曲率值对应计算得到的响应值。

EI_min：表示在选取的多组EI响应值中，按照曲率由小到大排序，最小的曲率值对应计算得到的响应值。

本构组一

给f_c、f_t一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，40]，对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取两个不同

值

对应的[EI_k,EI_k+1]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如下表20所示。

表20部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

本构组二

给f_c、f_t、f_s一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，40]、f_s∈[10，40]，对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取三个不同

值

对应的响应值[EI_k,EI_k+1,EI_k+2]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如下表所示。

表21部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

本构组三

给f_c、f_t一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，40]，对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取两个不同

值

对应的响应值[EI_k,EI_k+1]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如表22所示。

表22部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

本构组四

给f_c、f_t一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，40]，对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取两个不同

值

对应的响应值[EI_k,EI_k+1]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如表23所示。

表23部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

本构组五

给f_c、f_t、f_s一个包含理论值的区间为f_c∈[10，40]、f_t∈[10，40]、f_s∈[10，40]，对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取三个不同

值

对应的响应值[EI_k,EI_k+1,EI_k+2]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如表24所示。

表24部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

本构组六

给f_c、f_t、E_s1、E_s2一个包含理论值的区间为f_c∈[15，30]、f_t∈[1.75，2.5]、E_s1∈[160000，210000]、E_s2∈[160000，210000]，对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取四个不同

值

对应的[EI_k,EI_k+1,EI_k+2,EI_k+3]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如表25所示。

表25部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

3)结构受弯状态下材料本构模型参数修正的判别方法

前文列举了多组表达形式不同的本构模型，经过理论推导及分析计算，对各组本构模型不同测点工况下的本构参数进行修正并对修正效果进行评判，分析发现：

通过轴力、弯矩平衡关系计算可知，本构组一和二受压区高度z恒为定值。本构组二无论何种工况修正效果都极差，受压区高度z为恒定值C，且轴力平衡方程中存在线性相关项，进一步EI不随

值的变化而变化，因此该本构不能进行本构参数修正。本构组一只要满足EI不取相同的响应值各工况修正效果皆良好，受压区高度z恒等于H/2，不随

值的变化而变化，但EI随

值的变化而变化，采用判别式(44)进行分析，说明受压区高度z为恒值时，判别式满足EIB≠0，则该本构能够进行本构参数修正。该法则同样适用于本构二。

同样通过两平衡公式计算，其他本构组皆满足受压区高度z不为恒值，但修正效果参差不齐。经大量分析计算，对结构中各材料本构进行理论分析，若平衡方程不存在线性相关项时，且受压区高度z不为恒值时，若该本构可进行本构参数修正，且修正效果良好，则根据判别式(45)、(46)计算需满足以下条件：当待修正参数为两个时，满足判别式(阀值)EIA＞0.4％，如本构组三和本构组四；当待修正参数为三个时，满足判别式(阀值)EIR＞0.7％，如本构组五；当待修正参数为四个时，满足判别式(阀值)EIR＞1.5％，如本构组六。若确定本构能够进行参数修正时，解的特征是唯一解还是双解，需要由轴力、弯矩平衡方程及参数灵敏度进一步确定，如本构组三的修正解为双解。

根据以上分析对其他大量本构组进行理论分析计算并对参数进行修正，发现均能满足以上分析中给定的阀值，因此由以上的算例分析，针对结构受弯状态下材料本构模型能否进行参数修正提出的判别方法，具体流程如图7所示。

具体的，步骤6中，响应面模型还可采用克里格模型、径向基函数、BP神经网络等。

具体的，步骤6中，模型修正的目标响应值还可取梁的挠度w或者梁的挠度倒数w^-1等响应。

具体的，步骤6中，模型修正可采用单目标或多目标优化算法。

实施例一：素混凝土梁本构模型参数修正

混凝土强度等级为C45，截面为300mm×450mm的矩形截面，跨径为5100mm，结构构造设计如图8所示。进行EI响应计算时，测点选择曲率

为0.1×10^-6mm^-1、0.2×10^-6mm^-1、0.4×10^-6mm^-1时；进行挠度w响应计算时，加载方式为跨中集中力加载(P＝5kN)，测点位置为曲率对应的跨中位置。

本算例模拟素混凝土材料梁结构，先判别该材料本构能否进行参数修正，然后采用修正方法对强度参数(即本构模型参数)进行修正以及评判其修正效果。

步骤1：选择材料本构方程，在本算例中，混凝土受压本构采用Rüsch建议的非线性本构公式，混凝土受拉本构采用规范中建议的非线性本构公式，如式(47)。

式中

步骤2：选择f_c和f_t两个参数作为待修正参数，以EI_i和w_i作为响应值，混凝土峰值应变ε_c0＝0.002、ε_t0＝0.0001，给f_c、f_t一个包含理论组的区间为f_c∈[15，50]MPa、f_t∈[1.50，3.10]MPa。

步骤3：根据平截面假定，设截面受压区高度为z，列轴力平衡方程：

求解受压区高度：

步骤4：由式(48)可以判断待修正参数f_c、f_t线性无关。

步骤5：由力矩平衡，列弯矩平衡方程：

可求得截面抗弯刚度：

根据C45混凝土拉、压强度标准值为f_c＝29.6MPa、f_t＝2.51MPa，模拟理论组(模拟假设结构使用若干年后混凝土强度值下降10％的工况)拉、压强度取值为f_c＝26.64MPa、f_t＝2.259MPa。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应的计算结果见表26。

表26试验组在不同

值时计算得出的z_i、M_i、EI_i

对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取两个不同

值

对应的响应值[EI_k,EI_k+1]作为响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经分析计算，选取部分修正效果区别显著的修正组，如下表27所示。

表27部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

该本构组属于两参数型且受压区高度z≠C，根据判别方法需满足判别式EIA值大于0.4％，以上分析计算皆满足。而本算例中，对选定的

工况进行EIA值计算，由式(45)计算可得：

在该算例的条件下满足判别式EIA值大于0.4％的要求，所以在该工况下能对本构进行参数修正。

步骤6：模型修正

①理论模型优化

以

为已知量，应用MATLAB或有限元分析软件进行结构计算，截面受压区高度z的计算结果见表28。

表28在不同

值时计算得出的受压区高度z(单位:mm)

弯矩M的计算结果见表29，响应值EI_i的计算结果见表30。

表29在不同

值时计算得出的弯矩M(单位:kN·m)

表30在不同

值时计算得出的EI_i响应值(单位:N·mm²)

根据无约束最优化问题，根据建立的理论模型，通过目标响应值优化拟合计算得到f_c、f_t的修正解，比较修正解与理论解的误差。输入理论值对应的响应EI_i，求解参数修正解，计算结果见表31，修正结果准确。

表31参数修正的计算结果表

②响应面方法修正(w响应)

试验设计及样本序列同前，根据拟定的结构，采用

非线性分析程序对结构进行挠度w响应计算。以

为已知量求出的M_i值映射到梁对应的测点位置为跨中0.128169L、0.25565L、0.5L处，计算得到三组挠度目标响应，响应计算结果如表32所示。

表32试验样本不同测点的w_i响应值(mm)

试验中f_c、f_t与w_i的二次多项式响应面方程形如式(53)，A、B、C、D、E表示各参数项对应的系数，三个响应面方程各参数项对应的系数值如表33所示：

对拟合得到的各响应面进行精度检测，如表34所示。

理论组取f_c＝26.64MPa、f_t＝2.259MPa，输入计算得到对应的w目标响应值，通过优化拟合计算得到修正解，修正计算结果见表35，计算可知参数f_t的修正误差较大有7％。

表33挠度w响应面方程各参数项对应的系数值

表34参数f_c、f_t与响w_i的响应面精度检测表

表35算例w响应参数修正计算结果表

由于直接采用w响应进行参数修正的误差较大，而采用EI响应进行参数修正的效果比较好，根据梁结构挠度w的计算方法式(54)可近似看出w与EI成反比，因此可以尝试着采用w^-1作为响应进行参数修正。

在本算例采用挠度w^-1响应进行修正计算，验证该猜想，选用跨中0.128169L、0.25565L、0.5L处的三组挠度的-1次方作为响应，响应计算结果如下表36所示：

表36试验样本不同测点的w_i ^-1响应值(mm)

试验中f_c、f_t与w_i ^-1的二次多项式响应面方程各参数项对应的系数值如表37所示：

表37挠度w^-1响应面方程各参数项对应的系数值

对拟合得到的各响应面进行精度检测，如表38所示：

表38参数f_c、f_t与响w_i ^-1的响应面精度检测表

表39算例w^-1响应参数修正计算结果表

理论组取f_c＝26.64MPa、f_t＝2.259MPa，输入计算得到对应的w^-1目标响应值，通过优化拟合计算得到修正解，修正计算结果见表39，通过修正计算可以观察到参数的修正误差都低于5％，参数修正效果较好，修正效果相比于直接用w响应有了较显著提高。

步骤7：模型修正结束。计算表明：采用EI目标响应值基于理论模型及响应面模型的参数修正方法能够对该材料本构模型进行参数修正且修正效果良好，基于理论模型的修正误差可控制在1％以内、基于响应面模型的修正误差可控制在2％以内；采用w、w^-1目标响应值基于响应面模型的参数修正方法同样能够对该材料本构模型进行参数修正，但采用w^-1响应比直接采用w响应的修正效果好，采用w响应的修正误差存在7％左右，采用w^-1响应的修正误差可控制在3％左右。

实施例二：钢筋混凝土梁本构模型参数修正算例

结构拟定为钢筋混凝土简支梁，混凝土强度等级为C45，截面为300mm×450mm的矩形截面，在受拉侧配两根受力筋，采用钢筋类型为HRB335带肋钢筋、公称直径为28mm，配筋面积As＝1231.504mm²，配筋率为9.122％，跨径为5100mm，结构构造设计如图9所示。进行EI响应计算时，测点选择曲率

为0.1×10^-6mm^-1、0.2×10^-6mm^-1、0.4×10^-6mm^-1。

步骤1：选择材料本构方程，在本算例中，混凝土本构方程的采用同算例一，钢筋本构方程采用上升段为线性的本构，如式(55)所示：

式中

步骤2：选择f_c、f_t和E_s三个参数作为待修正参数，以EI_i作为响应值对本构材料参数进行参数修正，混凝土峰值应变ε_c0＝0.002、ε_t0＝0.0001，给f_c、f_t和E_s一个包含理论组的区间为f_c∈[15，50]MPa、f_t∈[1.50，3.10]MPa、E_s∈[170000，220000]MPa。

步骤3：根据平截面假定，设截面受压区高度为z，列轴力平衡方程，并求解z：

步骤4：由式(56)可以判断待修正参数f_c、f_t和E_s线性无关。

步骤5：由力矩平衡，列弯矩平衡方程，并求得截面抗弯刚度：

根据C45混凝土拉、压强度标准值为f_c＝29.6MPa、f_t＝2.51MPa，HRB335普通钢筋的弹性模量规范值为E_s＝2.0×10⁵MPa，模拟理论组(模拟假设结构使用若干年后混凝土强度值和钢筋的弹性模量下降10％的工况)f_c、f_t和E_s取值分别为f_c＝26.64MPa、f_t＝2.259MPa、E_s＝1.8×10⁵MPa。

应用MATLAB对理论试验组进行计算分析，在不同

值时计算得出有效受压区高度z_i、弯矩M_i、EI_i响应的计算结果见表40。

表40试验组在不同

值时计算得出的z_i、M_i、EI_i

对各组本构模型不同测点位置下的工况进行本构参数修正，每个工况分别取三个不同

值

对应的响应值[EI_k,EI_k+1,EI_k+2]作为目标响应值，并对修正效果进行判别，同时计算相应判别式的指标值，经大量分析计算，选取部分修正效果区别显著的修正组，如表41。

表41部分修正组的效果表

注：表中●表示参数修正效果良好，○表示参数修正效果一般，×表示参数修正效果差。

该本构组属于三参数型且受压区高度z≠C，根据判别方法需满足判别式EIR值大于0.7％，以上分析计算皆满足。而本算例中，对选定的

工况进行EIR值计算，由式(46)计算可得：

在该算例的条件下满足判别式EIR值大于0.7％的要求，所以在该工况下能对本构进行参数修正。

步骤6：模型修正

以

为已知量，应用MATLAB或有限元分析软件进行结构计算，截面受压区高度z的计算结果见表42，弯矩M的计算结果见表43，目标响应值EI_i的计算结果见表44。

表42在不同

值时计算得出的受压区高度z(单位:mm)

表43在不同

值时计算得出的弯矩M(单位:kN·m)

表44在不同

值时计算得出的EI_i响应值(单位:N·mm²)

根据无约束最优化问题，通过理论模型计算的目标响应值优化拟合计算得到f_c、f_t和E_s的修正解，比较修正解与理论解的误差。输入理论值对应的目标响应EI_i，求解参数修正解，计算结果见表45，修正结果准确。

表45参数修正的计算结果表

步骤7：模型修正结束。

以上所述仅为本发明的2个实施例，凡依本发明申请专利范围所做的均等变化与修饰，皆属于本发明的涵盖范围。

Claims (8)

1.一种受弯梁材料本构参数模型修正及判别方法，其特征在于包括如下步骤：

(1)根据工程经验，选择合适的材料本构方程；

(2)选取材料本构方程中的待修正参数，并确定合理的待修正参数取值区间，a)仅一个待修正参数时，可进行材料本构参数模型修正，转步骤(6)；b)两个及以上待修正参数时，继续下述步骤；

(3)基于平截面假定，假设截面曲率为

由截面应力平衡，列轴力平衡方程，求解截面受压区高度

表示关于自变量

的函数；

(4)以待修正参数为变量，判断轴力平衡方程中是否存在线性相关项，若存在，则不能进行材料本构参数模型修正，转步骤(7)；若不存在，继续下述步骤；

(5)由力矩平衡，列弯矩平衡方程，求得弯矩M，计算截面抗弯刚度

若对多个不同截面曲率

截面受压区高度

为常数时：若判别式EIB＝0，不能进行材料本构参数模型修正，转步骤(7)；判别式EIB≠0，能进行材料本构参数模型修正，继续下述步骤；

b)截面受压区高度

不为常数时：两个待修正参数，判别式EIA＞0.4％，三个待修正参数，判别式EIR＞0.7％，四个待修正参数，判别式EIR＞1.5％，能进行材料本构参数模型修正，继续下述步骤；否则，不能进行材料本构参数模型修正，转步骤(7)；

(6)采用理论模型或者响应面模型，通过优化算法对待修正参数进行识别，完成材料本构参数的修正，模型修正时目标响应值的个数不小于待修正参数的个数；

(7)模型修正结束；

具体的，步骤(5)中，判别式分别为：

其中，带下标的EI表示对不同的

值，通过理论模型计算得到的截面抗弯刚度，

EI为响应值，n表示待修正参数的数量，下标k表示顺序号，k＝1,2,3...,n-1；

EI_k：表示在选取的多组EI响应值中，按照曲率由小到大排序，曲率为

时对应计算得到的截面抗弯刚度EI_k；

EI_k+1：表示在选取的多组EI响应值中，按照曲率由小到大排序，曲率为

时对应计算得到的截面抗弯刚度EI_k+1，与EI_k相邻；

EI_max：表示在选取的多组EI响应值中，按照曲率由小到大排序，最大的曲率值对应计算得到的响应值；

EI_min：表示在选取的多组EI响应值中，按照曲率由小到大排序，最小的曲率值对应计算得到的响应值。

2.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(1)中材料本构方程可根据相应的设计规范选择，如混凝土、普通钢筋、预应力钢筋本构模型可参照《JTG 3362-2018公路钢筋混凝土及预应力混凝土桥涵设计规范》选用。

3.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(2)中材料本构方程中的待修正参数取值区间要包含待修正参数的实际值。

4.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(3)中，轴力平衡方程与梁截面形状、组成截面材料的种类数量、材料本构方程相关，以单一材料矩形截面梁为例：

∑N＝0；

其中，∑N为截面轴力之和，N表示轴力，B、H分别为矩形截面的宽度和高度，z为截面受压区高度，

分别为材料的受压应力、受拉应力，由材料本构方程计算，

为距离截面中性轴y位置的应变，

为截面曲率，y为距离截面中性轴的位置；

假设截面曲率

已知，可求得应变

通过材料本构方程，可求得受压应力

受拉应力

轴力平衡方程中只有截面受压区高度z未知，故可求解截面受压区高度

表示关于自变量

的函数。

5.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(5)中，以单一材料矩形截面梁为例：

M_i＝∑M；

其中，∑M为截面弯矩之和，M表示弯矩，M_i为与截面曲率

对应的总弯矩。

6.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(6)中，理论模型指采用数学方法建立的精确理论数学模型，响应面模型指由试验设计方法生成样本，拟合得到响应与待修正参数之间的响应面模型，是近似模型，与理论模型相比，有一定的误差，响应面模型可采用多项式模型、克里格模型、径向基函数、BP神经网络等模型。

7.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(6)中，模型修正的目标响应值取截面抗弯刚度EI、梁的挠度w或者梁的挠度倒数w^-1等响应。

8.根据权利要求1所述受弯梁材料本构参数模型修正及判别方法，其特征在于：步骤(6)中，模型修正可采用单目标或多目标优化算法。

Images