Eigenvalue computations in the context of data-sparse approximations of integral operators

In this work, we consider the numerical solution of a large eigenvalue problem resulting from a finite rank discretization of an integral operator. We are interested in computing a few eigenpairs, with an iterative method, so a matrix representation that allows for fast matrix-vector products is required. Hierarchical matrices are appropriate for this setting, and also provide cheap LU decompositions required in the spectral transformation technique. We illustrate the use of freely available software tools to address the problem, in particular SLEPc for the eigensolvers and HLib for the construction of H-matrices. The numerical tests are performed using an astrophysics application. Results show the benefits of the data-sparse representation compared to standard storage schemes, in terms of computational cost as well as memory requirements.This work was partially supported by the Spanish Ministerio de Ciencia e Innovacion under projects TIN2009-07519, TIN2012-32846 and AIC10-D-000600 and by Fundacao para a Ciencia e a Tecnologia - FCT under project FCT/MICINN proc 441.00.Román Moltó, JE.; Vasconcelos, PB.; Nunes, AL. (2013). Eigenvalue computations in the context of data-sparse approximations of integral operators. Journal of Computational and Applied Mathematics. 237(1):171-181. doi:10.1016/j.cam.2012.07.021S171181237

Approximation spectrale de matrices issues d opérateurs discrétisés

Cette thèse considère la solution numérique d'un problème aux valeurs propres de grandes dimensions, dans lequel l'opérateur est dérivé d'un problème de transfert radiatif. Ainsi, cette thèse étudie l'utilisation de matrices hiérarchiques, une représentation efficace de tableaux, très intéressante pour une utilisation avec des problèmes de grandes dimensions. Les matrices sont des représentations hiérarchiques de structures de données efficaces pour les matrices denses, l'idée de base étant la division d'une matrice en une hiérarchie de blocs et l approximation de certains blocs par une matrice de petite caractéristique. Son utilisation permet de diminuer la mémoire nécessaire tout en réduisant les coûts informatiques. L'application de l'utilisation de matrices hiérarchique est analysée dans le contexte de la solution numérique d'un problème aux valeurs propres de grandes dimensions résultant de la discrétisation d'un opérateur intégral. L'opérateur est de convolution et est défini par la première fonction exponentielle intégrale, donc faiblement singulière. Pour le calcul informatique, nous avons accès à HLIB (Hierarchical matrices LIBrary) qui fournit des routines pour la construction de la structure hiérarchique des matrices et des algorithmes pour les opérations approximative avec ces matrices. Nous incorporons certaines routines comme la multiplication matrice-vecteur ou la decomposition LU, en SLEPc (Hierarchical matrices LIBrary) pour explorer les algorithmes existants afin de résoudre les problèmes de valeur propre. Nous développons aussi des expressions analytiques pour l'approximation des noyaux dégénérés utilisés dans la thèse et déduire ainsi les limites supérieures d'erreur pour ces approximations. Les résultats numériques obtenus avec d'autres techniques pour résoudre le problème en question sont utilisés pour la comparaison avec ceux obtenus avec la nouvelle technique, illustrant l'efficacité de ce dernierIn this thesis, we consider the numerical solution of a large eigenvalue problem in which the integral operator comes from a radiative transfer problem. It is considered the use of hierarchical matrices, an efficient data-sparse representation of matrices, especially useful for large dimensional problems. It consists on low-rank subblocks leading to low memory requirements as well as cheap computational costs. We discuss the use of the hierarchical matrix technique in the numerical solution of a large scale eigenvalue problem arising from a finite rank discretization of an integral operator. The operator is of convolution type, it is defined through the first exponential-integral function and hence it is weakly singular. We access HLIB (Hierarchical matrices LIBrary) that provides, among others, routines for the construction of hierarchical matrix structures and arithmetic algorithms to perform approximative matrix operations. Moreover, it is incorporated the matrix-vector multiply routines from HLIB, as well as LU factorization for preconditioning, into SLEPc (Scalable Library for Eigenvalue Problem Computations) in order to exploit the available algorithms to solve eigenvalue problems. It is also developed analytical expressions for the approximate degenerate kernels and deducted error upper bounds for these approximations. The numerical results obtained with other approaches to solve the problem are used to compare with the ones obtained with this technique, illustrating the efficiency of the techniques developed and implemented in this workST ETIENNE-Bib. électronique (422189901) / SudocSudocFranceF

MPI-CUDA parallel linear solvers for block-tridiagonal matrices in the context of SLEPc's eigensolvers

