# [SSS08-P04] O(N) methods for spatiotemporal BIEM

Keywords:Spatiotemporal BIEM, O(N) method, Dynamic Rupture of Faults, H-matrix

Recent progressions of seismic inversions has shown off the numerous results which cannot be explained by the ordinary source modelings, as seen in Tohoku 2011 and Kumamoto 2016. Particularly, the hierarchical features of fault sources have been clearly detected, e.g. hierarchical asperity distribution of Tohoku-Oki [Ide Aochi 2013] and the scaling of critical slip distance [Mikumo et al. 2003]. To get the whole descriptions, i.e. the theory of earthquakes, we must exceedingly develop the source modelings including those hierarchical behaviors.

The disturbance is the numerical cost. Boundary integral equation method (BIEM) has been widely used for the source modelings to resolve the highly nonlinear boundary conditions of the fault, but it is quite time-consuming method [Ando 2016]. The original BIEM needs the cost O(N^2L) (N:fault unit number, L fault length). For example, to resolve M4 (km) events on the M9 fault (400km), we need a Peta-scale simulater! Our aim is O(N) (theoretical fastest) algorithm to solve this problematic situation. We combine the fast domain partitioning method (FDPM)[Ando et al. 2009, Ando2016] and H-matrix method, then construct O(N) method for hyperbolic equations (FDPM=H-matrix). Furthermore, this O(N) method is principally faster than the corresponding finite element modelings.

In the presentation, we talk about how to apply the H-matrix method to elastic equations and the performance evaluations of the implemented algorithm. (i) How to achieve O(N); The elastic equations have the singular wave front, thus it is known that the naive applications of H-matrix method cannot achieve O(N) [Yoshikawa and Yamamoto 2015]. FDPM provides the key idea to solve this problem. The implementation is based on the physical fact that the kernel is regular along the ray, although the kernel is singular across the ray [Aki and Richards 1980]. Thus we can safely apply the H-matrix method if we can define the ray coordinate on the kernels. This abstract idea can be mostly implemented by the Adaptive cross approximations on Front domain (Domain F in the terminology of [Ando 2016]) and Tensor Cross Approximations on Near-field and Static regions (Domain I and S). Some approximations of causality is also required and we clear the problem analytically. (ii) Performance evaluations; We show the accuracy and cost of the proposed algorithm based on some case studies. FDPM=H-matrix has some crucial approximations of time-directions, thus we carefully review it and discuss the accuracy. The accuracy is quite good (with 0.3 percent error in some cases) if we sufficiently resolve the dynamic rupture. However, because of the causality conditions, O(N) looks not to be achieved rigorously in some cases and the cost reduction looks to draw back to O(N^1.5), thus we also discuss the alternative approximations to achieve O(N) even in these situations. The result of cost reduction is also discussed in the presentation.

If we use this FDP=H-matrix method, we can extensively study the complex fault modelings and examine the theoretical hypothesis of the source parameters as already referred. Because of the numerical costs, nonplanar features of faults(for example, damage zones, branchings and fractal fault roughness) are not so studied in the current stages compared with planar fault modelings. Those features are crucially related to the source parameters [Andrews 1976], thus this method will finally contributes to the study of observed source parameters of faults.

The disturbance is the numerical cost. Boundary integral equation method (BIEM) has been widely used for the source modelings to resolve the highly nonlinear boundary conditions of the fault, but it is quite time-consuming method [Ando 2016]. The original BIEM needs the cost O(N^2L) (N:fault unit number, L fault length). For example, to resolve M4 (km) events on the M9 fault (400km), we need a Peta-scale simulater! Our aim is O(N) (theoretical fastest) algorithm to solve this problematic situation. We combine the fast domain partitioning method (FDPM)[Ando et al. 2009, Ando2016] and H-matrix method, then construct O(N) method for hyperbolic equations (FDPM=H-matrix). Furthermore, this O(N) method is principally faster than the corresponding finite element modelings.

In the presentation, we talk about how to apply the H-matrix method to elastic equations and the performance evaluations of the implemented algorithm. (i) How to achieve O(N); The elastic equations have the singular wave front, thus it is known that the naive applications of H-matrix method cannot achieve O(N) [Yoshikawa and Yamamoto 2015]. FDPM provides the key idea to solve this problem. The implementation is based on the physical fact that the kernel is regular along the ray, although the kernel is singular across the ray [Aki and Richards 1980]. Thus we can safely apply the H-matrix method if we can define the ray coordinate on the kernels. This abstract idea can be mostly implemented by the Adaptive cross approximations on Front domain (Domain F in the terminology of [Ando 2016]) and Tensor Cross Approximations on Near-field and Static regions (Domain I and S). Some approximations of causality is also required and we clear the problem analytically. (ii) Performance evaluations; We show the accuracy and cost of the proposed algorithm based on some case studies. FDPM=H-matrix has some crucial approximations of time-directions, thus we carefully review it and discuss the accuracy. The accuracy is quite good (with 0.3 percent error in some cases) if we sufficiently resolve the dynamic rupture. However, because of the causality conditions, O(N) looks not to be achieved rigorously in some cases and the cost reduction looks to draw back to O(N^1.5), thus we also discuss the alternative approximations to achieve O(N) even in these situations. The result of cost reduction is also discussed in the presentation.

If we use this FDP=H-matrix method, we can extensively study the complex fault modelings and examine the theoretical hypothesis of the source parameters as already referred. Because of the numerical costs, nonplanar features of faults(for example, damage zones, branchings and fractal fault roughness) are not so studied in the current stages compared with planar fault modelings. Those features are crucially related to the source parameters [Andrews 1976], thus this method will finally contributes to the study of observed source parameters of faults.