Applications of -Matrix Techniques in Micromagnetics

The variational model by Landau and Lifshitz is frequently used in the simulation of stationary micromagnetic phenomena. We consider the limit case of large and soft magnetic bodies, treating the associated Maxwell equation exactly via an integral operator . In numerical simulations of the resulting minimization problem, difficulties arise due to the imposed side-constraint and the unboundedness of the domain. We introduce a possible discretization by a penalization strategy. Here the computation of is numerically the most challenging issue, as it leads to densely populated matrices. We show how an efficient treatment of both and the corresponding bilinear form can be achieved by application of -matrix techniques.     

Adaptive Recompression of -Matrices for BEM

The efficient treatment of dense matrices arising, e.g., from the finite element discretisation of integral operators requires special compression techniques. In this article, we use a hierarchical low-rank approximation, the so-called -matrix, that approximates the dense stiffness matrix in admissible blocks (corresponding to subdomains where the underlying kernel function is smooth) by low rank matrices. The low rank matrices are assembled by the ACA+ algorithm, an improved variant of the well-known ACA method. We present an algorithm that can determine a coarser block structure that minimises the storage requirements (enhanced compression) and speeds up the arithmetic (e.g., inversion) in the -matrix format. This coarse approximation is done adaptively and on-the-fly to a given accuracy such that the matrix is assembled with minimal storage requirements while keeping the desired approximation quality. The benefits of this new recompression technique are demonstrated by numerical examples.     

Adaptive Geometrically Balanced Clustering of <Emphasis Type="Italic">H</Emphasis>-Matrices

In [8], a class of (data-sparse) hierarchical (-) matrices is introduced that can be used to efficiently assemble and store stiffness matrices arising in boundary element applications. In this paper, we develop and analyse modifications in the construction of an -matrix that will allow an efficient application to problems involving adaptive mesh refinement. In particular, we present a new clustering algorithm such that, when an -matrix has to be updated due to some adaptive grid refinement, the majority of the previously assembled matrix entries can be kept whereas only a few new entries resulting from the refinement have to be computed. We provide an efficient implementation of the necessary updates and prove for the resulting -matrix that the storage requirements as well as the complexity of the matrix-vector multiplication are almost linear, i.e., AMS Subject Classifications: 65F05, 65F30, 65N38, 65N50.     

-Matrix Arithmetics in Linear Complexity

For hierarchical matrices, approximations of the matrix-matrix sum and product can be computed in almost linear complexity, and using these matrix operations it is possible to construct the matrix inverse, efficient preconditioners based on approximate factorizations or solutions of certain matrix equations. -matrices are a variant of hierarchical matrices which allow us to perform certain operations, like the matrix-vector product, in ``true' linear complexity, but until now it was not clear whether matrix arithmetic operations could also reach this, in some sense optimal, complexity. We present algorithms that compute the best-approximation of the sum and product of two -matrices in a prescribed -matrix format, and we prove that these computations can be accomplished in linear complexity. Numerical experiments demonstrate that the new algorithms are more efficient than the well-known methods for hierarchical matrices.     

The Fourier Finite-element Approximation of the Lamé Equations in Axisymmetric Domains with Edges

This paper is concerned with a priori error estimates and convergence analysis of the Fourier-finite-element solutions of the Neumann problem for the Lamé equations in axisymmetric domains with reentrant edges. The Fourier-FEM combines the approximating Fourier method with respect to the rotational angle using trigonometric polynomials of degree N (N→∞), with the finite element method on the plane meridian domain of with mesh size h (h→0) for approximating the Fourier coefficients. The asymptotic behavior of the solution near reentrant edges is described by singularity functions in non-tensor product form and treated numerically by means of finite element method on locally graded meshes. For the rate of convergence of the combined approximations in is proved to be of the order      

Approximation of Integral Operators by -Matrices with Adaptive Bases

-matrices can be used to construct efficient approximations of discretized integral operators. The -matrix approximation can be constructed efficiently by interpolation, Taylor or multipole expansion of the integral kernel function, but the resulting representation requires a large amount of storage.In order to improve the efficiency, local Schur decompositions can be used to eliminate redundant functions from an original approximation, which leads to a significant reduction of storage requirements and algorithmic complexity.     

Hierarchical Tensor-Product Approximation to the Inverse and Related Operators for High-Dimensional Elliptic Problems