[EN] We consider the computation of a few eigenpairs of a generalized eigenvalue problem Ax = lambda Bx with block-tridiagonal matrices, not necessarily symmetric, in the context of Krylov methods. In this kind of computation, it is often necessary to solve a linear system of equations in each iteration of the eigensolver, for instance when B is not the identity matrix or when computing interior eigenvalues with the shift-and-invert spectral transformation. In this work, we aim to compare different direct linear solvers that can exploit the block-tridiagonal structure. Block cyclic reduction and the Spike algorithm are considered. A parallel implementation based on MPI is developed in the context of the SLEPc library. The use of GPU devices to accelerate local computations shows to be competitive for large block sizes.This work was supported by Agencia Estatal de Investigacion (AEI) under grant TIN2016-75985-P, which includes European Commission ERDF funds. Alejandro Lamas Davina was supported by the Spanish Ministry of Education, Culture and Sport through a grant with reference FPU13-06655.Lamas Daviña, A.; Roman, JE. (2018). MPI-CUDA parallel linear solvers for block-tridiagonal matrices in the context of SLEPc's eigensolvers. Parallel Computing. 74:118-135. https://doi.org/10.1016/j.parco.2017.11.006S1181357

High Resolution Numerical Methods for Coupled Non-linear Multi-physics Simulations with Applications in Reactor Analysis

The modeling of nuclear reactors involves the solution of a multi-physics problem with widely varying time and length scales. This translates mathematically to solving a system of coupled, non-linear, and stiff partial differential equations (PDEs). Multi-physics applications possess the added complexity that most of the solution fields participate in various physics components, potentially yielding spatial and/or temporal coupling errors. This dissertation deals with the verification aspects associated with such a multi-physics code, i.e., the substantiation that the mathematical description of the multi-physics equations are solved correctly (both in time and space). Conventional paradigms used in reactor analysis problems employed to couple various physics components are often non-iterative and can be inconsistent in their treatment of the non-linear terms. This leads to the usage of smaller time steps to maintain stability and accuracy requirements, thereby increasing the overall computational time for simulation. The inconsistencies of these weakly coupled solution methods can be overcome using tighter coupling strategies and yield a better approximation to the coupled non-linear operator, by resolving the dominant spatial and temporal scales involved in the multi-physics simulation. A multi-physics framework, KARMA (K(c)ode for Analysis of Reactor and other Multi-physics Applications), is presented. KARMA uses tight coupling strategies for various physical models based on a Matrix-free Nonlinear-Krylov (MFNK) framework in order to attain high-order spatio-temporal accuracy for all solution fields in amenable wall clock times, for various test problems. The framework also utilizes traditional loosely coupled methods as lower-order solvers, which serve as efficient preconditioners for the tightly coupled solution. Since the software platform employs both lower and higher-order coupling strategies, it can easily be used to test and evaluate different coupling strategies and numerical methods and to compare their efficiency for problems of interest. Multi-physics code verification efforts pertaining to reactor applications are described and associated numerical results obtained using the developed multi-physics framework are provided. The versatility of numerical methods used here for coupled problems and feasibility of general non-linear solvers with appropriate physics-based preconditioners in the KARMA framework offer significantly efficient techniques to solve multi-physics problems in reactor analysis

Meshless Methods for the Neutron Transport Equation

