Source author record

Leslie Greengard

Leslie Greengard appears in the imported research catalog. Authorship, coauthor and topic links are available while profile ownership is still unclaimed.

ResearcherUnclaimed source record

Catalog footprint

What is connected

32works
14topics
4close collaborators

Actions

Connect this record

Log in to claim

Research graph

See the researcher in context

Open full explorer

Inspect adjacent papers, topics, institutions and collaborators without losing the researcher page.

Building this map preview

BZPEER is loading the nearby papers, people, topics and institutions for this page.

Published work

32 published item(s)

preprint2023arXiv

FMM-LU: A fast direct solver for multiscale boundary integral equations in three dimensions

We present a fast direct solver for boundary integral equations on complex surfaces in three dimensions using an extension of the recently introduced recursive strong skeletonization scheme. For problems that are not highly oscillatory, our algorithm computes an ${LU}$-like hierarchical factorization of the dense system matrix, permitting application of the inverse in $\mathcal O(n)$ time, where $n$ is the number of unknowns on the surface. The factorization itself also scales linearly with the system size, albeit with a somewhat larger constant. The scheme is built on a level-restricted adaptive octree data structure, and therefore it is compatible with highly nonuniform discretizations. Furthermore, the scheme is coupled with high-order accurate locally-corrected Nyström quadrature methods to integrate the singular and weakly-singular Green's functions used in the integral representations. Our method has immediate applications to a variety of problems in computational physics. We concentrate here on studying its performance in acoustic scattering (governed by the Helmholtz equation) at low to moderate frequencies, and provide rigorous justification for compression of submatrices via proxy surfaces.

preprint2020arXiv

Geometry of the Phase Retrieval Problem

One of the most powerful approaches to imaging at the nanometer or subnanometer length scale is coherent diffraction imaging using X-ray sources. For amorphous (non-crystalline) samples, the raw data can be interpreted as the modulus of the continuous Fourier transform of the unknown object. Making use of prior information about the sample (such as its support), a natural goal is to recover the phase through computational means, after which the unknown object can be visualized at high resolution. While many algorithms have been proposed for this phase retrieval problem, careful analysis of its well-posedness has received relatively little attention. In this paper, we show that the problem is, in general, not well-posed and describe some of the underlying issues that are responsible for the ill-posedness. We then show how this analysis can be used to develop experimental protocols that lead to better conditioned inverse problems.

preprint2019arXiv

High-order discretization of a stable time-domain integral equation for 3D acoustic scattering

We develop a high-order, explicit method for acoustic scattering in three space dimensions based on a combined-field time-domain integral equation. The spatial discretization, of Nyström type, uses Gaussian quadrature on panels combined with a special treatment of the weakly singular kernels arising in near-neighbor interactions. In time, a new class of convolution splines is used in a predictor-corrector algorithm. Experiments on a torus and a perturbed torus are used to explore the stability and accuracy of the proposed scheme. This involved around one thousand solver runs, at up to 8th order and up to around 20,000 spatial unknowns, demonstrating 5-9 digits of accuracy. In addition we show that parameters in the combined field formulation, chosen on the basis of analysis for the sphere and other convex scatterers, work well in these cases.

preprint2016arXiv

Fast convolution with free-space Green's functions

We introduce a fast algorithm for computing volume potentials - that is, the convolution of a translation invariant, free-space Green's function with a compactly supported source distribution defined on a uniform grid. The algorithm relies on regularizing the Fourier transform of the Green's function by cutting off the interaction in physical space beyond the domain of interest. This permits the straightforward application of trapezoidal quadrature and the standard FFT, with superalgebraic convergence for smooth data. Moreover, the method can be interpreted as employing a Nystrom discretization of the corresponding integral operator, with matrix entries which can be obtained explicitly and rapidly. This is of use in the design of preconditioners or fast direct solvers for a variety of volume integral equations. The method proposed permits the computation of any derivative of the potential, at the cost of an additional FFT.