The class of -matrices allows an approximate matrix arithmetic with almost linear complexity. In the present paper, we apply the -matrix technique combined with the Kronecker tensor-product approximation (cf. [2, 20]) to represent the inverse of a discrete elliptic operator in a hypercube (0, 1)^dd in the case of a high spatial dimension d. In this data-sparse format, we also represent the operator exponential, the fractional power of an elliptic operator as well as the solution operator of the matrix Lyapunov-Sylvester equation. The complexity of our approximations can be estimated by (dn log ^qn), where N=n^d is the discrete problem size.     

Approximation of 1/||x−y|| by Exponentials for Wavelet Applications (Short Communication)

We discuss the approximation of by exponentials in order to apply it to the treatmentof 1/||x-y||. In the case of a wavelet basis, one has in addition the vanishing moment property, which allows to add polynomials without increasing the computational effort. This leads to the question whether an approximation of by the sum of a polynomial and an exponential part yields an improvement. We show that indeed the approximation error is remarkably reduced. The improvement depends on the interval on which is approximated.     

On the Crossing Number of Complete Graphs

Let (G) denote the rectilinear crossing number of a graph G. We determine (K₁₁)=102 and (K₁₂)=153. Despite the remarkable hunt for crossing numbers of the complete graph K_n – initiated by R. Guy in the 1960s – these quantities have been unknown forn>10 to date. Our solution mainly relies on a tailor-made method for enumerating all inequivalent sets of points (order types) of size 11. Based on these findings, we establish a new upper bound on (K_n) for general n. The bound stems from a novel construction of drawings of K_n with few crossings.     

Verification of Reduced Convergence Rates

In this short article, we recalculate the numerical example in Kíek and Neittaanmäki (1987) for the Poisson solution u=x(1–x)siny in the unit square S as . By the finite difference method, an error analysis for such a problem is given from our previous study by where h is the meshspacing of the uniform square grids used, and C₁ and C₂ are two positive constants. Let =u–u_h, where u_h is the finite difference solution, and is the discrete H¹ norm. Several techniques are employed to confirm the reduced rate of convergence, and to give the constants, C₁=0.09034 and C₂=0.002275 for a stripe domain. The better performance for arises from the fact that the constant C₁ is much large than C₂, and the h in computation is not small enough.     

Decidability in Analysis

In this paper, we introduce the notion of almost decidable predicates of real variables. The notion comes from the concept of decidability of number theoretical predicates together with the idea of effective convergence in computable analysis. The weakness of traditional definitions of decidability is discussed. The definition of almost decidable predicates of real variables is given. It is proved that some commonly used predicates on such as between computable reals and recursive open/closed sets are almost decidable, which justifies our definition of decidability in Euclidean spaces. Additionally, the relation between almost decidability and computability of reals and recursiveness of subsets of is considered, which provides a bridge to include our works here to the earlier literature on computability in analysis.     

Hierarchical Quadrature for Singular Integrals

We introduce a method for the computation of singular integrals arising in the discretization of integral equations. The basic method is based on the concept of admissible subdomains, known, e.g., from panel clustering techniques and -matrices: We split the domain of integration into a hierarchy of subdomains and perform standard quadrature on those subdomains that are amenable to it. By using additional properties of the integrand, we can significantly reduce the algorithmic complexity of our approach. The method works also well for hypersingular integrals.     

Parallel -Matrix Arithmetics on Shared Memory Systems

-matrices, as they were introduced in previous papers, allow the usage of the common matrix arithmetic in an efficient, almost optimal way. This article is concerned with the parallelisation of this arithmetics, in particular matrix building, matrix-vector multiplication, matrix multiplication and matrix inversion.Of special interest is the design of algorithms, which reuse as much as possible of the corresponding sequential methods, thereby keeping the effort to update an existing implementation at a minimum. This could be achieved by making use of the properties of shared memory systems as they are widely available in the form of workstations or compute servers. These systems provide a simple and commonly supported programming interface in the form of POSIX-threads.The theoretical results for the parallel algorithms are confirmed with numerical examples from BEM and FEM applications.     

The Singular Sources Method for an Inverse Transmission Problem