Mesh-based methods for the numerical solution of partial differential equations (PDEs) require the division of the problem domain into non-overlapping, contiguous subdomains that conform to the problem geometry. The mesh constrains the placement and connectivity of the solution nodes over which the PDE is solved. In meshless methods, the solution nodes are independent of the problem geometry and do not require a mesh to determine connectivity. This allows the solution of PDEs on geometries that would be difficult to represent using even unstructured meshes. The ability to represent difficult geometries and place solution nodes independent of a mesh motivates the use of meshless methods for the neutron transport equation, which often includes spatially-dependent PDE coefficients and strong localized gradients. The meshless local Petrov-Galerkin (MLPG) method is applied to the steady-state and k-eigenvalue neutron transport equations, which are discretized in energy using the multigroup approximation and in angle using the discrete ordinates approximation. The MLPG method uses weighted residuals of the transport equation to solve for basis function expansion coefficients of the neutron angular flux. Connectivity of the solution nodes is determined by the shared support domain of overlapping meshless functions, such as radial basis functions (RBFs) and moving least squares (MLS) functions. To prevent oscillations in the neutron flux, the MLPG transport equation is stabilized by the streamline upwind Petrov-Galerkin (SUPG) method, which adds numerical diffusion to the streaming term. Global neutron conservation is enforced by using MLS basis and weight functions and appropriate SUPG parameters. The cross sections in the transport equation are approximated in accordance with global particle balance and without constraint on their spatial dependence or the location of the basis and weight functions. The equations for the strong-form meshless collocation approach are derived for comparison to the MLPG equations. Two integration schemes for the basis and weight functions in the MLPG method are presented, including a background mesh integration and a fully meshless integration approach. The method of manufactured solutions (MMS) is used to verify the resulting MLPG method in one, two and three dimensions. Results for realistic problems, including two-dimensional pincells, a reflected ellipsoid and a three-dimensional problem with voids, are verified by comparison to Monte Carlo simulations. Finally, meshless heat transfer equations are derived using a similar MLPG approach and verified using the MMS. These heat equation are coupled to the MLPG neutron transport equations, and results for a pincell are compared to values from a commercial pressurized water reactor.PHDNuclear Engineering & Radiological SciencesUniversity of Michigan, Horace H. Rackham School of Graduate Studieshttps://deepblue.lib.umich.edu/bitstream/2027.42/145796/1/brbass_1.pd

Numerical methods for computing Casimir interactions

We review several different approaches for computing Casimir forces and related fluctuation-induced interactions between bodies of arbitrary shapes and materials. The relationships between this problem and well known computational techniques from classical electromagnetism are emphasized. We also review the basic principles of standard computational methods, categorizing them according to three criteria---choice of problem, basis, and solution technique---that can be used to classify proposals for the Casimir problem as well. In this way, mature classical methods can be exploited to model Casimir physics, with a few important modifications.Comment: 46 pages, 142 references, 5 figures. To appear in upcoming Lecture Notes in Physics book on Casimir Physic

Riesz-projection-based methods for the numerical simulation of resonance phenomena in nanophotonics

