Mathematical and Computational Materials Science

We study the Representative Volume Element (RVE) method, which is a method to approximately infer the effective behavior $a_{\rm hom}$ of a stationary random medium, described by a coefficient field $a(x)$ and the corresponding linear elliptic operator $-\nabla\cdot a\nabla$. In line with the theory of homogenization, the method proceeds by computing $d=3$ ( denoting the space dimension) correctors, however on a “representative” volume element, i.e. box with, say, periodic boundary conditions. The main message of this article is: Periodize the ensemble instead of its realizations.

By this we mean that it is better to sample from a suitably periodized ensemble than to periodically extend the restriction $a(x)$ of a realization from the whole-space ensemble $\langle\cdot\rangle$. We make this point by investigating the bias (or systematic error), i.e. the difference between $a_{\rm hom}$ and the expected value of the RVE method, in terms of its scaling w.r.t. the lateral size $L$ of the box. In case of periodizing $a(x)$, we heuristically argue that this error is generically $O(L^{-1})$. In case of a suitable periodization of $\langle\cdot\rangle$, we rigorously show that it is $O(L^{-d})$. In fact, we give a characterization of the leading-order error term for both strategies.

We carry out the rigorous analysis in the convenient setting of ensembles $\langle\cdot\rangle$ of Gaussian type, which allow for a straightforward periodization, passing via the (integrable) covariance function. This setting has also the advantage of making Malliavin calculus available for optimal stochastic estimates of correctors. We actually need control of second-order correctors to capture the leading-order error term in the presence of cancellations due to point symmetry. This is joint work with Nicolas Clozeau, Marc Josien, and Qiang Xu.

Accurate molecular simulation requires computationally expensive quantum chemistry models that makes simulating complex material phenomena or large molecules intractable. However, if no chemistry is required, but only interatomic forces then it should in principle be possible to construct much cheaper surrogate models, interatomic potentials, that capture full QM accuracy. This talk will review recent attempts, with focus on my personal analysis and numerical analysis perspectives, to achieve this. Specifically, I will explore whether we can rigorously justify the extremely low-dimensional functional forms proposed for interatomic potentials, and whether we can construct practical approximation schemes that can, in principle at least, close the complexity gap between density functional theory and interatomic potentials.

The talk will describe some recent work on the application of accelerated computing and machine learning for multi-scale modeling of materials. We will demonstrate the ideas with illustrative applications and discuss open issues.

This work is concerned with the topology optimization of so-called lattice materials, i.e., porous structures made of periodically perforated material, where the microscopic periodic cell can be macroscopically modulated and oriented. Lattice materials are becoming increasingly popular since they can be built by additive manufacturing techniques. The main idea is to optimize the homogenized formulation of this problem, which is an easy task of parametric optimization, then to project the optimal microstructure at a desired length-scale, which is a delicate issue, albeit computationally cheap. The main novelty of our work is, in a plane setting, the conformal treatment of the optimal orientation of the microstructure. In other words, although the periodicity cell has varying parameters and orientation throughout the computational domain, the angles between its members or bars are conserved. Several numerical examples are presented for compliance minimization in 2-d. Extension to the 3-d case will also be discussed.

Electronic structure calculations, especially those using density functional theory (DFT), have been very useful in understanding and predicting a wide range of materials properties. The importance of DFT calculations to engineering and physical sciences is evident from the fact that ~20% of computational resources on some of the world’s largest public supercomputers are devoted to DFT calculations. Despite the wide adoption of DFT, and the tremendous progress in theory and numerical methods over the decades, the following challenges remain. Firstly, the state-of-the-art implementations of DFT suffer from cell-size and geometry limitations, with the widely used codes in solid state physics being limited to periodic geometries and typical simulation domains containing a few hundred atoms. This limits the complexity of materials systems that can be treated using DFT calculations. Secondly, there are many materials systems (such as strongly-correlated systems) where the widely used model exchange-correlation functionals, which account for the many-body quantum mechanical interactions between electrons, are inaccurate. Addressing these challenges will enable large-scale quantum-accuracy DFT calculations, and will significantly advance our predictive modeling capabilities to treat complex materials systems. This presentation will discuss our recent advances towards addressing the aforementioned challenges. In particular, the development of computational methods and numerical algorithms for conducting fast and accurate large-scale DFT calculations using adaptive finite-element discretization will be presented, which form the basis for the recently released DFT-FE open-source code. The computational efficiency, scalability and performance of DFT-FE will be presented, which will demonstrate a significant outperformance of widely used plane-wave DFT codes. DFT studies on dislocations and biomolecules will be presented to showcase the capability of DFT-FE in handling large-scale systems. Finally, recent efforts, and related thoughts, towards developing a framework for a data-driven approach to improve the exchange-correlation description in DFT will be also discussed.

