地震震源监测的数学建模分析

1、汶川大地震震源监测的数学建模分析1.摘要：汶川的地震给我们每一位中国人带来了巨大的伤痛，痛定思痛，我们应该为减少震后灾害做些事情。震中位置的快速确定，对第一时间展开抗震救灾起到非常重要的作用。本文应用非线性数学模型研究地震活动的非线性行为。使用数学中的非线性函数建模，模型中的参数通过多个地震观测站点接收到地震波的时间推算得到。并运用MATLAB非线性解方程求出震源;运用地形、地貌等高线图建立三维坐标图数学模型，运用MATLAB绘出地形三维立体图。根据地震方面的资料及汶川地理位置进行最佳逼近拟合。2.关键词：时间和速度尺寸 建模分析 非线性函数 最佳逼近 地震资料 3.正文1） 前言：2008年

2、5月12日发生在四川汶川的一场大地震给我们每一位中国人带来了巨大的伤痛，痛定思痛，我们应该为减少震后灾害做些事情。当地震发生时，震中位置的快速确定对第一时间展开抗震救灾起到非常重要的作用，而震中位置可以通过多个地震观测站点接收到地震波的时间推算得到。这里假定地面是一个平面，在这个平面上建立坐标系见图1。图1中给出了10个地震观测站点（AJ）的坐标位置。图1 地震观测站点示意图 2009年4月1日某时在某一地点发生了一次地震，图1中10个地震观测站点均接收到了地震波，观测数据见表1。地震灾区的地形、地貌对抗震救灾的进展会有很大影响，根据卫星遥感和飞机航拍得到的照片可以构建灾区地形、地貌

3、图。所构建的灾区局部地形、地貌等高线图见图2。注：（1）图2中等高线标注数字的单位为千米；（2）图1和图2不属于同一坐标系。图2 灾区局部地形、地貌等高线图 问题一：假定地震波在各种介质和各个方向的传播速度均相等，并且在传播过程中保持不变。请你根据表1中的数据确定这次地震的震中位置、震源深度以及地震发生的时间（不考虑时区因素，建议时间以分为单位）。问题二：请你根据图2中的局部地形、地貌等高线图，建立数学模型，绘制出相应的三维地形、地貌曲面图。问题三：查阅有关资料，了解地震波在各种介质和各个方向的传播速度问题，给出合理假设，根据表1中的数据确定这次地震的震中位置、震源深度以及地震发生的

4、时间表1 地震观测站坐标及接收地震波时间地震观测站横坐标x（千米）纵坐标y（千米）接收地震波时间A50033004月1日9时21分9秒B3002004月1日9时19分29秒C80016004月1日9时14分51秒D140022004月1日9时13分17秒E17007004月1日9时11分46秒F230028004月1日9时14分47秒G250019004月1日9时10分14秒H29009004月1日9时11分46秒I320031004月1日9时17分57秒J34001004月1日9时16分49秒2） 概述： 第一问解决方案：建立了XYZ坐标模型，得到了各检测站关于震源的时间和地点的数学表达式;

5、通过MATLAB算出震源的最佳估计值。第二问解决方案:通过地形、地貌等高线图建立XYZ坐标模型，通过MATLAB绘出三维地形、地貌曲面图。第三问解决方案：根据实际情况分析，建立数学模型。3) 模型建立与求解一问题一根据问题一所给的数据表建立XYZ坐标系，设震源坐标（x，y，z）地震波传播速度为V，各点传播时间为t-t0，建立十个点的非线性方程组如下： 利用MATLAB编写程序如下：h,yfun,key,c=fsolve('earthqu',1600;2000;100;1000;20)%非线性方程的求解主函数function q=earthqu(h)%非线性方程函数x1=500;

6、x2=300;x3=800;x4=1400;x5=1700;x6=2300;x7=2500;x8=2900;x9=3200;x10=3900;%十个点的x坐标y1=3300;y2=200;y3=1600;y4=2200;y5=700;y6=2800;y7=1900;y8=900;y9=3100;y10=100;%十个点的y坐标t1=1269;t2=1169;t3=891;t4=797;t5=706;t6=887;t7=614;t8=706;t9=1077;t10=1009;%以9：00为基准点的时间的变化量q=zeros(10,1); %H(1)=x H(2)=y p(3)=z H(4)=v