preprint2016arXiv

High resolution inverse scattering in two dimensions using recursive linearization

We describe a fast, stable algorithm for the solution of the inverse acoustic scattering problem in two dimensions. Given full aperture far field measurements of the scattered field for multiple angles of incidence, we use Chen's method of recursive linearization to reconstruct an unknown sound speed at resolutions of thousands of square wavelengths in a fully nonlinear regime. Despite the fact that the underlying optimization problem is formally ill-posed and non-convex, recursive linearization requires only the solution of a sequence of linear least squares problems at successively higher frequencies. By seeking a suitably band-limited approximation of the sound speed profile, each least squares calculation is well-conditioned and involves the solution of a large number of forward scattering problems, for which we employ a recently developed, spectrally accurate, fast direct solver. For the largest problems considered, involving 19,600 unknowns, approximately one million partial differential equations were solved, requiring approximately two days to compute using a parallel MATLAB implementation on a multi-core workstation.

preprint2015arXiv

A new hybrid integral representation for frequency domain scattering in layered media

A variety of problems in acoustic and electromagnetic scattering require the evaluation of impedance or layered media Green's functions. Given a point source located in an unbounded half-space or an infinitely extended layer, Sommerfeld and others showed that Fourier analysis combined with contour integration provides a systematic and broadly effective approach, leading to what is generally referred to as the Sommerfeld integral representation. When either the source or target is at some distance from an infinite boundary, the number of degrees of freedom needed to resolve the scattering response is very modest. When both are near an interface, however, the Sommerfeld integral involves a very large range of integration and its direct application becomes unwieldy. Historically, three schemes have been employed to overcome this difficulty: the method of images, contour deformation, and asymptotic methods of various kinds. None of these methods make use of classical layer potentials in physical space, despite their advantages in terms of adaptive resolution and high-order accuracy. The reason for this is simple: layer potentials are impractical in layered media or half-space geometries since they require the discretization of an infinite boundary. In this paper, we propose a hybrid method which combines layer potentials (physical-space) on a finite portion of the interface together with a Sommerfeld-type (Fourier) correction. We prove that our method is efficient and rapidly convergent for arbitrarily located sources and targets, and show that the scheme is particularly effective when solving scattering problems for objects which are close to the half-space boundary or even embedded across a layered media interface.

preprint2015arXiv

Fast Direct Methods for Gaussian Processes

A number of problems in probability and statistics can be addressed using the multivariate normal (Gaussian) distribution. In the one-dimensional case, computing the probability for a given mean and variance simply requires the evaluation of the corresponding Gaussian density. In the $n$-dimensional setting, however, it requires the inversion of an $n \times n$ covariance matrix, $C$, as well as the evaluation of its determinant, $\det(C)$. In many cases, such as regression using Gaussian processes, the covariance matrix is of the form $C = σ^2 I + K$, where $K$ is computed using a specified covariance kernel which depends on the data and additional parameters (hyperparameters). The matrix $C$ is typically dense, causing standard direct methods for inversion and determinant evaluation to require $\mathcal O(n^3)$ work. This cost is prohibitive for large-scale modeling. Here, we show that for the most commonly used covariance functions, the matrix $C$ can be hierarchically factored into a product of block low-rank updates of the identity matrix, yielding an $\mathcal O (n\log^2 n) $ algorithm for inversion. More importantly, we show that this factorization enables the evaluation of the determinant $\det(C)$, permitting the direct calculation of probabilities in high dimensions under fairly broad assumptions on the kernel defining $K$. Our fast algorithm brings many problems in marginalization and the adaptation of hyperparameters within practical reach using a single CPU core. The combination of nearly optimal scaling in terms of problem size with high-performance computing resources will permit the modeling of previously intractable problems. We illustrate the performance of the scheme on standard covariance kernels.

preprint2015arXiv

Fast, adaptive, high order accurate discretization of the Lippmann-Schwinger equation in two dimension