Nonlocal continuum models are in general integro-differential equations in place of the conventional partial differential equations. While nonlocal models show their effectiveness in modeling a number of anomalous and singular processes in physics and material sciences, for example, the peridynamics model of fracture mechanics, they also come with increased difficulty in computation with nonlocality involved. In this talk, we will give a review of the asymptotically compatible schemes for nonlocal models with a parameter dependence. Such numerical schemes are robust under the change of the nonlocal length parameter and are suitable for multiscale simulations where nonlocal and local models are coupled. We will discuss finite difference, finite element and collocation methods for nonlocal models as well as the related open questions for each type of the numerical methods.

We will discuss models used in classical molecular dynamics, and some mathematical questions raised by their simulations. In particular, we will present recent results on the connection between a metastable Markov process with values in a continuous state space (satisfying e.g. the Langevin or overdamped Langevin equation) and a jump Markov process with values in a discrete state space. This is useful to analyze and justify numerical methods which use the jump Markov process underlying a metastable dynamics as a support to efficiently sample the state-to-state dynamics (accelerated dynamics techniques à la A.F. Voter). It also provides a mathematical framework to justify the use of transition state theory and the Eyring-Kramers formula to build kinetic Monte Carlo or Markov state models.

Particle systems subject to long-range interactions can be described, for large numbers of particles, in terms of continuum models involving nonlocal energies. For radially symmetric interaction kernels, several authors have established qualitative properties of minimizers for this kind of energies. But what can be said for anisotropic kernels? Starting from an example describing dislocation interactions in metals, I will discuss how the anisotropy may affect the equilibrium measure and, in particular, its dimensionality.

A natural property to demand from discrete in time approximations to gradient flows is energy stability: Just like the exact evolution, the approximate evolution should decrease the cost function from one time step to the next. Often, approximation schemes with desirable (e.g. unconditional) energy stability, such as minimizing movements, are only first order accurate in time. We will discuss general (problem independent) procedures for boosting the order of accuracy of existing implicit and semi-implicit variational schemes for gradient flows while preserving their desirable stability properties, such as unconditional energy stability. The resulting high order versions are formulated only in terms of multiple calls of the original scheme per time step, and therefore also essentially preserve their per time step complexity.

Green’s functions play a central role in describing excited state electronic structures in quantum chemistry and materials science. At the heart of Green’s function computation is the solution of a linear system of size , where is the number of spin-orbitals in the quantum system. We will discuss how to use quantum linear system solvers to compute Green’s functions, and how to accelerate such computations using a new quantum primitive called the fast inversion, as well as preconditioning techniques.

Layers of two-dimensional materials stacked with a small twist-angle give rise to periodic beating patterns on a “moiré superlattice” scale much larger than the original lattice. This effective large-scale fundamental domain allows phenomena such as the fractal Hofstadter butterfly to be observed in crystalline materials at experimental magnetic fields. More recently, this new length scale has allowed experimentalists to observe new correlated electronic phases such as superconductivity at a lower electron density than previously accessible and has motivated an intense focus by theorists to develop models for this correlated behavior. We will present some mathematical and computational models for these experimental platforms and theoretical models.

We discuss continuum Schroedinger operators which are basic models of 2D-materials, like graphene; in its bulk form or deformed by edges (sharp terminations or domain walls). For non-magnetic and strongly non-magnetic systems we discuss the relationship to effective tight binding (discrete) Hamiltonians through a result on strong resolvent convergence. An application of this convergence is a result on the equality of topological (Fredholm) indices associated with continuum and discrete models (for bulk and edge systems). Finally, we discuss the construction of edge states in continuum systems with domain walls. Away from the tight binding regime there are resonant phenomena, and we conjecture that there are meta-stable (finite lifetime, but long-lived) edge states which slowly diffract into the bulk.

Simulations can be used to measure the properties of interfaces in materials. The central role of quantitative phase field simulations in this effort is illustrated by a rapid throughput method to determine grain boundary properties. By comparing the evolution of experimentally determined three-dimensional grain structures to those derived from simulation, we measure the reduced mobilities of thousands of grain boundaries. Using a time step from the experiment as an initial condition in a phase-field simulation, the computed structure is compared to that measured experimentally at a later time. An optimization technique is then used to find the reduced grain boundary mobilities of over 1300 grain boundaries in iron that yield the best match to the simulated microstructure. We find that the grain boundary mobilities are largely independent of the five macroscopic degrees of freedom given by the misorientation of the grains and the inclination of the grain boundary. The challenge of developing quantitatively accurate phase field simulations of grain growth will be highlighted, with an emphasis on the novel PDE’s suggested by methods that can account for the five degrees of freedom of the grain boundary energy.

