一种基于地层信息的场地卓越周期的精简计算方法与流程

本发明属于岩土工程勘察技术领域，具体涉及一种基于地层信息的场地卓越周期的精简计算方法。

背景技术：

卓越周期是指与地基土层发生共振作用的地震波的一个谐波分量。目前，确定场地卓越周期的方法大致分为直接测定法和波速法两大类。无论哪一类方法，都需要输入除地层信息之外的计算参数，甚至需要通过调试这些计算参数，才能得到合适的计算结果；且计算精度不高，或计算过程繁杂。

技术实现要素：

本发明的目的是提供一种基于地层信息的场地卓越周期的精简计算方法，计算过程简单、计算精度高。

本发明所采用的技术方案是，一种基于地层信息的场地卓越周期的精简计算方法，具体按照以下步骤实施：

步骤1、根据岩土工程的勘查结果，确定场地地层数n、地层密度ρi、剪切速度vi及地层厚度hi；其中，i＝1，2，......，n；

步骤2、计算场地卓越周期估计值t0；

步骤3、计算地震波持时tl和时间采样间隔δt；

步骤4、计算地面地震波持时的采样点数nt、地震波在地层i中传播时间的采样点数nti、地面地震波的傅立叶周期谱的采样点数nt；

步骤5、计算地震波在地层i与地层i+1交界界面处的反射系数ri；

步骤6、计算地面地震波时程u0^d；

步骤7、对地震波时程u0^d按照时间进行衰减处理；

步骤8、对地震波时程u进行离散傅立叶变换(dft)，得到fu(k)；

步骤9、判断振幅谱|fu(k)|是否有极大值，若否，将nt扩大20％，返回步骤8；若是，确定场地卓越周期t0。

本发明的特点还在于：

步骤2中，场地卓越周期估计值t0，由如下公式计算得到：

t0＝4(h1/v1+h2/v2+……+hn/vn)(1)。

步骤3中，地震波持时tl，由如下公式计算得到：

tl≥10t0(2)

时间采样间隔δt，由如下公式(3)计算得到：

δt<0.02min(h1/v1，h2/v2，……，hn/vn)(3)。

步骤4中，地面地震波持时的采样点数nt、地震波在地层i中传播时间的采样点数nti、地面地震波的傅立叶周期谱的采样点数nt，由如下公式计算得到：

nt＝int(tl/δt)(4)

nti＝int(hi/vi/δt)(5)

nt＝int(t0/δt)(6)

其中，int是取整数函数，i＝1，2，......，n。

步骤5中，反射系数ri，由如下公式计算得到：

ri＝(ρi+1vi+1-ρivi)/(ρi+1vi+1+ρivi)，其中，i＝1，2，......，n-1；(7)

rn＝1，其中，i＝n；(8)。

步骤6中，地震波时程u0^d，由如下公式计算得到：

u0^d(j)＝-u1^p(j-nt1)

u1^d(j)＝(1+r1)u0^d(j-nt1)-r1u2^p(j-nt2)

u1^p(j)＝r1u0^d(j-nt1)+(1-r1)u2^p(j-nt2)

u2^d(j)＝(1+r2)u1^d(j-nt2)-r2u3^p(j-nt3)

u2^p(j)＝r2u1^d(j-nt2)+(1-r2)u3^p(j-nt3)

......

ui^d(j)＝(1+ri)ui-1^d(j-nti)-riui+1^p(j-nti+1)

ui^p(j)＝riui-1^d(j-nti)+(1-ri)ui+1^p(j-nti+1)

......

un-1^d(j)＝(1+rn-1)un-2^d(j-ntn-1)-rn-1un^p(j-ntn)

un-1^p(j)＝rn-1un-2^d(j-ntn-1)+(1-rn-1)un^p(j-ntn)

un^p(j)＝un^d(j-ntn)

其中，u0^d是地面地震波时程；ui^d、分别是地层i下、上界面的地震波时程，i＝1，2，......，n-1；u0^d、ui^d、均是时间序列，j表示时间样点序号，j＝0，1，......，nt-1，j＝0表示时间等于零。

步骤7中，按照如下公式对地震波时程u0^d按照时间进行衰减处理：