We present a fast direct solver for two dimensional scattering problems, where an incident wave impinges on a penetrable medium with compact support. We represent the scattered field using a volume potential whose kernel is the outgoing Green's function for the exterior domain. Inserting this representation into the governing partial differential equation, we obtain an integral equation of the Lippmann-Schwinger type. The principal contribution here is the development of an automatically adaptive, high-order accurate discretization based on a quad tree data structure which provides rapid access to arbitrary elements of the discretized system matrix. This permits the straightforward application of state-of-the-art algorithms for constructing compressed versions of the solution operator. These solvers typically require $O(N^{3/2})$ work, where $N$ denotes the number of degrees of freedom. We demonstrate the performance of the method for a variety of problems in both the low and high frequency regimes.

preprint2014arXiv

A fast multipole method for the evaluation of elastostatic fields in a half-space with zero normal stress

In this paper, we present a fast multipole method (FMM) for the half-space Green's function in a homogeneous elastic half-space subject to zero normal stress, for which an explicit solution was given by Mindlin (1936). The image structure of this Green's function is unbounded, so that standard outgoing representations are not easily available. We introduce two such representations here, one involving an expansion in plane waves and one involving a modified multipole expansion. Both play a role in the FMM implementation.

preprint2014arXiv

A fast semi-direct least squares algorithm for hierarchically block separable matrices

We present a fast algorithm for linear least squares problems governed by hierarchically block separable (HBS) matrices. Such matrices are generally dense but data-sparse and can describe many important operators including those derived from asymptotically smooth radial kernels that are not too oscillatory. The algorithm is based on a recursive skeletonization procedure that exposes this sparsity and solves the dense least squares problem as a larger, equality-constrained, sparse one. It relies on a sparse QR factorization coupled with iterative weighted least squares methods. In essence, our scheme consists of a direct component, comprised of matrix compression and factorization, followed by an iterative component to enforce certain equality constraints. At most two iterations are typically required for problems that are not too ill-conditioned. For an $M \times N$ HBS matrix with $M \geq N$ having bounded off-diagonal block rank, the algorithm has optimal $\mathcal{O} (M + N)$ complexity. If the rank increases with the spatial dimension as is common for operators that are singular at the origin, then this becomes $\mathcal{O} (M + N)$ in 1D, $\mathcal{O} (M + N^{3/2})$ in 2D, and $\mathcal{O} (M + N^{2})$ in 3D. We illustrate the performance of the method on both over- and underdetermined systems in a variety of settings, with an emphasis on radial basis function approximation and efficient updating and downdating.

preprint2014arXiv

A fast solver for multi-particle scattering in a layered medium

In this paper, we consider acoustic or electromagnetic scattering in two dimensions from an infinite three-layer medium with thousands of wavelength-size dielectric particles embedded in the middle layer. Such geometries are typical of microstructured composite materials, and the evaluation of the scattered field requires a suitable fast solver for either a single configuration or for a sequence of configurations as part of a design or optimization process. We have developed an algorithm for problems of this type by combining the Sommerfeld integral representation, high order integral equation discretization, the fast multipole method and classical multiple scattering theory. The efficiency of the solver is illustrated with several numerical experiments.

preprint2014arXiv

Fast and Accurate Evaluation of Nonlocal Coulomb and Dipole-Dipole Interactions via the Nonuniform FFT