Particles with diameters of nanometres to micrometres form the building blocks of many of the materials around us, and can be designed in a multitude of ways to form new ones. One challenge in simulating such particles is that the range over which they interact attractively, is often much shorter than their diameters, so the SDEs describing the particles’ dynamics are stiff, requiring timesteps much smaller than the timescales of interest. I will introduce methods aimed at accelerating these simulations, which simulate instead the limiting equations as the range of the attractive interaction goes to zero. In this limit a system of particles is described by a diffusion process on a collection of manifolds of different dimensions, connected by “sticky” boundary conditions. I will introduce methods that simulate such sticky diffusion processes directly, and discuss some ongoing challenges to extending these methods to high dimensions.

The work aims at developing new numerical tools in order to permit real-time and robust data assimilation that could then be used in various engineering activities. A specific targeted activity is the implementation of applications in which a continuous exchange between simulation tools and experimental measurements is envisioned to the end of creating retroactive control loops and online health monitoring on mechanical systems. In this context, and in order to take various uncertainty sources (modeling error, measurement noise,..) into account, a general stochastic methodology with Bayesian inference is considered. However, a well-known drawback of such an approach is the computational complexity which makes real-time simulation a difficult task. The research work thus proposes to couple Bayesian inference with attractive and advanced numerical techniques so that real-time and sequential assimilation can be envisioned. First, PGD model reduction [1] is introduced to facilitate the computation of the likelihood function, the uncertainty propagation through complex models, and the sampling of the posterior density. PGD builds a multi-parametric solution in an offline phase and leads to cost effective evaluation of the numerical model depending on parameters in the online inversion phase. Second, Transport Map sampling [2] is investigated as a substitute to classical MCMC procedures for posterior sampling. It is shown that this technique leads to deterministic computations, with clear convergence criteria, and that it is particularly suited to sequential data assimilation. Here again, the use of PGD model reduction highly facilitates the process by recovering gradient and Hessian information in a straightforward manner [3]. Third, and to increase robustness, on-the-fly correction of model bias is addressed in a stochastic context using data-based enrichment terms. The overall cost-effective methodology is applied and illustrated on a specific test-case dealing with real-time model updating for the control of a mechanical test involving damageable concrete structures with full-field measurements [4].

We investigate the emergence of rigid polycrystalline structures from atomistic particle systems. The atomic interaction is governed by a suitably normalized pair interaction energy, where the `sticky disk’ interaction potential models the atoms as hard spheres that interact when they are tangential. The discrete energy is frame invariant and no underlying reference lattice on the atomistic configurations is assumed. By means of Gamma-convergence, we characterize the asymptotic behavior of configurations with finite surface energy scaling in the infinite particle limit. The effective continuum theory is described in terms of a piecewise constant field delineating the local orientation and micro-translation of the configuration. The limiting energy is local and concentrated on the grain boundaries, i.e., on the boundaries of the zones where the underlying microscopic configuration has constant parameters. The corresponding surface energy density depends on the relative orientation of the two grains, their microscopic translation misfit, and the normal to the interface.

We characterize skyrmions in ultrathin ferromagnetic films as local minimizers of a reduced micromagnetic energy appropriate for quasi two-dimensional materials with perpendicular magnetic anisotropy and interfacial Dzyaloshinskii-Moriya interaction. The minimization is carried out in a suitable class of two-dimensional magnetization configurations that prevents the energy from going to negative infinity, while not imposing any restrictions on the spatial scale of the configuration. We first demonstrate existence of minimizers for an explicit range of the model parameters when the energy is dominated by the exchange energy. We then investigate the conformal limit, in which only the exchange energy survives and identify the asymptotic profiles of the skyrmions as degree 1 harmonic maps from the plane to the sphere, together with their radii, angles and energies. A byproduct of our analysis is a quantitative rigidity result for degree ±1 harmonic maps from the two-dimensional sphere to itself.

We will derive and explore the mass action kinetics of chemical reactions by employing an energetic variational approach. The dynamics of the system is determined through the choice of the free energy, the dissipation (the entropy production), as well as the kimenatics (conservation of species). The method enables us to capture the coupling and competition of various mechanisms, including mechanical effects such as diffusion, viscoelasticity in polymerical fluids and muscle contraction, as well as the thermal effects. This is a joint work with Bob Eisenberg, Pei Liu, Jan-Eric Sulzbach, Yiwei Wang and Tengfei Zhang.

We report on joint work with Carlos Roman and Etienne Sandier in which we study the onset of vortex lines in the 3D Ginzburg-Landau model of superconductivity with magnetic field and derive an interaction energy for them.