简介概要

基于人工神经网络的岩土流变本构模型辨识

来源期刊：中国有色金属学报2002年第5期

论文作者：陈沅江潘长良曹平王文星

文章页码：1027 - 1034

关键词：岩土材料; 流变; 本构模型; 系统辨识; 人工神经网络

Key words：rock and soil material; rheology; constitutive model; system identification; artificial neural network

摘要：岩土流变是岩土工程失稳破坏的重要原因之一。从系统辨识的角度，首先将岩土材料流变本构的一般微分方程通式按实际的采样周期转化为线性时不变SISO系统的离散差分方程格式，构建了用于岩土流变本构模型辨识的BP神经网络模型；然后探讨了该神经网络模型用于岩土流变本构模型辨识的基本步骤以及其网络结构参数（输入层神经元数和网络连接权值）与SISO流变系统差分方程模型参数间相互转化的算法原理，并据此在Matlab软件平台中编制了BP网络辨识算法的相应程序CYJ1.M；最后，采用有关的考题验证证明该辨识算法是成功可信的。

Abstract: The rheology of rock and soil is one of the important reasons why geotechnical engineering is often apt to lose its stability and be damaged. From the view of system identification, the author first changes the general constitutive differential equation of the rheology of rock and soil into the discrete difference equation of the linear time—invariant single—input single—output system using the actual sampling period and builds the BP neural network model to identify the rheological constitutive model of rock and soil, and then probes into the basic steps of the identification of the rheological constitutive model of rock and soil using the BP model and studies the algorithmic principle to translate the BP network structure parameters (the number of the network’s neurons and the values of the network’s connections weights) into the model parameters of the difference equation of the SISO rheological system. Based on the principle the author also compiles the application program CYJ1.M for the identification arithmetic using the BP neural network on the base of MATLAB software .At last, by using an examination question, the identification arithmetic is testified to be successful and reliable.

DOI:10.19476/j.ysxb.1004.0609.2002.05.034

基于人工神经网络的岩土流变本构模型辨识

陈沅江潘长良曹平王文星

中南大学资源环境与建筑工程学院采矿与岩土工程研究所

中南大学资源环境与建筑工程学院采矿与岩土工程研究所长沙410083

摘要：

岩土流变是岩土工程失稳破坏的重要原因之一。从系统辨识的角度 , 首先将岩土材料流变本构的一般微分方程通式按实际的采样周期转化为线性时不变SISO系统的离散差分方程格式 , 构建了用于岩土流变本构模型辨识的BP神经网络模型 ;然后探讨了该神经网络模型用于岩土流变本构模型辨识的基本步骤以及其网络结构参数 (输入层神经元数和网络连接权值 ) 与SISO流变系统差分方程模型参数间相互转化的算法原理 , 并据此在Matlab软件平台中编制了BP网络辨识算法的相应程序CYJ1.M ;最后 , 采用有关的考题验证证明该辨识算法是成功可信的

关键词：

岩土材料;流变;本构模型;系统辨识;人工神经网络;

中图分类号： TP183

收稿日期：2001-03-12

基金：国家九五科技攻关资助项目 ( 96 -116 -0 1-0 1-0 2 );

Identification of rheological constitution model of rock and soil based on artificial neural network

Abstract：

The rheology of rock and soil is one of the important reasons why geotechnical engineering is often apt to lose its stability and be damaged. From the view of system identification, the author first changes the general constitutive differential equation of the rheology of rock and soil into the discrete difference equation of the linear time-invariant single-input single-output system using the actual sampling period and builds the BP neural network model to identify the rheological constitutive model of rock and soil, and then probes into the basic steps of the identification of the rheological constitutive model of rock and soil using the BP model and studies the algorithmic principle to translate the BP network structure parameters (the number of the network's neurons and the values of the network's connections weights) into the model parameters of the difference equation of the SISO rheological system. Based on the principle the author also compiles the application program CYJ1.M for the identification arithmetic using the BP neural network on the base of MATLAB software .At last, by using an examination question, the identification arithmetic is testified to be successful and reliable.

Keyword：

rock and soil material; rheology; constitutive model; system identification; artificial neural network;

