应用Laplace变换实现沥青混合料材料参数快速转换的方法与流程

本发明属于沥青混合料材料参数换算领域，具体涉及一种基于laplace变换实现沥青混合料材料参数快速转换的方法。

背景技术：

沥青混合料的松弛模量是开展沥青路面结构力学计算的材料参数之一。然而，目前并没有直接测试沥青混合料松弛模量仪器，而是采用沥青混合料静态蠕变试验获取沥青混合料的蠕变柔量，然后利用黏弹性力学中松弛模量与蠕变柔量的卷积公式进行松弛模量的计算。该方法首先是将卷积公式离散成有限个无穷积分式，然后将区段中值时刻的模量值假定为该区段的模量值，即将模量假定为一个定值，这样的简化处理使得卷积公式变成一个简单的积分式，通过数值积分的方法获取每个区段的积分值，然后通过等式反算松弛模量值。但该方法存在如下不足：

(1)每一个区段的模量值假定为一个定值，这样的处理方法需要将区间划分的很小，才能保证计算的精度，而多小的区段才能保证计算精度并没有明确处理方法。

(2)该种方法只能计算有限个时间点下的松弛模量值，不能计算任意时刻的松弛模量值。

(3)该方法需要处理复杂的数值积分，需要编写复杂的计算程序，计算效率低。

技术实现要素：

本发明的目的是为了解决现有计算沥青混合料的松弛模量的方法精度难以保证，算法复杂的问题，而提供一种应用laplace变换实现沥青混合料材料参数快速转换的方法。

本发明应用laplace变换实现沥青混合料材料参数快速转换的方法按下列步骤实现：

一、沥青混合料蠕变柔量表征

测试固定荷载下沥青混合料的蠕变变形量△l，然后采用方程(1)计算沥青混合料的蠕变柔量；

方程(1)中△l是加载前与加载结束后高度差，h是试件的高度，σ0是荷载大小；

采用prony级数模型对沥青混合料的蠕变柔量进行拟合，prony级数的表达式如方程(2)所示：

方程(2)中j0和ji是prony级数的参数(其中j0是t＝0时的prony参数，ji通过matlab数据拟合获得)，τi为延迟时间，t是时间，然后曲线拟合得到沥青混合料prony级数模型；

二、laplace域松弛模量表征

对沥青混合料prony级数模型施加laplace正变换，prony级数模型的laplace正变换如方程(3)所示：

方程(3)中s是laplace变量；

沥青混合料为(典型的)黏弹性材料，黏弹性材料中松弛模量与蠕变柔量的laplace域换算关系如方程(4)所示：

方程(4)中e(s)是laplace域的松弛模量；

根据方程(4)得到laplace域的松弛模量计算公式如下：

方程(5)中l^-1()是laplace逆变换，e(t)是时域中的松弛模量；

将方程(3)代入到方程(5)中，解得：

即得到laplace域松弛模量；

三、松弛模量时域值求解

采用ft法(fixedtalbot算法)对方程(6)施加laplace逆变换，完成沥青混合料松弛模量时域解的求解，从而完成沥青混合料材料的参数快速转换。

本发明先在室内开展沥青混合料蠕变试验，并采用prony级数来表征沥青混合料的蠕变柔量获取模型参数，然后对prony级数施加laplace正变换；之后基于黏弹性力学当中松弛模量与蠕变柔量的laplace转换公式，将蠕变柔量的laplace表达式代入转换公式当中，完成laplace域松弛模量表达式的求解；最后对表达式施加laplace数值逆变换，完成松弛模量时域值的求解。

本发明应用laplace变换实现沥青混合料材料参数快速转换的方法包含以下有益效果：

1、laplace域松弛模量表达式的推导简单易行，概念清晰；

2、避免了复杂的数值积分运算及程序编写；

3、可快速、精确的获取松弛模量值；

4、通过laplace积分逆变换，可完成任意时刻下松弛模量的求解。

附图说明

图1是实施例一中20℃下四种沥青混合料蠕变柔量，其中■代表type-a，●代表type-b，▲代表type-c，▼代表type-d；

图2是实施例一中广义kelvin模型的结构示意图；

图3是实施例一中采用卷积方法和ft法计算得到的松弛模量值(r²＝0.999)，其中▼代表卷积方法，▲代表ft法；

图4是实施例一中m＝5和m＝8的沥青混合料松弛模量值(r²＝1.00)，其中▼代表m＝5，▲代表m＝8；

图5是实施例一中m＝5和m＝10的沥青混合料松弛模量值(r²＝1.00)，其中▼代表m＝5，▲代表m＝10；

图6是实施例一中m＝5和m＝15的沥青混合料松弛模量值(r²＝1.00)，其中▼代表m＝5，▲代表m＝15；

