Search

Low-rank approximation (LRA) techniques have become crucial tools in scientific computing in order to reduce the cost of storing matrices and compute usual matrix operations. Since standard techniques like the SVD do not scale well with the problem size N, there has been recently a growing interest for alternative methods like randomized LRAs. These methods are usually cheap, easy to implement and optimize, since they involve only very basic operations like Matrix Vector Products (MVPs) or orthogonalizations. More precisely, randomization allows for reducing the cubic cost required to perform a standard matrix factorization to the quadratic cost required to apply a few MVPs, namely O(r × N^2) operations where r is the numerical rank of the matrix. First of all, we present a new efficient algorithm for performing MVPs in O(N) operations called the Uniform FMM (ufmm). It is based on a hierarchical (data sparse) representation of a kernel matrix combined with polynomial interpolation of the kernel on equispaced grids. The latter feature allows for FFT-acceleration and consequently reduce both running time and memory footprint but has implications on accuracy and stability. Then, the ufmm is used to speed-up the MVPs involved in the randomized SVD, thus reducing its cost to O(r^2 × N) and exhibiting very competitive performance when the distribution of points is large and highly heterogeneous. Finally, we make use of this algorithm to efficiently generate spatially correlated multivariate Gaussian random variables.; Les approximations de rang faible (LRA) sont devenus des outils fondamentaux en calcul scientifique en vue de réduire les coûts liés au stockage et aux opérations matricielles. Le coût des méthodes standards comme la SVD croît très rapidement avec la taille du problème N, c'est pourquoi des méthodes alternatives comme les approches aléatoires (i.e. basées sur la projection ou la sélection de colonnes/l'échantillonnage aléatoire) se popularisent. Ces méthodes sont en général peu coûteuses et facile à implémenter et optimiser, car elles ne mettent en oeuvre que des opérations matricielles simples comme des produits ou des orthogonalisations. Plus précisemment, les LRA aléatoires permettent de réduire le coût cubique en N des méthodes standards de factorisation au coût quadratique nécessaire à la réalisation de quelques produits matrices vecteurs, i.e., O(r × N^2) opérations où r est le rang numérique de la matrice. Dans un premier temps, nous présentons un algorithme efficace pour réaliser des MVPs en O(N) opérations, que nous appelons Uniform FMM (ufmm). Il est basé sur la combinaison d'une représentation hiérachique d'une matrice noyau et l'interpolation polynomiale du noyau associé sur une grille régulière (uniforme). Cette dernière propriété permet une accélération par FFT réduisant ainsi le temps de calcul et la consommation mémoire mais a des répercussions sur la précision et la stabilité de l'algorithme. Ensuite, la ufmm est utilisée pour accélerer les MVPs intervenants dans la SVD aléatoire (i.e., SVD par projection aléatoire) diminuant son coût asymptotique à O(r^2 × N). La méthode est particulièrement compétitive pour des distributions de points hétérogènes. Enfin, nous utilisons cet algorithme pour générer des champs de variables Gaussiens de manière efficace.