Received： 2001-03-12

当前, 随着数值计算方法的普及, 岩土工程施工设计安全可靠的重要保障归结到了岩土本构模型和输入参数的正确性和科学性上 ^[1,2] , 这在考虑岩土流变问题时尤以为甚。由此引出了岩土流变本构模型及其参数的辨识课题。文献 [ 3] 在岩石蠕变试验的基础上, 通过分离材料蠕变变形的各种分量, 采用“试凑法”对岩石流变本构模型及其参数进行了辨识。但“试凑法”在岩石蠕变曲线的处理过程中带有人为的主观性。文献 [ 4] 将n阶线性偏微分方程视为岩石流变本构的一般数学描述, 并将其引入围岩位移的观测方程中, 利用最小二乘法对本构模型结构及其参数进行了辨识。然而岩体系统的输出测试数据在实际中常受噪声污染, 这使最小二乘算法误差较大 ^[5] ; 文献 [ 6] 利用遗传算法并结合有限元计算对岩石流变本构模型及其参数进行了辨识, 由于遗传算法目前来讲其计算效率不是很高, 有许多需改进的地方, 且计算量大 ^[7] , 故此法目前来讲仍不成熟。为了探索岩土流变本构模型辨识的新的有效途径, 本文作者利用Matlab软件平台, 对现已成熟的人工神经网络方法用于岩土流变本构模型辨识的可行性及有关算法进行了研究, 通过有关考题验证表明, 该方法成功可信, 宜于实际推广。

1 人工神经网络模型及其辨识原理

人工神经网络是人脑神经网络的结构模型, 是一种非线性动力学系统, 它具有自学习、自组织、高容错性、并行处理、分布式知识存储、非线性动态处理等优点 ^[8] 。目前已有许多模型, 其中以BP网络模型理论最为成熟, 本文作者便采用这种网络来进行研究。 BP神经网络是以BP算法为基础建立的, 模型中主要有3类节点:输入节点、隐节点和输出节点。其算法的基本思想是:外界信息首先经输入节点, 传输到隐节点, 最后传送到输出节点, 完成网络的正向传播; 若输出不是期望输出, 则将实际输出与期望输出的误差, 经原连通路径返回, 修正各层间的连接权值, 即误差反传过程, 使误差变小。经反复地应用此算法, 使误差最小, 即使系统输出与实际输出最接近。从系统论的观点看, 人工神经网络的训练学习过程实际上是系统辨识的过程。系统辨识, 就是对一个未知的系统, 在通过大量试验获得其输入和输出数据的基础上, 从一组给定的模型中, 确定一个与所测系统等价的模型。所谓确定模型就是确定模型的各种参数。对于一个SISO系统来说, 其模型结构可用差分方程表达如下 ^[9] :

$y (k) = \sum_{i = 1}^{n} a_{i} y_{k - i} + \sum_{j = 0}^{m} b_{j} u_{k - j} ? ? ? (1)$