图7是本发明应用laplace变换实现沥青混合料材料参数快速转换的方法的流程图。

具体实施方式

具体实施方式一：本实施方式应用laplace变换实现沥青混合料材料参数快速转换的方法按下列步骤实施：

一、沥青混合料蠕变柔量表征

测试固定荷载下沥青混合料的蠕变变形量△l，然后采用方程(1)计算沥青混合料的蠕变柔量；

方程(1)中△l是加载前与加载结束后高度差，h是试件的高度，σ0是荷载大小；

采用prony级数模型对沥青混合料的蠕变柔量进行拟合，prony级数的表达式如方程(2)所示：

方程(2)中j0和ji是prony级数的参数(其中j0是t＝0时的prony参数，ji通过matlab数据拟合获得)，τi为延迟时间，t是时间，然后曲线拟合得到沥青混合料prony级数模型；

二、laplace域松弛模量表征

对沥青混合料prony级数模型施加laplace正变换，prony级数模型的laplace正变换如方程(3)所示：

方程(3)中s是laplace变量；

沥青混合料为(典型的)黏弹性材料，黏弹性材料中松弛模量与蠕变柔量的laplace域换算关系如方程(4)所示：

方程(4)中e(s)是laplace域的松弛模量；

根据方程(4)得到laplace域的松弛模量计算公式如下：

方程(5)中l^-1()是laplace逆变换，e(t)是时域中的松弛模量；

将方程(3)代入到方程(5)中，解得：

即得到laplace域松弛模量；

三、松弛模量时域值求解

采用ft法(fixedtalbot算法)对方程(6)施加laplace逆变换，完成沥青混合料松弛模量时域解的求解，从而完成沥青混合料材料的参数快速转换。

本实施方式步骤三中尝试采用durbin法和dubner法进行求解，其公式分别如下：

式中a，b，t——计算参数；

t——具体计算时间；

re，im——复数的实部和虚部。

在完成松弛模量时域值求解之前，需要选取计算效率高、精度高的laplace数值逆变换方法，采用如下原函数形式，函数的参数如下表1。

表1prony级数系数

表2不同数值反演方法计算效率

本实施方式发现随着时间的增加，durbin法和dubner法反演的精度不断下降，而ft法能够较好满足反演的计算精度。从表2三者计算效率可以看出，ft法不仅计算精度要高于durbin法和dubner法，而且其计算效率非常高，逆变换50个数据点只需要0.04ms。ft在prony级数laplace反演中不仅保证了计算精度，而且能够很好的保证计算效率。ft法能够较好的应用到沥青混合料黏弹参数的换算当中。

ft法中的参数m是指将laplace反演的路径划分的步数，这类似于数值积分方法中的积分步数。步长越小，步数越多，计算精度越高，但计算精度越低。kim在研究层状黏弹体系理论中，提出ft法中任意的参数m值不会影响反演精度。但其研究的是短时间内荷载函数的laplace逆变换，其研究结论不能应用于prony级数的laplace逆变换中。因此，研究参数m等于5、8、10、15时ft法对prony级数反演的精度和效率。

表3不同m值下ft反演效率

从表3可以看出，无论m值选为多少，对ft法的计算精度和效率的影响都很小。所以，采用ft法研究沥青混合料黏弹参数换算。

因此，对方程(6)施加laplace逆变换(ft法)，即可完成松弛模量时域解的求解。

具体实施方式二：本实施方式与具体实施方式一不同的是步骤一中沥青混合料制成直径100mm，高度150mm的标准试件。

具体实施方式三：本实施方式与具体实施方式一不同的是步骤一中蠕变变形量的测试过程是在静态加载下，控制荷载大小为0.3mpa，加载时间为240s。

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是步骤一中prony级数的表达式中当i＝1时，则τi＝0.001；当i＝2时，则τi＝0.01；当i＝3时，则τi＝0.1；当i＝4时，则τi＝1；当i＝5时，则τi＝10；当i＝6时，则τi＝100；当i＝7时，则τi＝1000。

具体实施方式五：本实施方式与具体实施方式一至四之一不同的是步骤一中利用maltab软件中曲线拟合得到沥青混合料prony级数模型。

具体实施方式六：本实施方式与具体实施方式一至五之一不同的是步骤三中所述的ft法的算法公式如方程(9)所示，方程(6)中laplace域松弛模量表达式e(s)为方程(9)中f()函数的表达式，即将f(r)和f(s(θk))中的自变量值r和s(θk)的值代入e(s)中的s变量，分别求取e(s)的值即为f(r)和f(s(θk))的值：