Resonance effects are ubiquitous in physics and essential for understanding wave propagation and interference. In the field of nanophotonics, devices are often based on the strong confinement of light by resonances. The numerical simulation of resonances plays a crucial role for the design and optimization of the devices. The resonances are electromagnetic field solutions to the time-harmonic source-free Maxwell's equations with loss mechanisms. The corresponding eigenproblems are non-Hermitian due to the losses leading to complex-valued eigenvalues. The material dispersion, which is typically significant in nanophotonics, results in nonlinear eigenproblems. In this thesis, we develop an approach based on Riesz projections for the expansion of electromagnetic fields caused by light sources into resonances. The Riesz projection expansion is computed by contour integration in the complex frequency plane. The numerical realization essentially relies on solving Maxwell's equations with a source term, meaning solving linear systems of equations. For this, Maxwell's equations are directly evaluated at the given frequencies on the integration contours, which implies that linearization of the corresponding nonlinear eigenproblems is not required. This makes Riesz-projection-based approaches a natural choice for dealing with eigenproblems from the field of nanophotonics. We further extend the Riesz projection expansion approach to optical far-field quantities, which is not straightforward due to the spatial divergence of the resonances with increasing distance from the underlying resonators. Based on the ideas of the Riesz projection expansion, we introduce approaches for the calculation of physically relevant eigenvalues and for computing eigenvalue sensitivities. Physically relevant means that the eigenvalues are significant with respect to the resonance expansion of the physical observable of interest. By using physical solutions to Maxwell's equations for the contour integration, the developed numerical methods have a strong relation to physics. The methods can be applied to any material system and to any measurable physical quantity that can be derived from the electric field. We apply the numerical methods to several recent nanophotonic applications, for example, single-photon sources from the field of quantum technology, plasmonic nanostructures characterized by nonlocal material properties, and nanoantennas based on bound states in the continuum. The approaches introduced in this thesis are developed for nanophotonic systems, but can be applied to any resonance problem.Resonanzeffekte treten in allen physikalischen Systemen auf, die durch Wellen beschrieben werden, und sie sind für die Beschreibung von Wellenausbreitung und Interferenz unerlässlich. Auf dem Gebiet der Nanophotonik basieren viele Geräte auf den durch Lichtquellen angeregten Resonanzen mit ihren stark erhöhten elektromagnetischen Feldern. Die numerische Simulation von Resonanzen ist ein wichtiges Hilfsmittel für die Entwicklung und Optimierung der Geräte. Die Resonanzen sind die Lösungen der zeitharmonischen quellenfreien Maxwell-Gleichungen mit Verlustmechanismen. Die entsprechenden Eigenwertprobleme sind aufgrund der Verluste nicht-Hermitesch, was zu komplexwertigen Eigenwerten führt. Die Materialdispersion, die in der Nanophotonik typischerweise signifikant ist, führt zu nichtlinearen Eigenwertproblemen. In dieser Dissertation entwickeln wir einen auf der Riesz-Projektion basierenden Ansatz für die Expansion von elektromagnetischen Feldern, die von Lichtquellen erzeugt werden, in Resonanzen. Wir berechnen die Riesz-Projektionen durch Konturintegration in der komplexen Frequenzebene. Die numerische Realisierung basiert im Wesentlichen auf der Lösung der Maxwell-Gleichungen mit einem Quellterm, das heißt der Lösung von linearen Gleichungssystemen. Dabei werden die Maxwell-Gleichungen direkt bei den gegebenen Frequenzen auf den Integrationskonturen ausgewertet, sodass eine Linearisierung der entsprechenden nichtlinearen Eigenwertprobleme nicht erforderlich ist. Das macht die auf der Riesz-Projektion basierenden Methoden zu einer natürlichen Wahl für die Behandlung von Eigenwertproblemen aus dem Bereich der Nanophotonik. Wir erweitern den Ansatz der Riesz-Projektions-Expansion auf optische Größen im Fernfeld, was aufgrund der räumlichen Divergenz der Resonanzen mit zunehmender Entfernung von den zugrunde liegenden Resonatoren problematisch ist. Basierend auf den Ideen der Riesz-Projektions-Expansion entwickeln wir außerdem Methoden zur Berechnung physikalisch relevanter Eigenwerte und zur Berechnung von Sensitivitäten von Eigenwerten. Physikalisch relevant bedeutet, dass die Eigenwerte in Bezug auf die Resonanzexpansion der interessierenden physikalischen Größe signifikant sind. Durch die Verwendung physikalischer Lösungen der Maxwell-Gleichungen für die Konturintegration haben die entwickelten numerischen Methoden einen starken Bezug zur zugrunde liegenden Physik. Die Methoden können auf jedes Materialsystem und auf jede messbare physikalische Größe angewendet werden, die sich aus dem elektrischen Feld herleiten lässt. Wir wenden die numerischen Methoden auf mehrere aktuelle nanophotonische Strukturen an, wie zum Beispiel Einzelphotonenquellen aus dem Bereich der Quantentechnologie, plasmonische Nanostrukturen, die sich durch nichtlokale Materialeigenschaften auszeichnen, und Nanoantennen, die auf gebundenen Zuständen im Kontinuum basieren. Die in dieser Dissertation vorgestellten Ansätze werden für nanophotonische Systeme entwickelt, lassen sich aber auf jedes Resonanzproblem anwenden

Numerical Methods for Computing Casimir Interactions

We review several different approaches for computing Casimir forces and related fluctuation-induced interactions between bodies of arbitrary shapes and materials. The relationships between this problem and well known computational techniques from classical electromagnetism are emphasized. We also review the basic principles of standard computational methods, categorizing them according to three criteria\-choice of problem, basis, and solution technique\-that can be used to classify proposals for the Casimir problem as well. In this way, mature classical methods can be exploited to model Casimir physics, with a few important modifications.Keywords: Imaginary Frequency, Perfectly Match Layer, Casimir Force, Casimir Energy, Perfect Electric ConductorUnited States. Army Research Office (Contract W911NF-07-D-0004)Massachusetts Institute of Technology. Ferry FundUnited States. Defense Advanced Research Projects Agency (Contract N66001-09-1-2070-DOD

Software for Exascale Computing - SPPEXA 2016-2019

This open access book summarizes the research done and results obtained in the second funding phase of the Priority Program 1648 "Software for Exascale Computing" (SPPEXA) of the German Research Foundation (DFG) presented at the SPPEXA Symposium in Dresden during October 21-23, 2019. In that respect, it both represents a continuation of Vol. 113 in Springer’s series Lecture Notes in Computational Science and Engineering, the corresponding report of SPPEXA’s first funding phase, and provides an overview of SPPEXA’s contributions towards exascale computing in today's sumpercomputer technology. The individual chapters address one or more of the research directions (1) computational algorithms, (2) system software, (3) application software, (4) data management and exploration, (5) programming, and (6) software tools. The book has an interdisciplinary appeal: scholars from computational sub-fields in computer science, mathematics, physics, or engineering will find it of particular interest