式中 a_i, b_j (i=1, 2, …, n; j=1, 2, …, m) 为系统模型参数; m, n为模型结构参数, 它们都有待辨识。若将该方程右端的输入输出延迟量作为神经网络的输入, 左端的待求输出作为神经网络的输出, 则差分模型的各参数便隐含在神经网络模型中, 修正网络权值使系统的误差函数达到极小值的过程便是系统辨识的过程。现从数学的角度来考察神经网络的这一辨识过程。据式 (1) , 可构造一具有n+m+1个输入节点、一个输出节点的神经网络模型如图1所示, 其中隐节点个数可据经验或自适用的方法确定 ^[10] , 这里采用经验方法取隐层节点个数为输入节点个数的两倍再加上1个节点以加快误差的下降速度。图中y_k-n, y_k-n+1, …, y_k-1, u_k, u_k-1, …, u_k-m为输入层神经元的输入; Z_h[h=1, 2, …, 2 (n+m+2) +1。]为隐含层神经元的输出; y_k为输出层神经元的输出。设输入层和隐含层间的连接权值为V_{i h} (i=1, 2, …, n+m+1; h=1, 2, …, 2 (n+m+1) +1。 ) ; 隐含层与输出层间的连接权值W_h0 (h=1, 2, …, 2 (n+m+1) +1) 。现假定有L个样本对, 即 $[y_{k - n}^{(p)}, y_{k - n + 1}^{(p)}, ?, y_{k - 1}^{(p)}, u_{k}^{(p)}, ? ? u_{k - m}^{(p)}, \hat{y}_{k}^{(p)}] (p = 1 ? 2 ? ? ? L)$ 用于神经网络的训练学习以确定其权值V_{i h}和W_h0。上述样本对中的 $\hat{y}_{k}^{(p)}$ 为网络的实际输出, 若设y $_{k}^{(p)}$ 为网络的计算输出, 则决定权值的标准是考虑实际的输出 $\hat{y}_{k}^{(p)}$ 与计算的输出y $_{k}^{(p)}$ 的差值的平方为最小, 即输出误差 $E^{(p)} = (y_{k}^{(p)} - \hat{y}_{k}^{(p)})^{2}$ 为最小。对全部样本值来说, 其总误差为 $E_{Ζ} = \sum_{Ρ = 1}^{L} E^{(p)} 。$ 将E_Z换成 $E = \frac{1}{2} E_{Ζ}$ , 则

将E作为全局代价函数 (即目标函数) , 对全部样本来说, E为最小时所对应的权值V_{i h}, W_h0即为辨识过程所求的参数 (即人工神经网络在学习或经过训练后的权值) 。按人工神经网络的计算原理, 输出层节点在第p次样本对时的计算输出为

$y_{k}^{(p)} = f_{1} (\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} ? Ζ_{h}^{(p)}) ? ? ? (3)$

同理, 隐含层节点在第p次样本对时的计算输出为

$\begin{array}{l} Ζ_{h}^{(p)} = f_{2} (\sum_{i = 1}^{n} V_{i h} ? y_{k - i}^{(p)} + \sum_{i = 0}^{m} V_{i h} ? u_{k - 1}^{(p)}) \\ (h = 1 ? 2 ? ? ? 2 (n + m + 1) + 1 。) ? ? ? (4) \end{array}$

于是, 神经网络的辨识过程可分为如下几个步骤:

1) W_h0的辨识

将系统隐含层和输出层之间的权值W_h0视为参数, 据LMS法, 当W_h0沿E的负梯度方向变化时, E将下降最快。取适当步长, 设迭代前的值为W_h0, 迭代后的值为W ′_h0, 于是有

$W^{'}_{h 0} = W_{h 0} - λ ? \frac{? E}{? W_{h 0}} ? ? ? (5)$

式中 λ>0, 一般取0<λ<1。现考察?E/?W_h0, 有下式:

$\begin{array}{l} \frac{? E}{? W_{h 0}} = ? [\frac{1}{2} \sum_{p = 1}^{L} E_{k}^{(p)}] / ? W_{h 0} \\ = \sum_{p = 1}^{L} (y_{k}^{(p)} - \hat{y}_{k}^{(p)}) ? \frac{? y_{k}^{(p)}}{? W_{h 0}} \\ = \sum_{p = 1}^{L} [f_{1} (\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} ? Ζ_{h}^{(p)}) - \hat{y}_{k}^{(p)}] ? \\ \frac{? [f_{1} (\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} ? Ζ_{h}^{(p)})]}{? W_{h 0}} \\ = \sum_{p = 1}^{L} [f_{1} (\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} ? Ζ_{h}^{(p)}) - \hat{y}_{k}^{(p)}] ? \\ f^{'}_{1} (\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} ? Ζ_{h}^{(p)}) ? Ζ_{h}^{(p)} \end{array}$

图1 神经网络辨识模型结构图

Fig.1 Structure schematic of neural network model

结合式 (4) , 有

?E/?W_h0=F₁ (y $_{k - i}^{(p)}$ , u $_{k - i}^{(p)}$ , V_ih, W_h0) (6)

将式 (6) 代入式 (5) , 有迭代式

W′_h0=W_h0-λ·F₁ (y $_{k - i}^{(p)}$ , u $_{k - i}^{(p)}$ , V_ih, W_h0) (7)

2) V_{i h}的辨识