We present a fast and accurate algorithm for the evaluation of nonlocal (long-range) Coulomb and dipole-dipole interactions in free space. The governing potential is simply the convolution of an interaction kernel $U(\bx)$ and a density function $ρ(\bx)=|ψ(\bx)|^2$, for some complex-valued wave function $ψ(\bx)$, permitting the formal use of Fourier methods. These are hampered by the fact that the Fourier transform of the interaction kernel $\widehat{U}(\bk)$ has a singularity at the origin $\bk={\bf 0}$ in Fourier (phase) space. Thus, accuracy is lost when using a uniform Cartesian grid in $\bk$ which would otherwise permit the use of the FFT for evaluating the convolution. Here, we make use of a high-order discretization of the Fourier integral, accelerated by the nonuniform fast Fourier transform (NUFFT). By adopting spherical and polar phase-space discretizations in three and two dimensions, respectively, the singularity in $\hat{U}(\bk)$ at the origin is canceled, so that only a modest number of degrees of freedom are required to evaluate the Fourier integral, assuming that the density function $ρ(\bx)$ is smooth and decays sufficiently fast as $\bx \rightarrow \infty$. More precisely, the calculation requires $O(N\log N)$ operations, where $N$ is the total number of discretization points in the computational domain. Numerical examples are presented to demonstrate the performance of the algorithm.

preprint2014arXiv

Inverse Obstacle scattering in two dimensions with multiple frequency data and multiple angles of incidence

We consider the problem of reconstructing the shape of an impenetrable sound-soft obstacle from scattering measurements. The input data is assumed to be the far-field pattern generated when a plane wave impinges on an unknown obstacle from one or more directions and at one or more frequencies. It is well known that this inverse scattering problem is both ill posed and nonlinear. It is common practice to overcome the ill posedness through the use of a penalty method or Tikhonov regularization. Here, we present a more physical regularization, based simply on restricting the unknown boundary to be band-limited in a suitable sense. To overcome the nonlinearity of the problem, we use a variant of Newton's method. When multiple frequency data is available, we supplement Newton's method with the recursive linearization approach due to Chen. During the course of solving the inverse problem, we need to compute the solution to a large number of forward scattering problems. For this, we use high-order accurate integral equation discretizations, coupled with fast direct solvers when the problem is sufficiently large.

preprint2014arXiv

Randomized methods for rank-deficient linear systems

We present a simple, accurate method for solving consistent, rank-deficient linear systems, with or without addi- tional rank-completing constraints. Such problems arise in a variety of applications, such as the computation of the eigenvectors of a matrix corresponding to a known eigenvalue. The method is based on elementary linear algebra combined with the ob- servation that if the matrix is rank-k deficient, then a random rank-k perturbation yields a nonsingular matrix with probability 1.

preprint2014arXiv

The Decoupled Potential Integral Equation for Time-Harmonic Electromagnetic Scattering

We present a new formulation for the problem of electromagnetic scattering from perfect electric conductors. While our representation for the electric and magnetic fields is based on the standard vector and scalar potentials ${\bf A},ϕ$ in the Lorenz gauge, we establish boundary conditions on the potentials themselves, rather than on the field quantities. This permits the development of a well-conditioned second kind Fredholm integral equation which has no spurious resonances, avoids low frequency breakdown, and is insensitive to the genus of the scatterer. The equations for the vector and scalar potentials are decoupled. That is, the unknown scalar potential defining the scattered field, $ϕ^{Sc}$, is determined entirely by the incident scalar potential $ϕ^{In}$. Likewise, the unknown vector potential defining the scattered field, ${\bf A}^{Sc}$, is determined entirely by the incident vector potential ${\bf A}^{In}$. This decoupled formulation is valid not only in the static limit but for arbitrary $ω\ge 0$.

preprint2014arXiv

Visualizing Skin Effects in Conductors with MRI: ${}^7$Li MRI Experiments and Calculations

While experiments on metals have been performed since the early days of NMR (and DNP), the use of bulk metal is normally avoided. Instead, often powders have been used in combination with low fields, so that skin depth effects could be neglected. Another complicating factor of acquiring NMR spectra or MRI images of bulk metal is the strong signal dependence on the orientation between the sample and the radio frequency (RF) coil, leading to non-intuitive image distortions and inaccurate quantification. Such factors are particularly important for NMR and MRI of batteries and other electrochemical devices. Here, we show results from a systematic study combining RF field calculations with experimental MRI of $^7$Li metal to visualize skin depth effects directly and to analyze the RF field orientation effect on MRI of bulk metal. It is shown that a certain degree of selectivity can be achieved for particular faces of the metal, simply based on the orientation of the sample. By combining RF field calculations with bulk magnetic susceptibility calculations accurate NMR spectra can be obtained from first principles. Such analyses will become valuable in many applications involving battery systems, but also metals, in general.

