A log-linear time algorithm for the elastodynamic boundary integral equation method

Abstract

We present a fast and memory-efficient algorithm for transient, space-time-domain, and elastodynamic boundary-integral analysis. Associated data-sparse approximations and operations are named fast domain partitioning hierarchical matrices (FDP=H-matrices). The fast domain partitioning method (the FDPM) solves a known problem of hierarchical matrices (H-matrices) in compressing discretized elastodynamic kernel functions. A novel set of plane-wave approximations then unites the FDPM and H-matrices in an accurate analytic manner. Memory usage is $\mathcal O(N \log N)$ and computation time $\mathcal O(NM \log N)$ in our algorithm throughout one run with $N$ boundary elements and $M$ time steps. The amount of associated cost reduction is remarkable, as the memory usage and computational time have been originally $\mathcal O(N^2M)$ and $\mathcal O(N^2M^2)$, respectively, to run the orthodox time-marching implementation. Numerical experiments indicate that FDP=H-matrices achieve $\mathcal O(NM/\log N)$ times smaller memory and computation time while ensuring the accuracy of the analyses.

Figures (29)

• Figure 1: Schematic of the FDPM. A 3D elastodynamic example problem of a linear boundary is considered in the figure. a, Schematic of the domain partitioning. The panel depicts a spatiotemporal BIE that convolves $K$ and $D$ over sources $j=1,...,N$ and time $t\in (0,M\Delta t)$ for evaluating $T$ of respective receivers $i=1,...,N$. The domain of kernel $K$ is partitioned into subdomains. Domain F (the red parts) fully encloses the wavefronts of the P- and S-waves. (Fp and Fs, respectively). The separators of the subdomains are the propagation times (the travel times) of the P- and S-waves ($t_{ij}^\alpha$ and $t_{ij}^\beta$, respectively) assigned to the collocation points of receiver $i$ and source $j$. Domain I (the orange part) is in-between Fp and Fs (the P- and S-wave parts of Domain F, respectively). Domain S (the ivory part) is after Fs. b, Schematic of the separation of variables. The kernel tensor $K^I$ in Domain I separates into the spatially-varying part and time-dependent part, expressed by ($i,j$)-dependent matrices $\hat{K}^I$ and ($t$-dependent) vectors $h^I$. The kernel tensor $K^S$ in Domain S is time-invariant, expressed by an ($i,j$)-dependent matrix $\hat{K}^S$.
• Figure 1: $N$ versus the number of partitions (corresponding to the cost of Quantization) per receiver, when Quantization is singly used with $\epsilon_Q=0.1$. Lines show some asymptotic scalings.
• Figure 1: Error and rank distributions of submatrices approximated by the partially-pivoting ACA. The other settings are the same as in Fig. \ref{['FDPHfig:12']}. (Left) Error distribution in $\hat{K}^F$ for constant $\eta$. (Right) Error distribution in $\hat{K}^S$ for constant $\eta$.
• Figure 1: Error distributions in approximate kernels $K^{\rm approx}$ of different temporal ranks, in a 2D planar boundary case. Errors are quantified by the difference between $K^{\rm approx}$ and original kernel $K^{\rm original}$. Used values of the approximation parameters are $\epsilon_{Q} =\epsilon_{ACA} =\epsilon_{st}= 10^{-3}$, $l_{min} / \Delta x = 5$, and $\eta = 5.67$ and the temporal distance between the travel time and the end of Domain F is enlarged by $3\Delta x / \beta$. (Top left) Relative error of the asymptotic kernel, being one example case of the temporally first rank. (Top right) Relative error, for the case of the temporally first rank, where the temporal pivot point is set at the start of Domain S. (Bottom left) Relative error, for the case of the temporally second rank. (Bottom right) Absolute error, for the case of the temporally second rank, normalized by the radiation damping term.
• Figure 2: Schematic of H-matrices, illustrating an example case of linearly aligned structured boundary elements in a static problem, convolving $K$ and $E$ to evaluate $T$. Kernel matrix $K$ is subdivided into submatrices as the associated pairs of source clusters and receiver clusters are divided. The levels of the source-receiver clusters represent their number of divisions. The figure also shows the two division-stopping conditions: $diam<\eta dist$ for admissibly distant source- and receiver-cluster pairs and $diam<l_{min}$ for inadmissibly small ones under given parameters $\eta$ and $l_{min}$. The size $diam$ and distance $dist$ of clusters are indicated in the matrix, particularly for the above-mentioned boundary geometry, after divided by element length $\Delta x$. The low rank approximation of the kernel for an admissible leaf is also described.