将系统输入层和隐含层之间的权V_{i h}视为参数, 同样有

$V^{'}_{i h} = V_{i h} - β ? \frac{? E}{? V_{i h}} ? ? ? (8)$

式中 β>0, 一般取0<β<1。考察?E/?V_ih, 有

$\begin{array}{l} \frac{? E}{? V_{i h}} = ? [\frac{1}{2} \sum_{p = 1}^{L} E_{k}^{(p)}] / ? V_{i h} \\ = ? [\frac{1}{2} \sum_{p = 1}^{L} (y_{k}^{(p)} - \hat{y} (_{k}^{(p)})^{2}] / ? V_{i h} \\ = \sum_{p = 1}^{L} (y_{k}^{(p)} - \hat{y}_{k}^{(p)}) \frac{? y_{k}^{(p)}}{? V_{i h}} \\ = \sum_{p = 1}^{L} (y_{k}^{(p)} - \hat{y}_{k}^{(p)}) ? \\ f^{'}_{1} (\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} z_{h}^{(p)}) ? \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} ? \frac{? Ζ_{h}^{(p)}}{? V_{i h}} ? ? ? (9) \end{array}$

再考察式 (9) 中的?Z $_{h}^{(p)}$ /?V_ih, 有

$\begin{array}{l} \frac{? Ζ_{h}^{(p)}}{? V_{i h}} = f^{'}_{2} (\sum_{i = 1}^{n} V_{i h} y_{k - i}^{(p)} + \sum_{i = 0}^{m} V_{i h} u_{k - i}^{(p)}) ? \\ {\begin{cases} y_{k - i}^{(p)} (当 i ＜ n) \\ u_{k - i}^{(p)} (当 n ＜ i ＜ n + m + 1) \end{cases} ? ? ? (10) \end{array}$

综合 (9) , (10) 两式, 有

?E/ ?V_ih=F₂ (y $_{k - i}^{(p)}$ , u $_{k - i}^{(p)}$ , V_ih, W_h0) (11)

将式 (11) 代入式 (8) 得

V′_ih=V_ih-β·F₂ (y $_{k - i}^{(p)}$ , u $_{k - i}^{(p)}$ , V_ih, W_h0) (12)

3) min E

据式 (7) , (12) 计算得W′_h0, V′_{i h}后, 再利用样本 $[y_{k - n}^{(p)}, y_{k - n + 1}^{(p)}, ?, y_{k - 1}^{(p)}, u_{k}^{(p)}, u_{k - 1}^{(p)}, ?, u_{k - m}^{(p)}, \hat{y}_{k}^{(p)}] (p = 1 ? 2 ? ? ? L)$ 求出新的输出误差E′, 对?ε>0, 使E′<ε或|E′-E|<ε作为循环结束条件, 当条件满足, 停止循环, 否则, 继续循环。

4) 确定系统参数值

当循环结束时, 神经网络模型便满足了系统的输入输出条件, 且对应于差分方程 (1) , 网络模型中输入层的个数便是系统差分方程的结构参数, 而所确定的权值W′_h0, V′_{i h}便隐含着其模型参数a_i, b_j。

以上便是神经网络用于系统辨识时自身训练学习的4个步骤。当这一训练过程结束后, 网络中各神经元之间的连接权值便已确定, BP网络的输出与用差分方程表达的系统的输出也相吻合。但此时BP网络的连接权值只是一些数值, 尽管它隐含了系统差分方程的参数, 却并不能从其中得出被辨识系统的任何特性, 因此, BP网络的辨识结果很不直观, 需将其辨识所得的权值转换成系统差分方程中的参数。设神经网络模型中隐含层各单元的处理函数为S型函数, 即

$f_{1} (x) = \frac{1}{1 + e^{- x}} ? ? ? (13)$

而输出层单元处理函数为线性函数:

f₂ (x) = x (14)

设神经网络经训练学习后所确定的各层权值为V_ih, W_h0 (其中:i=1, 2, …, n, n+1, …, n+m+1;

h=, 1, 2, …, 2 (n+m+1) +1) , 于是有以下表达式:

$\begin{array}{l} y_{k} = \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} f_{1} (\sum_{i = 1}^{n} V_{i h} y_{k - i} + \sum_{j = 0}^{m} V_{j h} u_{k - j} + \\ b_{h}) + b_{0} ? ? ? (15) \end{array}$

式中 b_h为隐含层各单元的偏差, b₀为输出层单元的偏差。将f₁函数用麦克劳林级数展开并甩去高阶项代入上式得:

$\begin{array}{l} y_{k} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} (V_{1 h} y_{k - 1} + V_{2 h} y_{k - 2} + ? + V_{n h} y_{k - n} + V_{(n + 1) h} u_{k} + V_{(n + 2) h} u_{k - 1} + ? + V_{(n + m + 1) h} u_{k - m}) + \\ {\sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 0} [\frac{1}{2} + \frac{1}{4} b_{h} + R_{3} (ξ)] + b_{0}} ? ? ? (16) \end{array}$

