计算科学 - 在阻抗边界处模拟压力波 - 吾爱随笔录

我试图在 Matlab 中模拟从一种介质到另一种介质（例如，水到空气）的边界的压力波。到目前为止，我得到的代码主要取自维基百科的部分传输，它只会在波穿过边界时改变速度。该代码使用递归形式或使用中心差分方法找到的波动方程工作。

我对如何将密度对介质阻抗的影响结合起来感到困惑。据我了解：

Z = ρ C

$Z=\rho C$ 在哪里

Z

$Z$ 是阻抗，

ρ

$\rho$ 是密度，

C

$C$ 是声速

在边界处，反射和透射的波振幅的比例由反射和透射系数给出，它们是：

R = \frac{Z_{1} - Z_{2}}{Z_{1} + Z_{2}}, T = R = \frac{2 Z_{1}}{Z_{1} + Z_{2}}

$R=\frac{Z_1-Z_2}{Z_1+Z_2},\\ T=R=\frac{2Z_1}{Z_1+Z_2}$

目前，模拟将正确计算两种密度相同但波速不同的介质的反射波和透射波幅度；但是，我无法改变介质的密度。

鉴于密度或阻抗在波动方程中都不相关，我对它们最终如何确定反射和透射系数以及如何将它们实施到我的模拟中感到困惑。

以下是我的代码的相关部分：

   % length of the string and the grid
   L = 5;
   N = 151;
   X=linspace(0, L, N);

   h = X(2)-X(1); % space grid size
   c = 0.01; % speed of the wave for visualisation 
   tau = 0.25*h/c; % time grid size

   % form a medium with a discontinuous wave speed
   C = 0*X+c;  %this has formed a vector the same dimension as X with every entry =c 

   D=0.5*L;
   c_right = 2*c; % speed to the right of the disc
   for i=1:N
      if X(i) > D
         C(i) = c_right;
      end
   end
   % Now C = c fo x < D, and C=c_right for x > D

   K = 10; % steepness of the bump
   S = 0; % shift the wave
   f=inline('exp(-K*(x-S).^2)', 'x', 'S', 'K'); % a gaussian as an initial wave
   df=inline('-2*K*(x-S).*exp(-K*(x-S).^2)', 'x', 'S', 'K'); % derivative of f

   % wave at time 0 and tau
   U0 = 0*f(X, S, K);
   U1 = U0 - 2*tau*c*df(X, S, K);

   U = 0*U0; % current U

   % plot between Start and End
   Start=0; End=1500;

   for j=1:End

      %  fixed end points
       U(1)=0; U(N)=0;

      % finite difference discretization in time
      for i=2:(N-1)
          %this is the wave equation written algebraically with second
          %order central difference theorem and then rearranged for U(i)
         U(i) = (C(i)*tau/h)^2*(U1(i+1)-2*U1(i)+U1(i-1)) + 2*U1(i) - U0(i);
      end

      % update info, for the next iteration
      U0 = U1; U1 = U;