We employ the singular sources method introduced in [4] to solve an inverse transmission scattering problem for the Helmholtz equation or D, respectively, where the total field u satisfies the transmission conditions on the boundary of some domain D with some constant β. The main idea of the singular sources scheme is to reconstruct the scattered field of point sources or higher multipoles (·, z) with source point z in its source point from the far field pattern of scattered plane waves. The function (z, z) is shown to become singular for z→∂D. This can be used to detect the shape D of the scattering object.Here, we will show how in addition to reconstructions of the shape D of the scattering object, the constant β can be reconstructed without solving the direct scattering problem. This extends the singular sources method from the reconstruction of geometric properties of an object to the reconstruction of physical quantities.     

Low-rank Kronecker-product Approximation to Multi-dimensional Nonlocal Operators. Part I. Separable Approximation of Multi-variate Functions

The Kronecker tensor-product approximation combined with the -matrix techniques provides an efficient tool to represent integral operators as well as certain functions F(A) of a discrete elliptic operator A in ℝ ^d with a high spatial dimension d. In particular, we approximate the functions A ⁻¹ and sign(A) of a finite difference discretisation A∈ℝ ^{N × N} with a rather general location of the spectrum. The asymptotic complexity of our data-sparse representations can be estimated by (n ^p log ^q n), p = 1, 2, with q independent of d, where n=N ^{1/ d} is the dimension of the discrete problem in one space direction. In this paper (Part I), we discuss several methods of a separable approximation of multi-variate functions. Such approximations provide the base for a tensor-product representation of operators. We discuss the asymptotically optimal sinc quadratures and sinc interpolation methods as well as the best approximations by exponential sums. These tools will be applied in Part II continuing this paper to the problems mentioned above.     

Fast Summation of Radial Functions on the Sphere

Radial functions are a powerful tool in many areas of multi-dimensional approximation, especially when dealing with scattered data. We present a fast approximate algorithm for the evaluation of linear combinations of radial functions on the sphere . The approach is based on a particular rank approximation of the corresponding Gram matrix and fast algorithms for spherical Fourier transforms. The proposed method takes (L) arithmetic operations for L arbitrarily distributed nodes on the sphere. In contrast to other methods, we do not require the nodes to be sorted or pre-processed in any way, thus the pre-computation effort only depends on the particular radial function and the desired accuracy. We establish explicit error bounds for a range of radial functions and provide numerical examples covering approximation quality, speed measurements, and a comparison of our particular matrix approximation with a truncated singular value decomposition.     

A Generalized Kahan-Babuška-Summation-Algorithm

In this article, we combine recursive summation techniques with Kahan-Babuška type balancing strategies [1], [7] to get highly accurate summation formulas. An i-th algorithm have only error beyond 1upl and thus allows to sum many millions of numbers with high accuracy. The additional afford is a small multiple of the naive summation. In addition we show that these algorithms could be modified to provide tight upper and lower bounds for use with interval arithmetic.     

Low-rank Kronecker-product Approximation to Multi-dimensional Nonlocal Operators. Part II. HKT Representation of Certain Operators

This article is the second part continuing Part I [16]. We apply the -matrix techniques combined with the Kronecker tensor-product approximation to represent integral operators as well as certain functions F(A) of a discrete elliptic operator A in a hypercube (0,1) ^d ∈ ℝ ^d in the case of a high spatial dimension d. We focus on the approximation of the operator-valued functions A ^{− σ} , σ>0, and sign (A) for a class of finite difference discretisations A ∈ ℝ ^{N × N} . The asymptotic complexity of our data-sparse representations can be estimated by (n ^p log ^q n), p = 1, 2, with q independent of d, where n=N ^{1/ d} is the dimension of the discrete problem in one space direction.     

A Simplified Approach to the Order Conditions of Integration Methods

We present an approach to the numerical integration of ordinary differential equations based on the algebraic theory of Butcher (Math. Comp. 26, 79–106, 1972) and the -series theory of Hairer and Wanner (Computing 13, 1–15, 1974). We clarify the differences of these two approaches by equating the elementary weight functions and showing the differences of the composition rules. By interpreting the elementary weight function as a mapping from input values to output values and introducing some special mappings, we are able to derive the order conditions of several types of integration methods in a straight-forward way. The simplicity of the derivation is illustrated by linear multistep methods that use the second derivative as an input value, Runge-Kutta type methods that use the second as well as first derivatives, and general two-step Runge-Kutta methods. We derive new high stage-order methods in each example. In particular, we found a symmetric and stiffly-accurate method of order eight in the second example.