式中 R₃ (ξ) 为麦克劳林级数的误差项, 比较式 (16) 和差分方程式 (1) , 可近似认为差分方程的待辨识模型参数为

$a_{1} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 1} V_{1 h}$ ; $a_{2} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 1} V_{2 h};$

…; $a_{n} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 1} V_{n h};$

$b_{0} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 1} V_{(n + 1) h};$

$b_{1} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 1} V_{(n + 2) h}$ ;

…; $b_{m} = \frac{1}{2} \sum_{h = 1}^{2 (n + m + 1) + 1} W_{h 1} V_{(n + m + 1) h} ? ? ? (17)$

此时误差为式 (16) 式中最后一项的值, 它就是将非线性系统近似看作线性系统时带来的误差 ^[11,12] 。

2 岩土流变本构方程的辨识

岩土材料在外界荷载下常表现为粘弹、粘塑、弹粘塑等复杂特性, 可用弹簧、滑块和粘壶等3类元件通过串、并联组合成的复合元件模型来对其进行反映 ^[13] , 其应力—应变—时间关系即本构模型的通式可用线性常微分方程表达如下:

$\begin{array}{l} c_{1} \frac{d^{n} ε (t)}{d t^{n}} + c_{2} \frac{d^{n - 1} ε (t)}{d t^{(n - 1}} + ? + c_{n} \frac{d ε (t)}{d t} + \\ c_{n + 1} ε (t) = d_{1} \frac{d^{m} σ (t)}{d t^{m}} + d_{2} \frac{d^{m - 1} σ (t)}{d t^{m - 1}} + ? + \\ d_{m} \frac{d σ (t)}{d t} + d_{m + 1} σ (t) ? ? ? (18) \end{array}$

其表达的是一以应力σ (t) 为输入、应变ε (t) 为输出的单输入单输出线性时不变连续时间模型 ^[14] 。对其两边进行Laplace变换, 并假定初始条件为零, 则可得其传递函数表达为:

$G (s) = \frac{\overset{—}{σ} (s)}{\overset{—}{ε} (s)} = \frac{d_{1} s^{m} + d_{2} s^{m - 1} + ? + d_{m + 1}}{c_{1} s^{n} + c_{2} s^{n - 1} + ? + c_{n + 1}} ? ? ? (19)$

然而在岩土力学试验或工程监测中, 得到的往往是关于岩土流变过程之应力、应变时间序列的离散数据信息, 依此信息由上节所述的神经网络辨识算法可辨识得到系统的离散化差分方程模型式 (1) 。假定在试验或监测中测定不同载荷σ (t_i) 作用下的应变ε (t_i) 时间序列时, 所用的采样周期为T, 则采用差分方程式 (1) 表达的岩石流变系统的相应力学模型可化为 ^[9]

$\begin{array}{l} ε [(i + n) Τ] + (- a_{1}) ε [(i + n - 1) Τ] + \\ (- a_{2}) ε [(i + n - 2) Τ] + ? + \\ (- a_{n - 1}) ε [(i + 1) Τ] + (- a_{n}) ε (i Τ) = \\ b_{0} σ [(i + m + 1) Τ] + b_{1} σ [(i + m) Τ] + ? + \\ b_{m} σ (i Τ) ? ? ? (20) \end{array}$

对上式进行z变换, 则可得到其离散传递函数为:

$\begin{array}{l} G (z) = \\ \frac{b_{0} z^{m} + b_{1} z^{m - 1} + ? + b_{m - 1} z + b_{m}}{z^{n} + (- a_{1}) z^{n - 1} + ? + (- a_{n - 1}) z + (- a_{n})} ? ? ? (21) \end{array}$

按双线形变换法 ^[15] , 传递函数 (19) 和 (21) 中的变量s与z存在如下关系:

$s = \frac{2}{Τ} ? \frac{1 - z^{- 1}}{1 + z^{- 1}} ? ? ? (22)$

由此可实现连续模型 (18) 与离散模型 (20) 之间的相互转化。值得指出的是, 在该转化过程中, 系统的采样周期应取得足够小, 以使离散化模型与连续性模型的等效性足够好。

至此, 可总结岩土材料流变本构模型的人工神经网络辨识步骤如下:

步1. 据试验或监测中所获得的应力应变时间序列数据[σ (t_i) , ε (t_i) , t_i], 按采样周期T对其进行适当的离散处理, 得到用于神经网络训练学习的输入和系统实际输出数据。

步2. 据岩土流变的有关先验知识, 确定流变模型通式中的模型阶次n, m, 由此按采样周期T写出流变岩土系统的差分方程式 (20) 或 (1) 。

步3. 据差分方程式 (20) 或 (1) 构建人工神经网络模型, 按2节所述利用步1中数据对其进行训练。

步4. 利用网络训练的结果, 根据式 (17) 计算系统差分方程 (20) 或 (1) 中的待辨识参数, 至此系统的差分方程已完全确定。

步5. 按双线性变化法利用式 (19) , (21) , (22) 将确定的差分方程转换为系统连续模型 (18) 的形式, 至此, 模型的各力学参数已经确定。

步6. 据辨识的最终结果进行仿真计算, 并按有关误差准则考察所得模型是否满足精度要求, 若是, 则其为最佳的模型; 否则, 改变模型阶次数n, m, 返回步2重新进行辨识。这里, 误差准则函数选用 $J (y, y_{m}) = [\sum_{i = 1}^{n} (y_{i} - y_{m i})^{2})]^{\frac{1}{2}}$ 。上述各步计算可在Matlab软件包中轻松实现, 笔者 ^[16,17] 由此编制了用于岩土流变参数辨识的程序CYJ1.M, 图2所示为其有关的算法流程。