u(j)＝u0^d(j)exp[-βjδt](9)

其中，衰减系数β＝-lnδ/tl；δ是控制衰减系数的计算精度值，满足0<δ<0.01；j＝0，1，......，nt-1。

步骤8中，按照如下公式对地震波时程u进行离散傅立叶变换(dft)：

其中，i是虚数单位；k＝0，1，......，nt-1。

步骤9中，按照如下公式判断振幅谱|fu(k)|是否有极大值：

|fu(k)|>|fu(k-1)|且|fu(k)|>|fu(k+1)|(11)

|fu(k)|>0.5max(|fu(1)|，|fu(2)|，......，|fu(nt-1)|)(12)

如果|fu(k)|满足上述条件，则场地卓越周期由如下公式计算得到：

t0＝0.5t0+kδt(13)。

本发明的有益效果是：

本发明应用卓越周期简便算法得到卓越周期的估计值，根据该值大小设定除场地地层信息以外的所有计算参数；把场地土视为地震波信号滤波器，计算脉冲波作用下的地面地震时程，并对其进行傅立叶变换得到振幅谱；通过判断振幅谱极大值与时间(或频率)的关系，获得高精度的卓越周期计算值。

附图说明

图1是本发明基于地层信息的场地卓越周期的精简计算方法中场地地层模型；

图2是本发明基于地层信息的场地卓越周期的精简计算方法实施例1原始地面地震时程；

图3是本发明基于地层信息的场地卓越周期的精简计算方法实施例1衰减处理后的地震时程；

图4是本发明基于地层信息的场地卓越周期的精简计算方法实施例1地震时程的傅立叶振幅谱；

图5是本发明基于地层信息的场地卓越周期的精简计算方法实施例2原始地面地震时程；

图6是本发明基于地层信息的场地卓越周期的精简计算方法实施例2衰减处理后的地震时程；

图7是本发明基于地层信息的场地卓越周期的精简计算方法实施例2地震时程的傅立叶振幅谱；

图8是本发明基于地层信息的场地卓越周期的精简计算方法实施例3原始地面地震时程；

图9是本发明基于地层信息的场地卓越周期的精简计算方法实施例3衰减处理后的地震时程；

图10是本发明基于地层信息的场地卓越周期的精简计算方法实施例3地震时程的傅立叶振幅谱；

图11是本发明基于地层信息的场地卓越周期的精简计算方法实施例4原始地面地震时程；

图12是本发明基于地层信息的场地卓越周期的精简计算方法实施例4衰减处理后的地震时程；

图13是本发明基于地层信息的场地卓越周期的精简计算方法实施例4地震时程的傅立叶振幅谱；

图14是本发明基于地层信息的场地卓越周期的精简计算方法实施例5原始地面地震时程；

图15是本发明基于地层信息的场地卓越周期的精简计算方法实施例5衰减处理后的地震时程；

图16是本发明基于地层信息的场地卓越周期的精简计算方法实施例5地震时程的傅立叶振幅谱。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

本发明一种基于地层信息的场地卓越周期的精简计算方法，具体按照以下步骤实施：

步骤1、根据岩土工程的勘查结果，确定场地地层数n、地层密度ρi、剪切速度vi及地层厚度hi；其中，i＝1，2，......，n。

步骤2、计算场地卓越周期估计值t0；

场地卓越周期估计值t0，由如下公式计算得到：

t0＝4(h1/v1+h2/v2+……+hn/vn)(1)。

步骤3、计算地震波持时tl和时间采样间隔δt；

地震波持时tl，由如下公式计算得到：

tl≥10t0(2)

时间采样间隔δt，由如下公式(3)计算得到：

δt<0.02min(h1/v1，h2/v2，……，hn/vn)(3)。

步骤4、计算地面地震波持时的采样点数nt、地震波在地层i中传播时间的采样点数nti、地面地震波的傅立叶周期谱的采样点数nt；

地面地震波持时的采样点数nt、地震波在地层i中传播时间的采样点数nti、地面地震波的傅立叶周期谱的采样点数nt，由如下公式计算得到：

nt＝int(tl/δt)(4)