Low-rank approximation (LRA) techniques have become crucial tools in scientific computing in order to reduce the cost of storing matrices and compute usual matrix operations. Since standard techniques like the SVD do not scale well with the problem size N, there has been recently a growing interest for alternative methods like randomized LRAs. These methods are usually cheap, easy to implement and optimize, since they involve only very basic operations like Matrix Vector Products (MVPs) or orthogonalizations. More precisely, randomization allows for reducing the cubic cost required to perform a standard matrix factorization to the quadratic cost required to apply a few MVPs, namely O(r × N^2) operations where r is the numerical rank of the matrix. First of all, we present a new efficient algorithm for performing MVPs in O(N) operations called the Uniform FMM (ufmm). It is based on a hierarchical (data sparse) representation of a kernel matrix combined with polynomial interpolation of the kernel on equispaced grids. The latter feature allows for FFT-acceleration and consequently reduce both running time and memory footprint but has implications on accuracy and stability. Then, the ufmm is used to speed-up the MVPs involved in the randomized SVD, thus reducing its cost to O(r^2 × N) and exhibiting very competitive performance when the distribution of points is large and highly heterogeneous. Finally, we make use of this algorithm to efficiently generate spatially correlated multivariate Gaussian random variables.; Les approximations de rang faible (LRA) sont devenus des outils fondamentaux en calcul scientifique en vue de réduire les coûts liés au stockage et aux opérations matricielles. Le coût des méthodes standards comme la SVD croît très rapidement avec la taille du problème N, c'est pourquoi des méthodes alternatives comme les approches aléatoires (i.e. basées sur la projection ou la sélection de colonnes/l'échantillonnage aléatoire) se popularisent. Ces méthodes sont en général peu coûteuses et facile à implémenter et optimiser, car elles ne mettent en oeuvre que des opérations matricielles simples comme des produits ou des orthogonalisations. Plus précisemment, les LRA aléatoires permettent de réduire le coût cubique en N des méthodes standards de factorisation au coût quadratique nécessaire à la réalisation de quelques produits matrices vecteurs, i.e., O(r × N^2) opérations où r est le rang numérique de la matrice. Dans un premier temps, nous présentons un algorithme efficace pour réaliser des MVPs en O(N) opérations, que nous appelons Uniform FMM (ufmm). Il est basé sur la combinaison d'une représentation hiérachique d'une matrice noyau et l'interpolation polynomiale du noyau associé sur une grille régulière (uniforme). Cette dernière propriété permet une accélération par FFT réduisant ainsi le temps de calcul et la consommation mémoire mais a des répercussions sur la précision et la stabilité de l'algorithme. Ensuite, la ufmm est utilisée pour accélerer les MVPs intervenants dans la SVD aléatoire (i.e., SVD par projection aléatoire) diminuant son coût asymptotique à O(r^2 × N). La méthode est particulièrement compétitive pour des distributions de points hétérogènes. Enfin, nous utilisons cet algorithme pour générer des champs de variables Gaussiens de manière efficace.

International audience

International audience; In this talk we will describe how H-matrix data sparse techniques can be implemented in a parallel hybrid sparse linear solver based on algebraic non overlapping domain decomposition approach.Strong-hierarchical matrix arithmetic and various clustering techniques to approximate the local Schur complements will be investigated, aiming at reducing workload and memory consumption while complying with structures of the local interfaces of the sub-domains. Utilization of these techniques to form effective global preconditioner will be presented.

International audience; ScalFMM (Parallel Fast Multipole Library for Large Scale Simulations) offers all the functionalities needed to perform large parallel simulations while enabling an easy customization of the simulation components: kernels, particles and cells. We will present how we use our library on two kinds of application involving boundary integral representations of physical fields. The first one implements isotropic dislocation kernels for Dislocation Dynamics and the second a time dependent kernel for acoustic problems.

Recently, graphics processors (GPUs) have been increasingly leveraged in a variety of scientific computing applications. However, architectural differences between CPUs and GPUs necessitate the development of algorithms that take advantage of GPU hardware. As sparse matrix vector multiplication (SPMV) operations are commonly used in finite element analysis, a new SPMV algorithm and several variations are developed for unstructured finite element meshes on GPUs. The effective bandwidth of current GPU algorithms and the newly proposed algorithms are measured and analyzed for 15 sparse matrices of varying sizes and varying sparsity structures. The effects of optimization and differences between the new GPU algorithm and its variants are then subsequently studied. Lastly, both new and current SPMV GPU algorithms are utilized in the GPU CG Solver in GPU finite element simulations of the heart. These results are then compared against parallel PETSc finite element implementation results. The effective bandwidth tests indicate that the new algorithms compare very favorably with current algorithms for a wide variety of sparse matrices and can yield very notable benefits. GPU finite element simulation results demonstrate the benefit of using GPUs for finite element analysis, and also show that the proposed algorithms can yield speedup factors up to 12-fold for real finite element applications., Comment: 35 pages, 22 figures