3 辨识算法的考题验证

所谓考题验证是利用已知的有明确答案的题目, 来考核所设计的辨识算法成功与否。首先据考题中已知的模型参数正分析得到“测试”数据, 然后据此数据利用所设计的算法对所选的模型参数进行辨识, 如果其能再回到原来已知的模型参数, 则辨识是成功的, 方法是可行的, 反之则反。所用考题如下:一标准长方体砂岩试件, 在恒定单轴压应力作用下发生蠕变变形, 据试验分析已知岩石的流变模型为Burgers模型, 其流变本构方程为 ^[18,19]

$E_{1} \overset{? ?}{ε} + \frac{E_{1} E_{2}}{η_{2}} \dot{ε} = \overset{? ?}{σ} + (\frac{E_{1}}{η_{1}} + \frac{E_{2}}{η_{2}} + \frac{E_{1}}{η_{2}}) \dot{σ} + \frac{E_{1} E_{2}}{η_{1} η_{2}} σ ? ? ? (23)$

且E₁=200 MPa, η₁=0.02 MPa·d, E₂=400 MPa, η₂=0.01 MPa·d。其连续性方程和在离散周期为T=1.0×10^-5d时离散方程的具体形式分别如下:

$\overset{? ?}{ε} + 40 000 \dot{ε} = 0.005 \overset{? ?}{σ} + 700 \dot{σ} + 2 000 000 σ ? ? ? (24)$

图2 岩土流变参数辨识程序设计流程图

Fig.2 Flow chart of program of parameter identification of rock and soil rheological model

