## 2 Stochastic Green's Function Method for 3D Ground Motion Simulation and Its Improvement
### Small Earthquakes in Subfaults Superimposed to Form the Large Earthquake  
After the P wave, SV wave, and SH wave of the small earthquakes in subfaults are independently generated, the time delay of their propagation in the crust and the rupture propagation on the fault plane is considered, and the three components of the waves are superimposed to generate the P wave, SV wave, and SH wave of the large earthquake, respectively. In this study, a ground motion synthesis method based on the scaling law is used [34; 35], and the ground motion of the large earthquake U\({}_{\text{K}}\)(t) at the bedrock is expressed as:  
\[U_{\text{K}}(t)=\sum_{i=1}^{NL}\sum_{j=1}^{NW}f(t)*u^{K}_{ij}\Big{(}t-t^{K}_{ ij}\Big{)} \tag{8}\]  
where\[f(t)=\delta(t)+\frac{a}{n!(1-\exp(-a))}\times\sum_{k=1}^{(ND-1)n!}\biggl{\{}\exp \biggl{[}-\frac{a(k-1)}{(ND-1)n!}\biggr{]}\cdot\delta\biggl{[}t-\frac{(k-1)\tau_ {ij}}{(ND-1)n!}\biggr{]}\biggr{\}} \tag{9}\]  
where \(\delta(t)\) is Dirac function; \(a\) is taken as 1.0; \(\tau_{ij}\) is the rise time of the _ij_ subfault; NL, _NW_, and _ND_ represent the division number of the length, width, and slip of the fault, respectively; \(n^{\prime}\) is the further subdivision number of the slip; and \(t_{ij}^{K}\) is the time delay when the ground motion of the _ij_ subfault reaches the observation point. In order to avoid the occurrence of dominant frequencies due to the close interval between the rupture start times of subfaults, a random value is set for the rupture start time of each subfault [36]. The time delay is calculated according to Equation (10):  
\[t_{ij}^{K}=\frac{r_{ij}}{v}+\frac{\xi_{ij}}{V_{\rm R}}+\frac{(\lambda-0.5) \Delta W}{V_{\rm R}} \tag{10}\]  
where, \(\xi_{ij}\) is the distance from the rupture starting point on the fault plane to the center of the _ij_ subfault; \(V_{\rm R}\) is the velocity of rupture propagation on the fault plane; \(\lambda\) is a random number between 0 and 1; \(\Delta W\) is the width of the subfault. The rise time \(\tau\) is generally calculated according to the commonly used Equation (11) [8]:  
\[\tau=0.5\frac{W}{V_{\rm R}} \tag{11}\]  
However, the rise time at a point on the fault plane is related to its position. Near the rupture starting point, the rise time is longer, while near the edge of the fault, the rise time is shorter. Equation (11) is only approximately true for the long, narrow faults [8]. According to the numerical simulation results of Day [8], for faults with various aspect ratios, there is a law that the longest rise time occurs at the rupture starting point, which equals approximately the width of the fault divided by the propagation velocity of the rupture. As the calculation point gradually moves away from the rupture starting point and approaches the edge of the fault, the rise time gradually shortens, and this attenuation law approximates a parabola. Therefore, in order to improve the accuracy of the stochastic Green's function method for ground motion simulation, based on the above laws, an improved approximate expression for the rise time on the rectangular fault considering that the rupture starting point may be at any position and the aspect ratio may be arbitrary is proposed in this study:  
\[\tau_{ij}=\frac{W}{V_{\rm R}}\biggl{[}1-\biggl{(}\frac{x_{ij}-x_{\rm R}}{C_{ 1}}\biggr{)}^{2}\biggr{]}\left[1-\biggl{(}\frac{y_{ij}-y_{\rm R}}{C_{2}} \biggr{)}^{2}\right] \tag{12}\]  
where, (\(x_{\rm R}\), \(y_{\rm R}\)) and (\(x_{ij}\), \(y_{ij}\)) stand for the local coordinates of the rupture starting point and the center of the _ij_ subfault on the fault plane, respectively. The local coordinate system of the fault plane is \(oxy\) presented in Figure 1, \(o\) is the upper left corner of the fault plane, \(x\) and \(y\) are the direction of fault length and width, respectively.