nti＝int(hi/vi/δt)(5)

nt＝int(t0/δt)(6)

其中，int是取整数函数，i＝1，2，......，n。

步骤5、计算地震波在地层i与地层i+1交界界面处的反射系数ri；

反射系数ri，由如下公式计算得到：

ri＝(ρi+1vi+1-ρivi)/(ρi+1vi+1+ρivi)，其中，i＝1，2，......，n-1；(7)

rn＝1，其中，i＝n；(8)。

步骤6、计算地面地震波时程u0^d；

地震波时程u0^d，由如下公式计算得到：

u0^d(j)＝-u1^p(j-nt1)

u1^d(j)＝(1+r1)u0^d(j-nt1)-r1u2^p(j-nt2)

u1^p(j)＝r1u0^d(j-nt1)+(1-r1)u2^p(j-nt2)

u2^d(j)＝(1+r2)u1^d(j-nt2)-r2u3^p(j-nt3)

u2^p(j)＝r2u1^d(j-nt2)+(1-r2)u3^p(j-nt3)

......

ui^d(j)＝(1+ri)ui-1^d(j-nti)-riui+1^p(j-nti+1)

ui^p(j)＝riui-1^d(j-nti)+(1-ri)ui+1^p(j-nti+1)

......

un-1^d(j)＝(1+rn-1)un-2^d(j-ntn-1)-rn-1un^p(j-ntn)

un-1^p(j)＝rn-1un-2^d(j-ntn-1)+(1-rn-1)un^p(j-ntn)

un^p(j)＝un^d(j-ntn)

其中，u0^d是地面地震波时程；ui^d、分别是地层i下、上界面的地震波时程，i＝1，2，......，n-1；u0^d、ui^d、均是时间序列，j表示时间样点序号，j＝0，1，......，nt-1，j＝0表示时间等于零。

步骤7、对地震波时程u0^d按照时间进行衰减处理；

按照如下公式对地震波时程u0^d按照时间进行衰减处理：

u(j)＝u0^d(j)exp[-βjδt](9)

其中，衰减系数β＝-lnδ/tl；δ是控制衰减系数的计算精度值，满足0<δ<0.01；j＝0，1，......，nt-1。

步骤8、对地震波时程u进行离散傅立叶变换(dft)，得到fu(k)；

按照如下公式对地震波时程u进行离散傅立叶变换(dft)：

其中，i是虚数单位；k＝0，1，......，nt-1。

步骤9、判断振幅谱|fu(k)|是否有极大值，若否，将nt扩大20％，返回步骤8；若是，确定场地卓越周期t0；

按照如下公式判断振幅谱|fu(k)|是否有极大值：

|fu(k)|>|fu(k-1)|且|fu(k)|>|fu(k+1)|(11)

|fu(k)|>0.5max(|fu(1)|，|fu(2)|，......，|fu(nt-1)|)(12)

如果|fu(k)|满足上述条件，则场地卓越周期由如下公式计算得到：

t0＝0.5t0+kδt(13)。

实施例1

步骤1、根据岩土工程的勘查结果，确定场地地层数n、地层密度ρi、剪切速度vi及地层厚度hi；其中，i＝1，2，......，n；场地地层模型如图1所示；地层参数如表1所示；

n＝3

ρ1＝2041v1＝200h1＝4

ρ2＝2143v2＝300h2＝4

ρ3＝2041v3＝200h3＝12

步骤2、计算场地卓越周期估计值t0；

t0＝4×(4/200+4/300+12/200)＝0.3733

步骤3、计算地震波持时tl和时间采样间隔δt；

tl＝10t0＝10×0.3733＝3.7333

δt＝0.0001<0.02min(4/200，4/300，12/200)＝0.01333

步骤4、计算地面地震波持时的采样点数nt、地震波在地层i中传播时间的采样点数nti、地面地震波的傅立叶周期谱的采样点数nt；

nt＝int(tl/δt)＝int(3.733/0.0001)＝37333

nt1＝int(hi/vi/δt)＝int(4/200/0.0001)＝200

nt2＝int(hi/vi/δt)＝int(4/300/0.0001)＝133