preprint2013arXiv

A fast direct solver for scattering from periodic structures with multiple material interfaces in two dimensions

We present a new integral equation method for the calculation of two-dimensional scattering from periodic structures involving triple-points (multiple materials meeting at a single point). The combination of a robust and high-order accurate integral representation and a fast direct solver permits the efficient simulation of scattering from fixed structures at multiple angles of incidence. We demonstrate the performance of the scheme with several numerical examples.

preprint2013arXiv

Debye Sources, Beltrami Fields, and a Complex Structure on Maxwell Fields

The Debye source representation for solutions to the time harmonic Maxwell equations is extended to bounded domains with finitely many smooth boundary components. A strong uniqueness result is proved for this representation. Natural complex structures are identified on the vector spaces of time-harmonic Maxwell fields. It is shown that in terms of Debye source data, these complex structures are uniformized, that is, represented by a fixed linear map on a fixed vector space, independent of the frequency. This complex structure relates time-harmonic Maxwell fields to constant-k Beltrami fields, i.e. solutions of the equation curl(E) = kE. A family of self-adjoint boundary conditions are defined for the Beltrami operator. This leads to a proof of the existence of zero-flux, constant-k, force-free Beltrami fields for any bounded region in R^3, as well as a constructive method to find them. The family of self-adjoint boundary value problems defines a new spectral invariant for bounded domains in R^3.

preprint2013arXiv

Efficient sum-of-exponentials approximations for the heat kernel and their applications

In this paper, we show that efficient separated sum-of-exponentials approximations can be constructed for the heat kernel in any dimension. In one space dimension, the heat kernel admits an approximation involving a number of terms that is of the order $O(\log(\frac{T}δ) (\log(\frac{1}ε)+\log\log(\frac{T}δ)))$ for any $x\in\bbR$ and $δ\leq t \leq T$, where $ε$ is the desired precision. In all higher dimensions, the corresponding heat kernel admits an approximation involving only $O(\log^2(\frac{T}δ))$ terms for fixed accuracy $ε$. These approximations can be used to accelerate integral equation-based methods for boundary value problems governed by the heat equation in complex geometry. The resulting algorithms are nearly optimal. For $N_S$ points in the spatial discretization and $N_T$ time steps, the cost is $O(N_S N_T \log^2 \frac{T}δ)$ in terms of both memory and CPU time for fixed accuracy $ε$. The algorithms can be parallelized in a straightforward manner. Several numerical examples are presented to illustrate the accuracy and stability of these approximations.

preprint2013arXiv

Norm-preserving discretization of integral equations for elliptic PDEs with internal layers I: the one-dimensional case

We investigate the behavior of integral formulations of variable coefficient elliptic partial differential equations (PDEs) in the presence of steep internal layers. In one dimension, the equations that arise can be solved analytically and the condition numbers estimated in various L_p norms. We show that high-order accurate Nyström discretization leads to well-conditioned finite dimensional linear systems if and only if the discretization is both norm-preserving in a correctly chosen L_p space and adaptively refined in the internal layer.

preprint2013arXiv

On the convergence of local expansions of layer potentials

In a recently developed quadrature method (quadrature by expansion or QBX), it was demonstrated that weakly singular or singular layer potentials can be evaluated rapidly and accurately on surface by making use of local expansions about carefully chosen off-surface points. In this paper, we derive estimates for the rate of convergence of these local expansions, providing the analytic foundation for the QBX method. The estimates may also be of mathematical interest, particularly for microlocal or asymptotic analysis in potential theory.

