Split-Step Fourier Approximation

The extrapolation wavenumber defined in equations 14 and 15 cannot be implemented purely in the Fourier domain due to the presence of mixed-domain fields (i.e. a function of both $\xi_1$ and $k_\xi_1$simultaneously). This can be addressed using an extended version of the split-step Fourier approximation Stoffa et al. (1990), a popular approach that uses Taylor expansions to separate $k_\xi_3$ into two parts: $k_\xi_3\approx k_\xi_3^{PS} + k_\xi_3^{SSF}$. Wavenumbers $k_\xi_3^{PS}$and $k_\xi_3^{SSF}$ represent a pure Fourier domain phase-shift and a mixed $\omega-\mathbf{x}$ domain split-step correction, respectively.

The phase-shift term is given by,

$$k_\xi_3^{PS} = - b_1 k_\xi_1 - b_2 k_\xi_2 + i b_3 \pm \le... ... i\,k_\xi_1 + b_9 i\,k_\xi_2 - b_{10}^2 \right]^{\frac{1}{2}},$$ (16)

where $b_i=b_i(\xi_3)$ are reference values of $a_i = a_i(\xi_1,\xi_2,\xi_3)$. The split-step approximation is developed by performing a Taylor expansion about each coefficient a_i and evaluating the results at stationary reference values b_i. The stationary values of $k_\xi_1$ and $k_\xi_2$ are assumed to be zero. This leads to a split-step correction of,  

$$k_\xi_3^{SSF}= \left. \frac{\partial k_\xi_3}{\partial a_3} ... ...3}{\partial a_{10}} \right\vert _0 \left(a_{10}-b_{10} \right),$$ (17)

where ``'' denotes "with respect to a reference medium". The partial differential expressions in equation 17 are,

$$\left. \frac{\partial k_\xi_3}{\partial a_3} \right\vert _0 ... ...ac{b_{10}}{\sqrt{b_{10}^2 \, \omega^2 - b_{10}^2 }}, \;\;\;\;\;$$ (18)

resulting in the following split-step Fourier correction,

$$k_\xi_3^{SSF} = i\,b_3 \, \left(a_3-b_3 \right)+ \frac{b_4 \... ...eft(a_{10}-b_{10} \right)}{\sqrt{b_4\,\omega^2 - b_{10}^2 }}.$$ (19)

Note that we could use many reference media followed by interpolation similar to the phase-shift plus interpolation (PSPI) technique of Gazdag and Sguazzero (1984).

Importantly, even though there are additional a_i coefficients in the dispersion relationship, these can be made smooth through mesh regularization such that fewer sets of reference parameters are needed to accurately represent wavenumber $k_\xi_3$. In addition, situations exist where some coefficients are zero or negligible. For example, the coefficients for a weakly non-orthogonal coordinate system (i.e. ${\rm max} \vert\{g^{12},g^{13},g^{23}\}\vert << {\rm max} \{g^{11},g^{22},g^{33} \}$) within a kinematic approximation reduce to,  

$$\mathbf{a} \approx \left[ 0 \;\;\; 0 \;\;\; 0 \;\;\; \ss / ... ... \;\;\; 0 \;\;\; 0 \;\;\; n_3 / m^{33} \right]^{\mathbf{T}} \;.$$ (20)

Additional special cases are presented in Appendix A.