## 2 Semi-Analytical Solution for the Coupled Track-Ground-Tree Grid System
### Fundamental Solutions to the Dynamic Response of a Stratified Half-Space
\begin{array}{ccc}\mathbf{M_{1}}&\mathbf{O}&\mathbf{O}\\ \mathbf{O}&\mathbf{M_{2}}&\mathbf{O}\\ \mathbf{O}&\mathbf{O}&\mathbf{M_{3}}\end{array} \right]^{\mathbf{T}}\]  
By applying Equation (5), Equation (4) can be decoupled into two distinct eigenvalue problems: one describing the modes, \(\widehat{\overline{\Psi}}_{x}\) and \(\widehat{\overline{\Psi}}_{z}\), associated with normalized waves of the Rayleigh mode propagating in the \(x\)-\(z\) plane and the other describing the mode, \(\widehat{\overline{\Psi}}_{y}\), associated with normalized waves in Love mode along the \(y\) coordinate:  
\[\left(k_{R}^{2}\left[ \begin{array}{ccc}\mathbf{A_{1}}&\mathbf{O}\\ \mathbf{B_{1}^{T}}&\mathbf{A_{3}}\end{array} \right]+\left[ \begin{array}{ccc }\mathbf{G_{1}}-\omega\mathbf{M_{1}}&\mathbf{B_{1}}\\ \mathbf{O}&\mathbf{G_{3}}-\omega\mathbf{M_{3}}\end{array} \right]\right)\left\{ \begin{array}{c}\widehat{\overline{\Psi}}_{x}\\ k_{R}\widehat{\overline{\Psi}}_{z}\end{array} \right\}=\mathbf{0} \tag{6}\] \[\left(k_{L}^{2}\mathbf{A_{2}}+\mathbf{G_{2}}-\omega\mathbf{M_{2}} \right)\widehat{\overline{\Psi}}_{y}=0\]  
Herein the modes governed by Equation (6) are required to meet the orthogonal conditions [51]. The TLM is restricted to simulating wave propagation within a finite-depth domain. To address wave absorption in the underlying half-space, the perfectly matched layer (PML) technique is introduced. In the TLM formulation, the thickness \(h\) of each layer forming the PML is replaced by an adjusted thickness \(\overline{h}_{l}\)[56]. For the \(l\)-th thin layer in the PML region, the modified thickness \(\overline{h}_{l}\) is given by  
\[\overline{h}_{l}=\overline{z}_{l}-\overline{z}_{l-1}=H_{\text{PML}}\left\{ \frac{1}{N_{\text{PML}}}-\mathrm{i}\zeta\left[\left(\frac{l}{N_{\text{PML}}} \right)^{m+1}+\left(\frac{l-1}{N_{\text{PML}}}\right)^{m+1}\right]\right\} \tag{7}\]  
where \(H_{\text{PML}}\) denotes the overall thickness of the PML. The corresponding PML constants \(m\) and \(\zeta\) can be referenced from the literature [56]. Once the eigen values and associated modes specified in Equation (6) are identified, the basic solution set is constructed using a modal combination approach. The expression that describes the vertical component of displacement, expressed in the joint domain of frequency and dual wavenumbers, is detailed in Equation (8):  
\[\widetilde{\overline{G}}_{zz}^{mn}\left(k_{x},k_{y},\omega\right)=\sum_{j}^{ R-\text{mod}es}\frac{1}{k_{x}^{2}+k_{y}^{2}-k_{j}^{2}}\varphi_{zj}^{m}\varphi_{zj}^{n} \tag{8}\]  
The Green's function \(\widetilde{\overline{G}}_{zz}^{mn}\left(k_{x},k_{y},\omega\right)\) characterizes the out-of-plane response at the \(m\)-th layer induced by a harmonic point load applied at the \(n\)-th layer and is expressed in the joint frequency-wavenumber domain. For the case of \(m=n=0\), both the excitation and measurement locations are situated at the ground surface. The term \(\varphi_{zj}^{n}\) denotes the \(n\)-th entry of the vertical mode shape vector \(\hat{\overline{\psi}}_{z}\).