preprint2013arXiv

Quadrature by Expansion: A New Method for the Evaluation of Layer Potentials

Integral equation methods for the solution of partial differential equations, when coupled with suitable fast algorithms, yield geometrically flexible, asymptotically optimal and well-conditioned schemes in either interior or exterior domains. The practical application of these methods, however, requires the accurate evaluation of boundary integrals with singular, weakly singular or nearly singular kernels. Historically, these issues have been handled either by low-order product integration rules (computed semi-analytically), by singularity subtraction/cancellation, by kernel regularization and asymptotic analysis, or by the construction of special purpose "generalized Gaussian quadrature" rules. In this paper, we present a systematic, high-order approach that works for any singularity (including hypersingular kernels), based only on the assumption that the field induced by the integral operator is locally smooth when restricted to either the interior or the exterior. Discontinuities in the field across the boundary are permitted. The scheme, denoted QBX (quadrature by expansion), is easy to implement and compatible with fast hierarchical algorithms such as the fast multipole method. We include accuracy tests for a variety of integral operators in two dimensions on smooth and corner domains.

preprint2013arXiv

The solution of the scalar wave equation in the exterior of a sphere

We derive new, explicit representations for the solution to the scalar wave equation in the exterior of a sphere, subject to either Dirichlet or Robin boundary conditions. Our formula leads to a stable and high-order numerical scheme that permits the evaluation of the solution at an arbitrary target, without the use of a spatial grid and without numerical dispersion error. In the process, we correct some errors in the analytic literature concerning the asymptotic behavior of the logarithmic derivative of the spherical modified Hankel function. We illustrate the performance of the method with several numerical examples.

preprint2012arXiv

A consistency condition for the vector potential in multiply-connected domains

A classical problem in electromagnetics concerns the representation of the electric and magnetic fields in the low-frequency or static regime, where topology plays a fundamental role. For multiply connected conductors, at zero frequency the standard boundary conditions on the tangential components of the magnetic field do not uniquely determine the vector potential. We describe a (gauge-invariant) consistency condition that overcomes this non-uniqueness and resolves a longstanding difficulty in inverting the magnetic field integral equation.

preprint2012arXiv

A fast direct solver for structured linear systems by recursive skeletonization

We present a fast direct solver for structured linear systems based on multilevel matrix compression. Using the recently developed interpolative decomposition of a low-rank matrix in a recursive manner, we embed an approximation of the original matrix into a larger, but highly structured sparse one that allows fast factorization and application of the inverse. The algorithm extends the Martinsson/Rokhlin method developed for 2D boundary integral equations and proceeds in two phases: a precomputation phase, consisting of matrix compression and factorization, followed by a solution phase to apply the matrix inverse. For boundary integral equations which are not too oscillatory, e.g., based on the Green's functions for the Laplace or low-frequency Helmholtz equations, both phases typically have complexity O(N) in two dimensions, where $N$ is the number of discretization points. In our current implementation, the corresponding costs in three dimensions are $O(N^{3/2})$ and $O(N \log N)$ for precomputation and solution, respectively. Extensive numerical experiments show a speedup of $\sim 100$ for the solution phase over modern fast multipole methods; however, the cost of precomputation remains high. Thus, the solver is particularly suited to problems where large numbers of iterations would be required. Such is the case with ill-conditioned linear systems or when the same system is to be solved with multiple right-hand sides. Our algorithm is implemented in Fortran and freely available.

preprint2012arXiv

A fast, high-order solver for the Grad-Shafranov equation

We present a new fast solver to calculate fixed-boundary plasma equilibria in toroidally axisymmetric geometries. By combining conformal mapping with Fourier and integral equation methods on the unit disk, we show that high-order accuracy can be achieved for the solution of the equilibrium equation and its first and second derivatives. Smooth arbitrary plasma cross-sections as well as arbitrary pressure and poloidal current profiles are used as initial data for the solver. Equilibria with large Shafranov shifts can be computed without difficulty. Spectral convergence is demonstrated by comparing the numerical solution with a known exact analytic solution. A fusion-relevant example of an equilibrium with a pressure pedestal is also presented.

