地震波数值模拟

利用波动方程对地震波传播进行数值模拟，利用MATLAB实现。同时还可研究不同速度及地层密度时地震波的传播。

clc
clear
dx=10;                  %差分网格10m*10m
dt=1/1000;              %取样间隔1ms
f0=10;                  %频率
c=2000;
u_0=zeros（400）;         %前一时刻的波场，相当于0时刻波场
u_1=zeros（400）;         %当前时刻的波场，相当于1时刻波场，
u_x=zeros（400）;         %对当前时刻的波场x方向的二阶偏导矩阵
u_z=zeros（400）;         %对当前时刻的波场z方向的二阶偏导矩阵
u_2=zeros（400）;         %后一时刻的波场，待求波场
f=zeros（400）;           %震源矩阵
u_3=zeros（200）;         %无反射波场
u_4=zeros（1000400）;
u_5=[]; 
z=1;
for t=0:0.001:0.5         %时间共1s，每1ms取样一次，共计1000个样本     
   a=pi*f0*（t-0.1）;     %对震源进行离散采样 
   s=（1-2*a.^2）.*exp（-a.^2）*100;
   f（200200）=s;          %震源更新
   u_1（200200）=s;
   for n=1:398      
        for k=1:398
            u_x（nk+1）=（u_1（nk+2）-2*u_1（nk+1）+u_1（nk））/（dx^2）;
            %对x方向求偏导，行不变，列变 
            u_z（n+1k）=（u_1（n+2k）-2*u_1（n+1k）+u_1（nk））/（dx^2）;
            %对z方向求偏导，列不变，行变
        end
   end
   
   for n=1:400                  
       for k=1:400
           u_2（nk）=（c^2）*（dt^2）*（u_x（nk）+u_z（nk）-f（nk））+2*u_1（nk）-u_0（nk）;
       end
   end
   
   for n=1:200                     
       for k=1:200
           u_3（nk）=u_1（n+100k+100）;
       end
   end
%    for i=1:400
%        u_4（zi）=u_1（100i）;
%    end
%    z=z+1;
% u_4（z）=u_1（100200）;
% u_5（z）=u_1（150200）;
% z=z+1;
   imagesc（u_3）
   set（gca‘XAxisLocation‘‘top‘）
   set（gca‘xticklabel‘{‘0‘;‘500‘;‘1000‘;‘1500‘;‘2000‘}‘xtick‘[0 50 100 150 200]）;
   set（gca‘yticklabel‘{‘0‘;‘500‘;‘1000‘;‘1500‘;‘2000‘}‘ytick‘[0 50 100 150 200]）;
   xlabel（‘位移/m‘）
   ylabel（‘深度/m‘）
   colorbar;
   drawnow
   u_0=u_1;            %更新时刻
   u_1=u_2; 
end
%  imagesc（u_4）
%   colorbar;
%   drawnow
%   title（‘距离/10m‘）
%   ylabel（‘时间/ms‘）
% t=0:0.001:1;
% figure（1）
% plot（tu_4）;
% title（‘匀速模型点（100200）处波形‘）;
% xlabel（‘x轴‘）
% ylabel（‘z轴‘）
% figure（2）
% plot（tu_5）;
% title（‘匀速模型点（150200）处波形‘）;
% xlabel（‘x轴‘）
% ylabel（‘z轴‘）

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       2100  2018-06-08 09:57  modelling.m

     文件        239  2018-05-29 10:32  source.m

----------- ---------  ---------- -----  ----

                 2339                    2