用有限差分法在matlab中绘制0-100km内地震传播动图，一维波动传播模拟

以下是一个简单的基于有限差分法的地震波传播模拟的 MATLAB 代码，可以绘制地震波在1维介质中的传播动图，范围为0-100km。请注意，这只是一个简单的示例，并且可能需要进行进一步的调整和优化才能适用于实际的地震模拟。

% 定义模拟参数
dx = 1000; % 网格间距，单位为米
dt = 0.01; % 时间步长，单位为秒
nt = 1000; % 模拟时间步数
nx = 101; % 网格数

% 定义介质参数
vp = 6000*ones(nx,1); % 定义纵波速度，单位为米/秒
vs = 3500*ones(nx,1); % 定义横波速度，单位为米/秒
rho = 2500*ones(nx,1); % 定义密度，单位为千克/立方米

% 初始化波场
u = zeros(nx,nt);

% 定义震源
src = zeros(nt,1);
src(1) = 1;

% 循环计算波场
for it = 2:nt
    % 计算纵波和横波速度的平均值
    v_avg = sqrt((vp.^2+4*vs.^2)/5);
    
    % 计算应力和速度的时间差分
    stress = rho.*v_avg.^2.*(u(3:nx,it-1)-2*u(2:nx-1,it-1)+u(1:nx-2,it-1))/dx^2;
    velocity = u(2:nx-1,it-1)+dt*stress./rho(2:nx-1);
    
    % 添加震源
    velocity(51) = velocity(51)+src(it);
    
    % 计算速度的空间差分
    u(2:nx-1,it) = velocity+dt*vs(2:nx-1).^2.*(u(3:nx,it-1)-u(1:nx-2,it-1))/dx^2;
    
    % 绘制波场动图
    plot(u(:,it));
    xlim([1 nx]);
    ylim([-0.1 0.1]);
    drawnow;
end

该代码使用有限差分法来模拟纵波和横波在1维介质中的传播。在代码中，我们首先定义了模拟参数，包括网格间距、时间步长、模拟时间步数和网格数。然后定义了介质参数，包括纵波速度、横波速度和密度。接着初始化了波场，并定义了震源。最后，我们使用一个循环来计算每个时间步长的波场，并绘制波场动图。

需要注意的是，该代码中使用了简单的震源函数，仅在第一个时间步加入了一个单位震源。在实际应用中，需要根据具体情况选择合适的震源函数。此外，该代码还假设介质参数是均匀的，实际情况中介质参数可能存在空间变化，需要进行相应的调整。