preprint2012arXiv

On the efficient representation of the half-space impedance Green's function for the Helmholtz equation

A classical problem in acoustic (and electromagnetic) scattering concerns the evaluation of the Green's function for the Helmholtz equation subject to impedance boundary conditions on a half-space. The two principal approaches used for representing this Green's function are the Sommerfeld integral and the (closely related) method of complex images. The former is extremely efficient when the source is at some distance from the half-space boundary, but involves an unwieldy range of integration as the source gets closer and closer. Complex image-based methods, on the other hand, can be quite efficient when the source is close to the boundary, but they do not easily permit the use of the superposition principle since the selection of complex image locations depends on both the source and the target. We have developed a new, hybrid representation which uses a finite number of real images (dependent only on the source location) coupled with a rapidly converging Sommerfeld-like integral. While our method applies in both two and three dimensions, we restrict the detailed analysis and numerical experiments here to the two-dimensional case.

preprint2011arXiv

Debye Sources and the Numerical Solution of the Time Harmonic Maxwell Equations, II

In this paper, we develop a new integral representation for the solution of the time harmonic Maxwell equations in media with piecewise constant dielectric permittivity and magnetic permeability in R^3. This representation leads to a coupled system of Fredholm integral equations of the second kind for four scalar densities supported on the material interface. Like the classical Muller equation, it has no spurious resonances. Unlike the classical approach, however, the representation does not suffer from low frequency breakdown. We illustrate the performance of the method with numerical examples.

preprint2011arXiv

Fast elliptic solvers in cylindrical coordinates and the Coulomb collision operator

In this paper, we describe a new class of fast solvers for separable elliptic partial differential equations in cylindrical coordinates $(r,θ,z)$ with free-space radiation conditions. By combining integral equation methods in the radial variable $r$ with Fourier methods in $θ$ and $z$, we show that high-order accuracy can be achieved in both the governing potential and its derivatives. A weak singularity arises in the Fourier transform with respect to $z$ that is handled with special purpose quadratures. We show how these solvers can be applied to the evaluation of the Coulomb collision operator in kinetic models of ionized gases.

preprint2011arXiv

Fast multi-particle scattering: a hybrid solver for the Maxwell equations in microstructured materials

A variety of problems in device and materials design require the rapid forward modeling of Maxwell's equations in complex micro-structured materials. By combining high-order accurate integral equation methods with classical multiple scattering theory, we have created an effective simulation tool for materials consisting of an isotropic background in which are dispersed a large number of micro- or nano-scale metallic or dielectric inclusions.

preprint2010arXiv

A new integral representation for quasiperiodic fields and its application to two-dimensional band structure calculations

In this paper, we consider band-structure calculations governed by the Helmholtz or Maxwell equations in piecewise homogeneous periodic materials. Methods based on boundary integral equations are natural in this context, since they discretize the interface alone and can achieve high order accuracy in complicated geometries. In order to handle the quasi-periodic conditions which are imposed on the unit cell, the free-space Green's function is typically replaced by its quasi-periodic cousin. Unfortunately, the quasi-periodic Green's function diverges for families of parameter values that correspond to resonances of the empty unit cell. Here, we bypass this problem by means of a new integral representation that relies on the free-space Green's function alone, adding auxiliary layer potentials on the boundary of the unit cell itself. An important aspect of our method is that by carefully including a few neighboring images, the densities may be kept smooth and convergence rapid. This framework results in an integral equation of the second kind, avoids spurious resonances, and achieves spectral accuracy. Because of our image structure, inclusions which intersect the unit cell walls may be handled easily and automatically. Our approach is compatible with fast-multipole acceleration, generalizes easily to three dimensions, and avoids the complication of divergent lattice sums.