nt3＝int(hi/vi/δt)＝int(12/200/0.0001)＝600

nt＝int(t0/δt)＝int(0.3733/0.0001)＝3733

步骤5、计算地震波在地层i与地层i+1交界界面处的反射系数ri；

r1＝(2143×300－2041×200)/(2143×300+2041×200)＝0.2233

r2＝(2041×200－2143×300)/(2041×200+2143×300)＝-0.2233

r3＝1

步骤6、计算地面地震波时程u0^d；

将nt，nt1，nt2，nt3以及r1，r2，r3等相关参数代入时间递推方程式，可得到地面地震时程u0^d(j)，j＝0，1，......，37332(如图2所示)。

步骤7、对地震波时程u0^d按照时间进行衰减处理；

取δ＝0.005，计算衰减系数：

β＝－lnδ/tl＝－ln0.005/3.7333＝1.4192

按照式u0^d(j)exp[-βjδt]进行衰减处理，得到衰减处理后的地震时程u(j)，j＝0，1，......，37332。(如图3所示)

步骤8、对地震波时程u进行离散傅立叶变换(dft)，得到fu(k)；

对地震波时程u进行离散傅立叶变换(dft)，得到fourier谱fu(k)，k＝0，1，......，nt-1，振幅谱|fu|如图4所示。

步骤9、判断振幅谱|fu(k)|是否有极大值，若否，将nt扩大20％，返回步骤8；若是，确定场地卓越周期t0；

按照k＝3731，3730，3729，......递减次序，寻找振幅谱|fu(k)|的极大值。经过计算可知，k＝2069，计算可得场地卓越周期：

t0＝0.5t0+kδt＝0.5×0.3733+2069×0.0001＝0.3935

卓越周期估计值、卓越周期计算值及主要计算参数详见表2。

实施例2

某2层土场地20m，地层参数如表3所示。按照上述实施例1步骤计算所得的原始地面地震时程、衰减处理后的地震时程及其傅立叶振幅谱如图5、图6及图7所示，卓越周期估计值、卓越周期计算值及主要计算参数详见表4。由表可见，尽管卓越周期估计值、主要计算参数与实施例2相同，但是卓越周期计算值却不相同，分别为0.3935s和0.4027s。

实施例3

某5层土场地10.7m，地层参数如表5所示。按照上述实施例1步骤计算所得的原始地面地震时程、衰减处理后的地震时程及其傅立叶振幅谱如图8、图9及图10所示，卓越周期估计值、卓越周期计算值及主要计算参数详见表6。由表可见，卓越周期的计算值与估计值接近，分别为0.2092s和0.2345。

实施例4

某场地土厚126m，地层较多，共12层，各个土层参数见表7所示。按照上述实施例1步骤计算所得的原始地面地震时程、衰减处理后的地震时程及其傅立叶振幅谱如图11、图12及图13所示，卓越周期估计值、卓越周期计算值及主要计算参数详见表8。地面地震时程比较复杂(见图11)，卓越周期达1.3195s。

实施例5

某场地土厚250m，地层多达20层，各个土层参数见表9所示。按照上述实施例1步骤计算所得的原始地面地震时程、衰减处理后的地震时程及其傅立叶振幅谱如图14、图15及图16所示，卓越周期估计值、卓越周期计算值及主要计算参数详见表10。由于场地分层多且厚度大，因此地面地震时程非常复杂(见图14)，卓越周期大，其值达2.0319s。

本发明基于地层信息的场地卓越周期的精简计算方法的优点为：本发明应用卓越周期简便算法得到卓越周期的估计值，根据该值大小设定除场地地层信息以外的所有计算参数；把场地土视为地震波信号滤波器，计算脉冲波作用下的地面地震时程，并对其进行傅立叶变换得到振幅谱；通过判断振幅谱极大值与时间(或频率)的关系，获得高精度的卓越周期计算值。

表1实施例1地层参数

表2实施例1计算参数

表3实施例2地层参数

表4实施例2计算参数

表5实施例3地层参数

表6实施例3计算参数

表7实施例4地层参数

表8实施例4计算参数

表9实施例5地层参数

表10实施例5计算参数