7、H(5)=t %(x-xi).2+(y-yi).2+z.2-(ti-t).2.*v.*vq=(h(1)-x1).2+(h(2)-y1).2+h(3).2-(t1-h(5).2.*h(4).*h(4);%A点非线性方程 (h(1)-x2).2+(h(2)-y2).2+h(3).2-(t2-h(5).2.*h(4).*h(4);%B点非线性方程 (h(1)-x3).2+(h(2)-y3).2+h(3).2-(t3-h(5).2.*h(4).*h(4);%C点非线性方程 (h(1)-x4).2+(h(2)-y4).2+h(3).2-(t4-h(5).2.*h(4).*h(4);%D点非线性方程 (h

8、(1)-x5).2+(h(2)-y5).2+h(3).2-(t5-h(5).2.*h(4).*h(4);%E点非线性方程 (h(1)-x6).2+(h(2)-y6).2+h(3).2-(t6-h(5).2.*h(4).*h(4);%F点非线性方程 (h(1)-x7).2+(h(2)-y7).2+h(3).2-(t7-h(5).2.*h(4).*h(4);%G点非线性方程 (h(1)-x8).2+(h(2)-y8).2+h(3).2-(t8-h(5).2.*h(4).*h(4);%H点非线性方程 (h(1)-x9).2+(h(2)-y9).2+h(3).2-(t9-h(5).2.*h(4).*h

9、(4);%I点非线性方程 (h(1)-x10).2+(h(2)-y10).2+h(3).2-(t10-h(5).2.*h(4).*h(4);%J点非线性方程运行结果为：h = 1.0e+003 * 2.3235 1.3042 -0.6993 3.1521 0.4791通过上面调试运行的结果，得震源坐标（2.3235，1.3042，-0.6993）速度是3.1521，时间为0.4791。这样我们就可以利用数学建模迅速的找到震源的具体位置。这可以尽量减小地震所带来的灾害。并方便有关部门第一时间展开抗震救灾行动。二问题二通过第二问的地形、地貌等高线图建立XYZ坐标模型，得数据表3-1如下： 3-1通

10、过已知的数据编写MATLAB程序如下：x=0 :1:12;y=0 :1:12;z=0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.12 0.18 0.16 0.13 0 0 0 0 0 0.1 0.1 0.08 0.11 0.1 0.26 0.28 0.21 0.15 0.14 0 0 0.12 0.16 0.17 0.15 0.17 0.3 0.41 0.42 0.34 0.22 0.17 0 0 0.15 0.20 0.22 0.19 0.25 0.4 0.56 0.56 0.46 0.29 0.2 0.1 0.1 0.17 0.22 0.24 0.26 0

11、.29 0.5 0.67 0.68 0.52 0.35 0.22 0.11 0.11 0.18 0.21 0.23 0.27 0.3 0.51 0.71 0.71 0.57 0.37 0.2 0.1 0.1 0.17 0.22 0.24 0.26 0.29 0.5 0.67 0.65 0.52 0.35 0.17 0 0 0.15 0.2 0.22 0.18 0.25 0.4 0.45 0.56 0.44 0.28 0.14 0 0 0.1 0.15 0.17 0.15 0.17 0.3 0.41 0.42 0.32 0.22 0 0 0 0 0.1 0.1 0.08 0.11 0.19 0.

12、25 0.25 0.2 0.14 0 0 0 0 0 0 0 0 0.1 0.15 0.16 0.12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 mesh(x,y,z) %绘原始数据图xi=linspace(0,12,200); %加密横坐标数据到50个yi=linspace(0,12,200); %加密纵坐标数据到60个xii,yii=meshgrid(xi,yi); %生成网格数据zii=interp2(x,y,z,xii,yii,'cubic'); %插值mesh(xii,yii,zii) %加密后的地貌图hold on % 保持图形xx,yy=me

13、shgrid(x,y); %生成网格数据plot3(xx,yy,z+0.1,ob) %原始数据用O汇编出：正视图3-2，俯视图3-3，仰视图3-4三图如下，这样我们就可以正确而又准确的描绘出地震周围的地貌。更方便救援工作顺利而又高效的进行。（侧视图见附录）正视图3-2俯视图3-3仰视图3-4三问题三1）分析问题三：地震波的传播速度与传播介质有关系，表现在地底深度上，纵波是推进波，地壳中传播速度为78.5千米秒，最先到达震中，又称P波，它使地面发生上下振动，破坏性较弱。横波是剪切波：在地壳中的传播速度为3.24.5千米/秒，第二个到达震中，又称S波，它使地面发生前后、左右抖动，破坏性较强。面波又

14、称L波，是由纵波与横波在地表相遇后激发产生的混合波：面波的传播速度为約為 34 千米/秒。其波长大、振幅强，只能沿地表面传播，是造成建筑物强烈破坏的主要因素。所以地震波的传播速度与各种因素有关系，地震波在实际地球内部传播示意图3-53-5弹性模量和波速这一名词描述波传播中由于传播路径上条件变化产生波前方向变化的现象。反射和折射也是光线通过透镜和棱柱时人们熟知的性质。均质各向同性的固体可由两个常数：k和来描述其弹性，两常数都可表示为单位面积的力。k是体积模量，表示不可压缩性。花岗岩：k约为27×1010达因/厘米2；水：k约为2×1010达因/厘米2。 是剪切模量，

15、表示其刚性。花岗岩：约为1.6×1010达因/厘米2；水：为0。密度为的弹性固体内，可以传播两种弹性波。P波，速度vP=(k+3/4)/。花岗岩：vP=5.5千米/秒；水：vP=1.5千米/秒。S波，速度vS=/。花岗岩：vS=3.0千米/秒；水：vS=0千米/秒。综上所诉可以推出纵波在地球内部传播速度最快。所以我们用纵波的波速解决第三问。2）解决问题：我们取纵波为已知量。代入问题一的方程中。运用MATLAB不难解的h = 1.0e+003 * 2.3156 1.3654 -0.67560.3678的具体震源坐标为（2.3156，1.3654，-0.6756）传播时间为0.35783）应用问题：所以，从现实中考虑，当震源深度为-0.6756时，破坏性较强的地震传播速度最快为4.5，照此计算，汶川地震发生瞬间立即鸣响500公里范围内的警报系统，那么可以让北川提前19 秒逃离，青川提前50秒逃离，广元提前55秒逃离。45秒足以让人从7楼跑到一楼外空地上。更重要的是：地震区的人们不会浪费判断的时间！想想多少生命可以继续，多少孩子 会有更美的明天！ 就算我们无法预测地震，那么地震发生后的瞬间不作为，那将是不可饶恕的罪孽！看看日本的临震预报，那是我们学习的榜样，学习他们的技术并不难，难的是学习他们对生命的尊 重！强烈要求地震带安置可靠的报警系统！4、最后