ε (k) =1.67ε (k-1) +0.670 3ε (k-2) +0.005σk-0.004 143σ (k-1) -0.000 692 5σ (k-2) (25)

模型在两级蠕变荷载σ₁=1.0 MPa, σ₂=1.5 MPa下的连续和离散仿真输出与时间关系分别如图3 (a) 中实线和 (b) 所示。据模型的输入输出时间序列按上节所述人工神经网络辨识算法对该砂岩的流变模型进行辨识, 其主要结果如表1所示, 可见, 表中模型3为最优模型, 其在上述两级应力下蠕变连续仿真输出如图3 (a) 中的虚线所示。模型存在一定的误差, 这是由于将神经网络计算结果近似转化为模型参数时舍去高阶项引起的, 但误差在工程允许的范围内, 故认为计算方法是成功可信的。实际岩土工程中存在大量加载条件可控、位移输出可测的实例, 如不同荷载下岩基的位移随时间的变化、不同围压下隧道周边的位移变化等等 ^[20,21] , 它们亦可近似描述为线性时不变系统, 从而利用上述方法进行辨识, 这方面的内容笔者将另文讨论。

5 结论

岩土材料流变是岩土工程失稳破坏的重要原因之一, 对其流变本构规律进行研究具有十分重大的现实意义。从系统辨识的角度, 探讨了人工神经网络用于岩土材料流变本构模型辨识的可行性及其算法原理, 在此基础上, 基于matlab平台编制了相应辨识程序。有关的考题验证表明, 人工神经网络用于岩土材料流变本构模型的辨识是成功可信的, 不失为岩土流变模型辨识的又一条有效途径。由于当今人工神经网络理论较为成熟, 充分利用其有关成果于岩土流变本构模型辨识领域的研究必将具有广阔的发展前景。

图3 某砂岩的蠕变仿真输出

Fig.3 Simulating output of certain sandstone's creep (a) —Continuous model; (b) —discrete model

表1 某砂岩的模型辨识结果

Table 1 Results of identification of rheological model of certain sandstoneStructure of model

Structure of model			Discrete model identified	Continuous model identified	Value of error function
No.	Output order n	Input order m	Discrete model identified	Continuous model identified	Value of error function
1	1	1	ε (k) =0.820 6ε (k-1) +0.004 435σ (k) + 0.000 048 56σ (k-1)	$\dot{ε} (t) + 19 770 ε (t) = 0.004 435 \dot{σ} (t) +$ 494.1σ (t)	0.028
2	2	1	ε (k) =1.035ε (k-1) +0.149 5ε (k-2) = 0.003 528 6σ (k) -0.000 434 56σ (k-1)	$\begin{array}{l} \overset{? ?}{ε} (t) + 1.901 \times 10^{5} \dot{ε} (t) + 2.608 \times 10^{9} ε (t) = \\ 352.8 \dot{σ} (t) + 7.057 \times 10^{7} σ (t) \end{array}$	0.013 4
3	2	2	ε (k) =1.67ε (k-1) +0.67ε (k-2) + 0.004 6σ (k) -0.004 144 6σ (k-1) - 0.000 693σ (k-2)	$\begin{array}{l} \overset{? ?}{ε} (t) + 4.01 \times 10^{4} \dot{ε} (t) - 1.668 \times 10^{- 5} ε (t) = \\ 0.004 7 \overset{? ?}{σ} (t) + 699.7 \dot{σ} (t) + 2 \times 10^{6} σ (t) \end{array}$	3.907 2×10^-4
4	2	3	ε (k) =0.726 7ε (k-1) +0.101 1ε (k-2) + 0.004 845×σ (k) +4.265×10^-5σ (k-1) - 0.000 262 4σ (k-2)	$\begin{array}{l} \overset{?}{ε} + 2.291 \times 10^{5} \overset{? ?}{ε} (t) + 4.601 \times 10^{10} \dot{ε} (t) + \\ 6.332 \times 10^{14} ε (t) = 317.3 \overset{? ?}{σ} (t) + \\ 1.529 \times 10^{8} \dot{σ} (t) + 1.701 \times 10^{13} σ (t) \end{array}$	0.010 5