After some time now, I have the first results to present. First the bad news: We are (currently) don’t have a faster solver compared to previous methods, which is rather untypical for utilizing GPUs. On the bright side, we can save 20-26% of memory used by a $\mathcal{H}^2$-matrices at the expense of 20-55% of increased computation time in the $\mathcal{H}^2$-matrix vector multiplication, decreasing antiproportional with the problem size. This in term is heavily utilized in, e.g., Krylov subspace methods to solve system of equations like the discretized integrals we consider.

Read more...

With the adaptive cross approximation we have a capable tool for creating low-rank approximation of arbitary matrices, thus saving memory storage and computation time [4, 5]. We are going to use this method to reduce the rank of the left factor of the

approximation of Green’s second identity $A_{t} \in \mathbb{R}^{\hat{t} \times 2K}$, $t \in \mathcal{T}_{\mathcal{I}}$ and K given like in the post mentioned before, and construct an approximation

$$A_{t} \approx L_{t}U_{t} = V_{t}P_{t}A_{t}, \qquad L_{t} \in \mathbb{R}^{\hat{t} \times k}, U_{t} \in \mathbb{R}^{k \times 2K}, P_{t} \in \mathbb{R}^{k \times \hat{t}} \text{ and } V_{t} := L_{t}(P_{t}L_{t})^{-1},$$

with $k \in \{ 1, \cdots, 2K \}$. Embedding this approximation in Green’s second identity leaves us with

$$G|_{\hat{t} \times \hat{s}} \approx A_{t} B_{ts}^{*} \approx L_{t}U_{t}B_{ts}^{*} = L_{t}(B_{ts}U_{t}^{*})^{*}.$$

In other words, we have reduced the rank of the approximation of Green’s second identity from $2K$ to $k$ [1].

Read more...

The low-rank approximation we prepared in this section still offers room for improvment [1]: the matrices $A_{t, +}$ and $A_{t, -}$ describe the influence of Neumann and Dirichlet values of $g(\cdot, y)$ on $\partial \omega_{t}$, respectively. Using the Pioncaré-Steklov operator, we can construct the Neumann from the Dirichlet values. In other words we can gain $A_{t, +}$ from $A_{t, -}$ and thus the rank of $A_{t} = \left ( A_{t, +} A_{t, -} \right ) \in \mathbb{R}^{\hat{t} \times 2K}$ is $K$.

Though it’s possible to discard $A_{t, +}$ explicitly by approximating the Poincaré-Steklov operator, this approach needs to solve an integral equation on the boundary $\partial \omega_{t}$ and a high accuraccy for the sake of preserving the exponential convergence of the quadrature approximation [1].

Because of their implactions on performance, we choose an implicit approach: we assume that the properties mentioned above holds true, e.g. the rank of $A_{t}$ is lower than the number of columns, and we can use an algebraic procedure to approximate $A_{t}$ by a lower-rank matrix. The adaptive cross approximation [1, 5] is an interesting method for this situation since it allows us to construct an algebraic interpolation operator we can use to approximate matrix blocks based on a small number of their entries.

Read more...

To receive a data-sparse approximation of a discretized matrix, which is both fast and reliable, applying cross approximation to the matrix directly would lead either to high runtime complexity of full pivoting or to potentially unreliable results of for example heurisitc pivoting strategies and error estimators.

For this reason, in a first step, we use Green’s second identity to gain a degenerate approximation of the kernel function and thus a $\mathcal{H}$-matrix representation of the Galerkin discretization [1].

Read more...

As we already mentioned, matrices resulting from integral equations are typically non-sparse, thus we consider data-sparsity: We are interested in an alternative representation or approximation of the matrix requiring just a small number of coefficients. Hierarchical matrix methods [8], [9] achieve this by dividing the matrix into submatrices approximating those who hold an admissibility condition or storing them in standard two-dimensional arrays if they are small enough.

Read more...

In this blog I will collect ideas for my master thesis with the working title “An efficient solver for the Green cross approximation method on GPUs” based on the article Approximation of integral operators by Green quadrature and nested cross approximation [1].

We consider integral equations e.g. of the form

$$\int\limits_{\Omega} g(x, y)\ u(y)\ dy = f(x) \qquad \text{for all } x \in \Omega \subseteq \mathbb{R}^d.$$

with $d \in \{ 2, 3 \}$ and

$$g(x, y) := \begin{cases} -\frac{1}{2 \pi} \log{\left \Vert x - y \right \Vert_2} & \text{if } x \neq y, \ d = 2 \\ \\ \frac{1}{4 \pi \left \Vert x - y \right \Vert_2} & \text{if } x \neq y, \ d = 3 \\ \\ 0 & \text{otherwise} \end{cases} \qquad \text{for all } x, y \in \Omega$$

as an sample kernel function for a fundamental solution of the negative Laplace operator $-\Delta$, e.g. occurring in the interpretation of Poisson’s equations in integral formulation.

Read more...