Kernel methods are widespread in machine learning; however, they are limited by the quadratic complexity of the construction, application, and storage of kernel matrices. Low-rank matrix approximation algorithms are widely used to address this problem and reduce the arithmetic and storage cost. However, we observed that for some datasets with wide intra-class variability, the optimal kernel parameter for smaller classes yields a matrix that is less well approximated by low-rank methods. In this paper, we propose an efficient structured low-rank approximation method -- the Block Basis Factorization (BBF) -- and its fast construction algorithm to approximate radial basis function (RBF) kernel matrices. Our approach has linear memory cost and floating-point operations for many machine learning kernels. BBF works for a wide range of kernel bandwidth parameters and extends the domain of applicability of low-rank approximation methods significantly. Our empirical results demonstrate the stability and superiority over the state-of-art kernel approximation algorithms., Comment: 16 pages, 5 figures

This article introduces a new fast direct solver for linear systems arising out of wide range of applications, integral equations, multivariate statistics, radial basis interpolation, etc., to name a few. \emph{The highlight of this new fast direct solver is that the solver scales linearly in the number of unknowns in all dimensions.} The solver, termed as Inverse Fast Multipole Method (abbreviated as IFMM), works on the same data-structure as the Fast Multipole Method (abbreviated as FMM). More generally, the solver can be immediately extended to the class of hierarchical matrices, denoted as $\mathcal{H}^2$ matrices with strong admissibility criteria (weak low-rank structure), i.e., \emph{the interaction between neighboring cluster of particles is full-rank whereas the interaction between particles corresponding to well-separated clusters can be efficiently represented as a low-rank matrix}. The algorithm departs from existing approaches in the fact that throughout the algorithm the interaction corresponding to neighboring clusters are always treated as full-rank interactions. Our approach relies on two major ideas: (i) The $N \times N$ matrix arising out of FMM (from now on termed as FMM matrix) can be represented as an extended sparser matrix of size $M \times M$, where $M \approx 3N$. (ii) While solving the larger extended sparser matrix, \emph{the fill-in's that arise in the matrix blocks corresponding to well-separated clusters are hierarchically compressed}. The ordering of the equations and the unknowns in the extended sparser matrix is strongly related to the local and multipole coefficients in the FMM~\cite{greengard1987fast} and \emph{the order of elimination is different from the usual nested dissection approach}. Numerical benchmarks on $2$D manifold confirm the linear scaling of the algorithm., 25 pages, 28 figures

High performance \FMM is crucial for the numerical simulation of many physical problems. In a previous study~\cite{Agullo2013}, we have shown that task-based \FMM provides the flexibility required to process a wide spectrum of particle distributions efficiently on multicore architectures. In this paper, we now show how such an approach can be extended to fully exploit heterogeneous platforms. For that, we design highly tuned GPU versions of the two dominant operators (P2P and M2L) as well as a scheduling strategy that dynamically decides which proportion of subsequent tasks are processed on regular CPU cores and on GPU accelerators. We assess our method with the StarPU runtime system for executing the resulting task flow on an Intel X5650 Nehalem multicore processor possibly enhanced with one, two or three Nvidia Fermi M2070 or M2090 GPUs. A detailed experimental study on two 30 million particle distributions (a cube and an ellipsoid) shows that the resulting software consistently achieves high performance across architectures.; Développer une méthode des Multipôles Rapide (FMM) à haute performance est cruciale pour des simulations numériques dans beaucoup de problèmes physiques. Dans une étude précédente~\cite{Agullo2013}, nous avons montré que l'utilisation d'un paradigme à base de tâches fournit la flexibilité nécessaire pour traiter efficacement un large spectre de distributions de particules sur des architectures homogènes. Dans ce document, nous montrons maintenant comment une telle approche peut être étendue pour exploiter toutes les unités de calculs (CPU et GPU) des machines hétérogènes. Pour cela, nous présentons une version optimisée pour GPU des deux opérateurs dominants (P2P et M2L) de la FMM ainsi qu'une stratégie d'ordonnancement qui décide dynamiquement quelle proportion de tâches est traitée par les cœurs CPU et par des accélérateurs GPU. Nous évaluons notre méthode avec le moteur d'exécution StarPU pour exécuter le flot de tâches résultant sur un processeur Intel X5650 Nehalem augmenté avec un, deux ou trois Nvidia Fermi M2070 ou M2090. Une étude expérimentale détaillée sur deux distributions de 30 millions de particules (un cube et un ellipsoïde) montre que nous obtenons des résultats performants sur cette architecture.

International audience; Among all the steps involved in DD simulations, the computation of the internal elastic forces and energy are by far the most resources consuming. However, since these are long-ranged and fast decreasing interactions, hierarchical algorithms like the Fast Multipole Method (FMM) are well suited for their fast evaluation. The relatively low accuracy required for the interaction forces between dislocations brought us to develop a more efficient approximation method for the farfield. On the other hand, the nearfield interactions are still evaluated analytically, which required a rather performant implementation (AVX, GPU...). Regarding parallelism, our code benefits from a hybrid OpenMP/MPI paradigm and a cache-conscious data structure. Finally, an accurate handling of topological elements intersecting the structure of the octree was considered. The latter feature implied careful modifications of the P2M/L2P operators in order to deal with shared memory model of parallelism.

Preliminary version of a paper to appear in SIAM SISC; Fast Multipole Methods (FMM) are a fundamental operation for the simulation of many physical problems. The high performance design of such methods usually requires to carefully tune the algorithm for both the targeted physics and hardware. In this paper, we propose a new approach that achieves high performance across architectures. Our method consists of expressing the Fast Multipole Method (FMM) algorithm as a task flow and employing a state-of-the-art runtime system, StarPU, to process the tasks on the different computing units. We carefully design the task flow, the mathematical operators, their implementations as well as scheduling schemes. Potentials and forces on 200 million particles are computed in 48.7 seconds on a homogeneous 160 cores SGI Altix UV 100 and good scalability is shown.; Les méthodes multipôles rapides (FMM) constituent une opération fondamentale pour la simulation de nombreux problèmes physiques. Leur mise en œuvre haute performance requiert habituellement d'optimiser attentivement l'algorithme à la fois pour la physique visée et le matériel utilisé. Dans ce papier, nous proposons une nouvelle approche qui atteint une performance élevée et portable. Notre méthode consiste à exprimer l'algorithme FMM comme un flot de tâches et d'employer un moteur d'exécution, StarPU, afin de traiter les tâches sur les différentes unités d'exécution. Nous concevons précisément le flot de tâches, les opérateurs mathématiques, leurs implémentations sur unité centrale de traitement (CPU) ainsi que les schémas d'ordonnancement. Nous calculons les potentiels et forces pour des problèmes avec 200 millions de particules en 48.7 secondes sur une machine homogène SGI Altix UV 100 comportant 160 cœurs

The method and the advantages of an evolutionary computing based approach using a steady state genetic algorithm (GA) for the parameterization of interatomic potentials for metal oxides within the shell model framework are developed and described. We show that the GA based methodology for the parameterization of interatomic force field functions is capable of (a) simultaneous optimization of the multiple phases or properties of a material in a single run, (b) facilitates the incremental re-optimization of the whole system as more data is available for either additional phases or material properties not included in previous runs, and (c) successful global optimization in the presence of multiple local minima in the parameter space. As an example, we apply the method towards simultaneous optimization of four distinct crystalline phases of Barium Titanate (BaTiO3 or BTO) using an ab initio density functional theory (DFT) based reference dataset. We find that the optimized force field function is capable of the prediction of the two phases not used in the optimization procedure, and that many derived physical properties such as the equilibrium lattice constants, unit cell volume, elastic properties, coefficient of thermal expansion, and average electronic polarization are in good agreement with the experimental results available from the literature., Comment: 31 pages, 16 figures

Computing reaction rates in biomolecular systems is a common goal of molecular dynamics simulations. The reactions considered often involve conformational changes in the molecule, either changes in the structure of a protein or the relative position of two molecules, for example when modeling the binding of a protein and ligand. Here we will consider the general problem of computing the rate of transfer from a subset A of the conformational space Omega to a subset B of Omega. It is assumed that A and B are associated with minimum energy basins and are long-lived states. Rates can be obtained using many different methods. We review some of the most popular approaches. We organize the different approaches roughly in chronological order and under four main categories: reactive flux, transition path sampling, conformation dynamics. The fourth class of methods, to which we do not give any specific name, in some sense attempts to combine features from transition path sampling and conformation dynamics. They include non-equilibrium umbrella sampling (Warmflash et al. [2007], Dickson et al. [2009b]), and weighted ensemble dynamics (Huber and Kim [1996]).

A fast multipole method (FMM) for asymptotically smooth kernel functions (1/r, 1/r^4, Gauss and Stokes kernels, radial basis functions, etc.) based on a Chebyshev interpolation scheme has been introduced in [Fong et al., 2009]. The method has been extended to oscillatory kernels (e.g., Helmholtz kernel) in [Messner et al., 2012]. Beside its generality this FMM turns out to be favorable due to its easy implementation and its high performance based on intensive use of highly optimized BLAS libraries. However, one of its bottlenecks is the precomputation of the multiple-to-local (M2L) operator, and its higher number of floating point operations (flops) compared to other FMM formulations. Here, we present several optimizations for that operator, which is known to be the costliest FMM operator. The most efficient ones do not only reduce the precomputation time by a factor up to 340 but they also speed up the matrix-vector product. We conclude with comparisons and numerical validations of all presented optimizations.

Fast Multipole Methods (FMM) are a fundamental operation for the simulation of many physical problems. The high performance design of such methods usually requires to carefully tune the algorithm for both the targeted physics and the hardware. In this paper, we propose a new approach that achieves high performance across architectures. Our method consists of expressing the FMM algorithm as a task flow and employing a state-of-the-art runtime system, StarPU, in order to process the tasks on the different processing units. We carefully design the task flow, the mathematical operators, their Central Processing Unit (CPU) and Graphics Processing Unit (GPU) implementations, as well as scheduling schemes. We compute potentials and forces of 200 million particles in 48.7 seconds on a homogeneous 160 cores SGI Altix UV 100 and of 38 million particles in 13.34 seconds on a heterogeneous 12 cores Intel Nehalem processor enhanced with 3 Nvidia M2090 Fermi GPUs., No. RR-7981 (2012)

The FIND algorithm is a fast algorithm designed to calculate certain entries of the inverse of a sparse matrix. Such calculation is critical in many applications, e.g., quantum transport in nano-devices. We extended the algorithm to other matrix inverse related calculations. Those are required for example to calculate the less-than Green's function and the current density through the device. For a 2D device discretized as an N_x x N_y mesh, the best known algorithms have a running time of O(N_x^3 N_y), whereas FIND only requires O(N_x^2 N_y). Even though this complexity has been reduced by an order of magnitude, the matrix inverse calculation is still the most time consuming part in the simulation of transport problems. We could not reduce the order of complexity, but we were able to significantly reduce the constant factor involved in the computation cost. By exploiting the sparsity and symmetry, the size of the problem beyond which FIND is faster than other methods typically decreases from a 130x130 2D mesh down to a 40x40 mesh. These improvements make the optimized FIND algorithm even more competitive for real-life applications.

Commercial graphics processors (GPUs) have high compute capacity at very low cost, which makes them attractive for general purpose scientific computing. In this paper we show how graphics processors can be used for N-body simulations to obtain improvements in performance over current generation CPUs. We have developed a highly optimized algorithm for performing the O(N^2) force calculations that constitute the major part of stellar and molecular dynamics simulations. In some of the calculations, we achieve sustained performance of nearly 100 GFlops on an ATI X1900XTX. The performance on GPUs is comparable to specialized processors such as GRAPE-6A and MDGRAPE-3, but at a fraction of the cost. Furthermore, the wide availability of GPUs has significant implications for cluster computing and distributed computing efforts like Folding@Home.