其中，θk＝kπ/m，s(θ)＝rθ(cotθ+i)，σ(θ)＝θ+(θcotθ-1)cotθ，r＝2m/(5t)；

式中t——计算时间，r，m——模型参数。

具体实施方式七：本实施方式与具体实施方式六不同的是m的取值为5～15。

具体实施方式八：本实施方式与具体实施方式七不同的是m的取值为5～8。

实施例一：本实施例应用laplace变换实现沥青混合料材料参数快速转换的方法按下列步骤实施：

一、沥青混合料蠕变柔量表征

在实验室内成型4种直径为150mm，高度170mm的圆柱体沥青混合料试件，分别命名为type-a、type-b、type-c、type-d，然后钻芯切割成直径100mm，高度150mm的标准试件，采用多功能材料试验机(utm-250)在静态加载下，荷载大小为0.3mpa，加载时间为240s，通过公式(1)计算沥青混合料的蠕变柔量，如图1所示；

方程(1)中△l是加载前与加载结束后高度差，h是试件的高度，σ0是荷载大小；

将蠕变柔的试验数据导入到matlab中，通过自带工具箱curvefitting确定不同沥青混合料的prony参数如下；

表4四种沥青混合料prony级数参数值

二、laplace域松弛模量表征

线黏弹性材料的应力-应变-时间关系，主要有微分和积分型两大类。恒定荷载下黏弹材料的力学行为通过时间相关的材料函数来描述，而对于交变荷载作用下的黏弹材料性能，采用频率相关的复函数来描述。该部分内容主要基于广义kelvin模型微分型本构关系，介绍黏弹材料蠕变柔量、松弛模量及动态模量的计算过程；

经典的黏弹本构由弹簧和黏糊串/并联组合而成，由此组成不同的本构模型，如：maxwell、kelvin模型及burgers模型等，广义kelvin模型由n个kelvin链组成，设kelvin链的第i个单元的应变为εi，其弹簧弹性模量和阻尼器黏性系数分别为ei和ηi，用微分算子表示应变：

由n个kelvin单元组成的广义kelvin模型的总应变：

将此展开经整理后，得到一般模型的本构方程：

写作

或

pσ＝qε(22)

其中微分算子：

为求得式(22)的解，对公式施加laplace积分变换，则得到代数方程：

或

为了求得本构方程的蠕变柔量函数，将固定荷载σ＝σ0h(t)的laplace变换式代入(24)中，求得蠕变柔量的像函数：

通过laplace逆变换求得蠕变柔量的时域解：

同样，为了求得松弛模量的时域解，将应变函数的laplace积分变换式代入(24)中，求得松弛模量的像空间函数如公式(27)，然后通过laplace逆变换完成时域求解如式(28)。

由(25)和(28)可求得蠕变柔量与松弛模量像空间解的数学关系：

通过laplace逆变换，可求得蠕变柔量与松弛模量时域解的数学关系：

目前，针对沥青及沥青混合料黏弹参数的换算研究主要集中蠕变柔量和松弛模量的换算，首先通过试验测得沥青的蠕变柔量，然后将蠕变柔量在时域内离散，利用蠕变柔量与松弛模量的时域下的卷积关系如式(29)，计算得到松弛模量。首先是将卷积公式离散成有限个积分公式如下式：

将松弛模量e(τ)看成区间[ti-1，ti]之间的中值，即方程(24)即简化为下式：

然后应用数值求积公式，如gauss积分，求得每个区间的积分值，然后求解得到松弛模量值；

应用laplace数值反演方法是首先得到蠕变柔量prony级数的laplace变换式，然后通过方程(4)完成laplace域松弛模量表达式的求解，如下式：

其目标就是采用laplace数值逆变换的方法，完成上式的求解；

三、松弛模量时域值求解

这里，采用type-a的蠕变柔量数据，对比分析laplace方法(ft法)与卷积方法计算的精度和效率，卷积方法中的积分采用精度较高的gauss积分法，计算效率采用计算总时长来评价，计算的精度采用平均误差来评价，平均误差的计算公式如下：

方程中，n是数据总量，ya、yb分别为两种方法的计算结果；

如图3所示，两种方法计算结果相似，通过计算分析，两种方法的计算结果平均误差误差小于10^-5，在240个时间点的计算中，卷积方法的所用时长为0.204s，而ft方法的计算时长为0.008s，ft方法比卷积方法快近25倍，这为沥青路面层状黏弹体系力学长时间计算节约了大量的时间；

在ft方法中，参数m值会影响laplace反演的精度和效率，普遍认为m值取5～15，这里将m的值取为5、8、10及15。如图4-6所示和表5所示，m值的大小对计算精度影响很小，对计算效率存在一定的影响，在计算过程可选取较小的m值。

表5不同m值下反算时长