Wavefield Separation

The theory for this section is based on this paper:

Zhou, X.-Y., Chang, X., Wang, Y.-B., Wen, X.-T., You, J.-C., & Sun, C. (2022). Non-artifact vector P- and S-wave separation for elastic reverse time migration. Petroleum Science, 19(6), 2695–2710.

One can separate the P- and S-wave components of a wavefield by taking the divergence and the curl of the displacement field. An easy way to do this in Salvus is to output the gradient-of-displacement and then use the components from the resulting tensor to construct the divergence and curl manually.

The following is an example demonstrating such a wavefield separation:

2D Separation

The gradient of the displacement field u\nabla \vec{\mathbf{u}} can be computed using

u=[xuxyuxxuyyuy]=[u0u1u2u3] \begin{aligned} \nabla \vec{\mathbf{u}} &= \begin{bmatrix} \partial_x u_x & \partial_y u_x \\ \partial_x u_y & \partial_y u_y \end{bmatrix} \\ &= \begin{bmatrix} u_0 & u_1 \\ u_2 & u_3 \end{bmatrix} \end{aligned}

where the indices of uiu_i represent which component of the output field to select.

To get the divergence of the displacement field (P-wave component), simply apply

u=xux+yuy=u0+u3. \begin{aligned} \nabla \cdot \vec{\mathbf{u}} &= \partial_x u_x + \partial_y u_y \\ &= u_0 + u_3. \end{aligned}

Similarly, to get the curl of the displacement field (S-wave component), use

×u=xuyyux=u2u1. \begin{aligned} \nabla \times \vec{\mathbf{u}} &= \partial_x u_y - \partial_y u_x \\ &= u_2 - u_1. \end{aligned}
Practical Advice

The wavefield which is computed during a Salvus simulation can be opened using the WavefieldOutput.

Note that the regularly gridded wavefield output has dimensions [n_time_slices, n_components, n_x, n_y]. The subscripts for uu effectively represent which dimension in the wavefield output should be extracted along axis=1.

For example, the divergence of the wavefield output could be obtained using something anaologous to

div_u = wavefield_output[:, 0, :, :] + wavefield_output[:, 3, :, :]

3D Separation

The same general procedure applies as in the 2D case, but now considering the additional z-component.

u=[xuxyuxzuxxuyyuyzuyxuzyuzzuz]=[u0u1u2u3u4u5u6u7u8]. \begin{aligned} \nabla \vec{\mathbf{u}} &= \begin{bmatrix} \partial_x u_x & \partial_y u_x & \partial_z u_x \\ \partial_x u_y & \partial_y u_y & \partial_z u_y \\ \partial_x u_z & \partial_y u_z & \partial_z u_z \end{bmatrix} \\ &= \begin{bmatrix} u_0 & u_1 & u_2 \\ u_3 & u_4 & u_5 \\ u_6 & u_7 & u_8 \end{bmatrix}. \end{aligned}

To get the divergence of the displacement field:

u=xux+yuy+zuz=u0+u4+u8. \begin{aligned} \nabla \cdot \vec{\mathbf{u}} &= \partial_x u_x + \partial_y u_y + \partial_z u_z \\ &= u_0 + u_4 + u_8. \end{aligned}

Similarly, to get the curl:

×u=(yuzzuy)i^+(zuxxuz)j^+(xuyyux)k^=(u7u5)i^+(u2u6)j^+(u3u1)k^. \begin{aligned} \nabla \times \vec{\mathbf{u}} &= \left( \partial_y u_z - \partial_z u_y \right) \hat{\text{i}} + \left( \partial_z u_x - \partial_x u_z \right) \hat{\text{j}} + \left( \partial_x u_y - \partial_y u_x \right) \hat{\text{k}} \\ &= \left( u_7 - u_5 \right) \hat{\text{i}} + \left( u_2 - u_6 \right) \hat{\text{j}} + \left( u_3 - u_1 \right) \hat{\text{k}}. \end{aligned}