Roadmap on multiscale materials modeling

What is the point of molecular simulations?¹⁹ Broadly speaking there are two types of simulations: (1) simulations using 'toy models' designed to understand possible behaviors of classes of materials, and (2) quantitative simulations aimed at predicting the behavior of a specific technological material. There is a large gray area between these two extremes where the simulation is presented as describing a real material (such as pure single crystal copper), but typically this is unconfirmed. Molecular simulations suffer from the disadvantage that in many cases direct experimental validation is not possible due to the very small material systems that these simulations can model, and the very high loading rates necessitated by stable integration of Newton's equations of motion. The situation is improved in multiscale methods [3] which can reach longer length and time scales, but even there, direct experimental validation (or even comparisons among multiscale methods) are rarely done [4].

Peer review of articles reporting molecular simulations is based largely on evaluating whether the simulations appear to have been performed correctly based on the procedures reported by the authors and that the analysis and any theory developed to explain the results appear to be correct. It is not realistic to expect reviewers to redo simulations to verify correctness. This is analogous to peer review of experimental work. In both simulations and experiments, verification of the research is left to follow up work where other researchers attempt to reproduce the results and build on them. This requires that the readers of an article have all the information that they need to replicate the work. For molecular simulations, this means a complete characterization of the system, boundary conditions, and any simulation procedures used.

The ability to reproduce work is critical for the self-correcting mechanism of Science as explained above. It is also of great value to researchers themselves. Experimentalists are famous for maintaining meticulous lab notebooks that help them keep track of the large number of experiments (many unreported) that are necessary in order to understand a problem and obtain high-quality results. The field of molecular simulation (and simulations in general) do not have a similar culture. Students are typically not taught how to maintain order among the large numbers of preliminary simulations that they perform. Numerical simulations leave a wake of directories full of inputs and outputs with little or no documentation. Even the researcher who did the work (let alone other researchers) will find it difficult (and sometimes impossible) to go back to an earlier step, understand what was done, and reproduce the results. This culture is beginning to change with the emergence of workflow management tools such as AiidA [5] and Jupyter [6]. These tools make it possible to document a simulation and in principle reproduce it²⁰.

The key challenge for workflow management systems when it comes to classical molecular simulations is the interatomic model (IM). The IM is a computer program that receives as input a configuration of atoms (including information on coordinates, species, charges, etc) and outputs the energy and its derivatives (e.g. the negative derivative of the energy with respect to coordinates are the forces on the atoms). IMs have traditionally been implemented within specific simulation codes, such as the MD code LAMMPS [7]. However, this makes it impossible to ensure reproducibility (since simulation codes and IMs are continuously mutating) and very difficult and error prone to transfer IMs between different simulation platforms. A recent article [8] discussing software reuse and reproducibility used the materials simulation community use of LAMMPS as an example of dysfunction.

The issue of reliable IM archiving and portability is being addressed by the Open Knowledgebase of IMs (OpenKIM) project [9–11]. OpenKIM is a cyberinfrastructure funded by the National Science Foundation hosted at https://openkim.org. OpenKIM archives IMs, verifies their coding integrity through a series of 'Verification Checks' (e.g. the forces returned by the IM are checked against numerical differentiation of the energy), and tests them by computing their predictions for a variety of material properties using codes uploaded by the community. As a member of DataCite [12], OpenKIM issues a unique permanent digital object identifier (DOI) to every IM archived in openkim.org. Any modifications to an IM, such as parameter changes due to an improved fit, lead to a version change in openkim.org and a new DOI. The DOI can be cited in publications and incorporated into workflow managers to ensure that the exact IM cited in the work is downloaded and used, thereby ensuring reproducibility.

The issue of portability is addressed in OpenKIM through the development of an application programming interface (API) for communication between simulation codes (simulators) and IMs [13]. Simulators and IMs conforming to the KIM API work seamlessly together. The KIM API is cross-language compatible, currently supporting Fortran, C, C++ and Python, and is lightweight with negligible performance overhead in most cases. The KIM API is currently supported by a number of major simulators including ASAP, ASE, DL_POLY, GULP, LAMMPS, and the multiscale Quasicontinuum method [14]. The existence of this standard is important to ensure technology transfer throughout the community (by allowing an IM to be used in a range of codes) and encourages the development of new molecular simulation methods since by conforming to the KIM API new codes have instant access to a large pool of IMs.

Molecular simulations and multiscale methods are coming of age. Increasing computing power and the development of new highly-accurate IMs—and in particular machine-learning based IMs—is now making it possible to perform predictive simulations for real materials at meaningful length and time scales (see section 6). The potential inherent in these developments is being held back by current practices in the field. To advance, the molecular simulation community must embrace computing best practices, which include methods to ensure reproducibility and standards to allow for rapid sharing of new technologies.

The goals and expectations of computational materials science have evolved over the last few decades. Underlying this evolution is the range of, often unspoken, objectives of modeling. In some cases, the goal is to develop qualitative understanding of fundamental mechanisms and how those mechanisms interact to produce macroscopic behavior. The emerging goal is to inform materials design and qualification processes, where being quantitatively predictive is important. An increasing emphasis on computational materials science as a key component of the engineering process is exemplified by initiatives such as Integrated Computational Materials Engineering [15], the Materials Genome Initiative or NASA Vision 2040 [16]. In this role, modeling is employed to make decisions as opposed to expanding understanding. For modeling to be useful in a decision-making environment, an assessment of the reliability of the model predictions is essential. Thus, the growing interest in the development of UQ for materials modeling. Meaningful UQ is a prerequisite for meaningful Validation, the assessment of if—and how far—the results of a simulated model can be trusted to predict reality.

In thinking about UQ, it is helpful to identify the different types of uncertainty [17]. The simplest source of potential uncertainty is reliability of the numerical implementation. A portion of this is Verification, assessing that simulation codes correctly solve the equations that underlie the calculation. Associated with this is the practitioner's attention to detail in the application of the code, e.g. in such considerations as adequate k-point sample for electronic structure calculations or sufficient simulation time for MD simulations of statistical quantities such as correlation functions.

Another class of uncertainty in modeling predictions is aleatory uncertainty. This uncertainty arises due to inherent randomness in a physical process. An example of this is the variation in material properties for items that are processed in nominally the same manner. Though processing conditions are the same, each instantiation will differ microscopically. For example, while the grain structure may be similar, it will not be exactly the same in each case. This will lead to differences in the properties of each item. One way to treat aleatory uncertainty is to describe the response in terms of probability distributions. This is in contrast to traditional approach to materials modeling that focuses on prediction of average behavior.

The more challenging class of uncertainty is epistemic uncertainty, the prediction uncertainty that results from our incomplete knowledge. A simple aspect might be a lack or poor knowledge of key material input parameters, say the elastic constants of a new alloy. A more difficult aspect of epistemic uncertainty, fundamental to a multi-scale approach to modeling of materials behavior, is model form error. As reiterated throughout this Roadmap, behavior of the material at more detailed scale is synthesized into a coarser scale model. The use of a reduced, approximate model form clearly can lead to errors [18]. Model form errors are more difficult to assess than are parametric uncertainties. Parametric uncertainties can be estimated by sampling techniques. The error from neglect of detailed physics aspects usually cannot be directly quantified. Phase field models abstract constitutive relations from atomistic data (next section). The form of an interatomic potential determines how faithfully a MD calculation, neglecting detailed electronic structure, can replicate chemistry described by a density functional calculation (see section 6), Electronic structure calculations themselves introduce model form errors in choice of density functional that can only be crudely estimated [19]. The form of that model dictates the fidelity of the information transfer between scales.

A key challenge of UQ in materials modeling is that it requires a culture shift in the modeling community, especially at smaller, sub-continuum length scales. Historically, UQ issues have generally only received minimal attention. This is beginning to change. There is a small but growing body of literature addressing methods for UQ in materials modeling. For example, a recent study demonstrated that DFT codes have small numerical uncertainties by comparing results for a suite of model calculations [20]. Symposia at national meetings of major societies such as MRS and TMS as well as focused conferences are addressing the role of UQ and its future directions. Such signs are encouraging. Two changes will help drive this culture change. The first is in the education of material modelers. While attention to numerical issues is often discussed, the broader issues need to be incorporated into academic curricula. Scientific journals, such as Modeling and Simulation in Materials Science and Engineering, also have a role to play. Similar to how many journals require the inclusion of error bars on experimental data points, peer review criteria should be expanded to require a discussion of estimates of uncertainty in computational models or in their absence a justification for omitting them for that paper.

In considering a path forward for the development of UQ methodologies for materials, there are two types of challenges moving forward. The first assesses the impact of the approximations at a given length/time scale associated with a certain computational technique. For example, classical MD simulations are based on an assumed interatomic potential or force field [21, 22]. Quantifying the range of results from MD simulations from a sampling of similarly realistic interatomic potentials would be a key component of an overall UQ method. Based on recent UQ symposia, the majority of on-going efforts address this class of problems [23]. This is a sensible starting point for the field because the questions are more clearly defined.

The second broad challenge arises in the context of multi-scale materials modeling [24]. Conceptually, information obtained from simulation(s) at smaller length/time scales are synthesized and used to inform models at higher scales [25–27]. Inherently, there is a loss of information from the smaller scales. The challenge is quantifying or at least bounding the prediction uncertainty that results from a chain of modeling modalities moving from electronic scales up to engineering scales [28, 29]. While conceptual multi-scale modeling hierarchies exist, complete quantitative multi-scale modeling hierarchies are actually rare [30, 31]. A major challenge is the transformation of information between the scales. As a simple example of this, consider the treatment of temperature at different scales. At mesoscale and continuum scales, temperature is typically a scalar field variable. In atomistic simulations, thermal energy exists in the random kinetic energy of the particles and temperature is a derived quantity or a boundary condition. While the treatment of temperature is generally understood, this demonstrates the type of conceptual challenge that can exist moving between scales. Similarly, in looking at mechanical response, the dynamics of ensembles of atoms is mapped onto a set of characteristic defects (vacancies, interstitials, dislocations, grain boundaries, ...). The evolution of these defects is often further synthesized into higher level constructs like shear bands. Such sequences are useful if they capture the essential behavior, but can fail if the synthesized information fails to capture essential features. A UQ analysis of such a hierarchy is clearly a formidable challenge which requires sophisticated error propagation methods [32, 33]. Another potential use of UQ concepts is in the reverse direction, using UQ analysis and observations at higher length scales to pinpoint knowledge gaps at lower length scales or poor model form used for information transfer [34]. And then integrating this UQ into a meaningful system of Validation to certify the predictive accuracy of materials simulations in prescribed regimes. Addressing these questions represents the long-term research area for UQ in the context of materials modeling.

The phase-field method over the years has established as the method of choice for simulation of microstructure evolution at the nano- and mesoscopic scale. Nevertheless, the term 'phase-field' may provoke some misunderstanding. Traditionally it denotes the region in an alloy phase-diagram where an individual crystallographic phase is stable for a given composition, pressure and temperature. It had been firstly used by Langer in 1978 [35] for the solution of a nonlinear wave equation, known in physics as the 'soliton' [36]. This field theoretical solution was applied to dendritic solidification, a phase transformation problem, therefrom the name 'phase-field'. The soliton solution simply helps to propagate the solidification front in a numerical simulation, meaning the change of phase from solid to liquid over space and time. The width of the transformation front, the 'diffuseness' of the phase-field, in this regard has no physical meaning. The theory is agnostic of an intrinsic scale. It lives at the meso-scale, i.e. large compared to atoms, but small compared to the sample dimensions, the size of a casting in solidification. Even 40 years later, there is still a debate about the interpretation of 'phase-field' as a microscopic order parameter model, or an elegant numerical tool. For a more in-depth review of the history, see [37–39]. The future of 'phase-field' clearly lies in 'making true its promises'. It is considered a thermodynamically consistent theory in the tradition of variational approaches of classical mechanics. It offers a consistent framework to incorporate interfaces and kinetics into thermodynamics [40]. Augmented by most advanced models of diffusional and advective transport, micro-elasticity and plasticity, magnetism and ionic mass transport, it will bring the important phenomenon of 'evolving microstructures' into full-field models of materials behavior. A phase-field model transfers atomistic scale properties, like interface energy anisotropy, into mesoscopic scale microstructures. From here it transfers into macroscopic scale materials properties. The microstructures and their evolution during processing and service determine materials properties. They are evaluated by direct numerical simulation of materials behavior under load. In the following section, I will highlight some research issues and challenges for future developments. They are based on my own experience and interest. The applicability of 'phase-field' as a scale bridging approach is, however, much broader.

Phase-field, as discussed above, poses a promise: We have a thermodynamically consistent theory to simulate materials behavior by solving well-posed partial differential equations (PDEs) on computers. The solutions for 3D problems deserve huge computational resources. Strategies of massive parallelization, intelligent time stepping strategies as well as efficient adaptive meshing schemes are developed to a high level. I consider these issues as 'technical' and have no doubt that in the case of real application the necessary resources will be made available. The real challenge, the 'big research issues', lie in the integration of best available constitutive relations for bulk materials with most advanced models for interfaces, their static and kinetic properties. Integration in this respect means that 'bulk' and 'interface' must be considered in common! The best example is diffusion controlled dendritic solidification: The morphology of the dendritic structure is determined by diffusion in the bulk melt around the growing solid, but is intrinsically linked to interface energy anisotropy living at the atomistic scale of the solid liquid interface. Solid state interface properties are even more involved, and their impact on microstructure evolution is out of question. Martensitic and mixed mode transformations critically depend on the balance of bulk and interface related phenomena. 'Phase-field' offers a platform, but it does not provide the solution by itself. It requires good input from atomistically informed constitutive relations and experiments.

A big challenge, however, is the consideration of the whole life cycle of a material: from production through service to failure (see figure 1). Materials are not 'dead bodies'. Their microstructure evolves continuously during the whole lifetime cycle of the material including refurbishment and, in general, also recycling. So, in general, different routes of production, as sketched in the figure, will lead to different microstructures and to different properties. Phase-field simulations have been applied to investigate most individual steps in this cycle: from solidification to failure. Solidification from the homogeneous melt sets the initial structure, at least for almost all metallic materials. The microsegregation, created during solidification, will persist in most conventional heat treatments, even after rolling, if slow diffusing elements as Mn in steel are considered. Further transformation steps should consider this information as a starting configuration. In particular, predictive simulation of crack initiation and failure will only be possible, in general, if important microstructural information through the lifetime cycle of the material is considered. This consistent through-process simulation is still a challenge for future applications. Phase-field simulation in combination with most advanced micromechanical models offer the possibility to attack this challenge. First steps in this direction by atomistically informed full-field simulation of quenching, tempering and testing of tempered martensite have been published recently [41, 42].

Zoom In Zoom Out Reset image size

The phase-field method as discussed above incorporates interfaces and kinetics into thermodynamics. Since it is a continuum method formulated as partial differential equations and resting on sophisticated constitutive relations, it requires a maximum of input compared to other methods, as discussed in figure 2. For bulk thermodynamic properties well established CALPHAD databases exist. Here the challenge is to address also metastable regions in the phase diagram as well as new phases and exotic materials. First principles calculations can help to supplement additional information like chemical mobilities, activation barriers for nucleation and fault formation. Databases of interfacial properties, static and dynamic, are still rare and in the state of development. Databases for mechanical properties are diverse. Mostly they relate to 'materials' as a homogeneous medium, not specific to the microstructure and the properties of individual components of the microstructure. In all these cases, improved data and models are needed urgently. Also commonly accepted standards of constitutive relations and materials data have to be developed. Actual activities in data mining and materials informatics must be utilized. The output of phase-field simulations to macroscale simulations would be local constitutive relations considering microstructural information. A last and very important issue is numerical accuracy and benchmark problems to be accepted by the community [43]. In both cases, the phase-field community must team up with the communities of applied mathematics and continuum mechanics, in particular micromechanical modeling and simulation, see the following section 5, to realize the necessary demand of accuracy and efficiency in solving the coupled multi-physics problems or evolving microstructures in real materials.

Zoom In Zoom Out Reset image size

The phase-field method bridges scales in several respects. It bridges from atomistic ordering to long-range transport. It bridges from microstructures to macroscopic materials properties. It bridges from physics to continuum mechanics and to engineering applications. Future application of phase-field can be found in fundamental research as well as in applied research. Fundamental aspects relate to pattern formation in various classes of phase transformations. Applied aspects relate to everyday engineering problems in metallurgy, processing of ceramics or functional materials such as magnetic microstructures, or ferroelectrics. Also problems in geoscience and biology are in the range of applications. All of this in combination with best constitutive relations, best data and best numerics.

Multiscale models of plasticity aim to predict the deformation behavior of metals and alloys suitable for engineering applications in the plastic regime (e.g. yield surface, strain hardening, texture evolution, creep rate, ductility limit, fatigue life) based on fundamental physics of atoms and crystalline defects. In a broader context, the goal is to reach a quantitative understanding of the microstructure-property relations, sufficient to enable prediction of beneficial and detrimental microstructural features. Successfully developed multiscale models of plasticity will have profound impacts on a wide range of engineering applications and industries. For example, they can inform the design of manufacturing processes, such as extrusion and forming, leading to better predictions of margins and safety factors respective to different failure modes. Multiscale models of plasticity can accelerate the design of high-performance materials, such as those used in high temperature applications like gas turbines and lightweight materials used in aeronautic and automotive industries. They can predict the performance of materials under extreme environments (such as inside nuclear reactors) in which experiments are very difficult or impossible to perform. They are also expected to play a vital role in establishing metal additive manufacturing (3D printing) as a reliable process to produce parts within acceptable property tolerances.

While many different plastic deformation mechanisms exist in crystalline solids (e.g. twinning, phase transformations, grain boundary sliding), slip induced by dislocation motion is dominant under most conditions. Therefore, a predictive model of plasticity requires understanding fundamental dislocation physics and dislocation interactions with other defect microstructures in the material (e.g. other dislocations, solute atoms, point defects, radiation defects, precipitates, grain/twin boundaries). Because these microstructural processes span a wide range of length and time scales, they exceed the capacity of any single computational model. Many models have been developed, such as atomistic models based on first-principles (e.g. DFT) and empirical potentials (e.g. MD), discrete dislocation dynamics (DDD) and continuum dislocation dynamics at the mesoscale, and crystal plasticity (CP) models at the polycrystalline microstructural scale [44]. These models need to be meaningfully connected to each other to obtain a predictive framework for multiscale modeling of plasticity. The most outstanding problem today is the lack of quantitative connections between CP models with the lower-scale dislocation models. As a result, existing CP models used in engineering applications are still phenomenological, while evidence continues to mount that they can make inaccurate predictions under realistically complex scenarios [45, 46].

There are three major challenges that must be overcome, in order to establish a successful framework for multiscale modeling of plasticity. The first challenge is to connect computational models of defect dynamics and experimental measurements of plastic deformation. A direct comparison between predictions and experiments under identical conditions would not only provide much needed validation of theory but also calibrate model parameters that may be impossible to determine from first principles. A promising approach is to start from simpler cases (e.g. pure single crystals) and progress towards more challenging ones (e.g. alloys then polycrystals), as illustrated in figure 3. The connection may be first established at high strain rates (e.g. >10² s⁻¹), by directly comparing DDD simulations [47] with Kolsky bar or micro-pillar compression experiments, and then expanded towards lower strain rates (e.g. <10⁻¹ s⁻¹).

Zoom In Zoom Out Reset image size

The second challenge is to quantitatively connect computational models at different scales, starting from a pure single crystal, as shown in figure 4. As always, the purpose is to benchmark and calibrate an upper scale model against a more fundamental, lower scale model. The key premise of the DDD method is that the response to straining of a statistically representative ensemble of dislocations can be assembled from the motion of its individual constituent dislocations. It remains unclear if and how much is 'lost in translation' in transferring knowledge gained in atomistic simulations of individual dislocations to DDD simulations of CP. Recent advances in direct ultra-scale atomistic simulations of CP [48] has reached simulation cell sizes of ∼1 μm containing up to 10⁶ dislocation lines. These ultra-scale simulations, coupled with efficient and accurate methods to extract dislocations from MD snapshots [49], allow direct comparison with DDD simulations in terms of, e.g. dislocation network structure, dislocation mobility and multiplication rates. Such comparisons provide a critical test case and a useful proving ground for improving physical fidelity and predictive accuracy of DDD simulations.

Zoom In Zoom Out Reset image size

The third challenge is to embrace the complex nature of real engineering alloys, i.e. realistically dirty materials. Keeping track of all the interactions between various defects is a daunting task. Even if an atomistic simulation could be carried out for an arbitrary defect configuration (e.g. dislocation-GB interaction), the total number of distinct configurations is too large to be practically considered. Therefore, the fundamental physics concerning the rules for defect interaction need to be accounted for in a statistical rather than an exhaustive manner. A theoretical framework is still lacking for constructing sufficiently accurate statistical models.

The convergence of several breakthroughs in computational and experimental capabilities in recent years has moved us much closer toward the goal of physics-based plasticity models. Further advances along these lines are needed to realize the full potential of multiscale modeling of plasticity.

First, the emergence of new computing architectures, such as graphical processing units (GPUs), has made a major impact in many fields of science and technology. For example, the use of GPUs, together with advanced time-stepping algorithms, has resulted in orders-of-magnitude increase in the efficiency of DDD simulations. MD simulations of metal plasticity [48] on million-core CPUs has enabled statistically meaningful comparisons between MD and DDD models. Therefore, it is essential for the developers of plasticity models across all scales to promptly take advantage of new computational architectures as they emerge, and to develop new algorithms that scale well on these new platforms. For example, further efficiency gains by several orders of magnitudes are needed for DDD simulations to reach quasi-static strain rates for fcc metals.

Second, the recent breakthroughs in microscopy have revealed microstructural details about materials that were previously unavailable during plastic deformation. For example, the near-field high-energy x-ray diffraction microscopy technique [50] at the advanced photon source has allowed the local lattice orientation of polycrystals to be followed as a function of plastic strain, and be compared with CP predictions. Time is ripe for quantitative comparisons between predictions from defect dynamics models and the wealth of microstructural information revealed by modern experimental techniques, e.g. 3D transmission electron microscopy (TEM), Laue micro-diffraction, Bragg coherent diffraction imaging (BCDI), high resolution digital image correlation, etc. Such comparisons would be greatly facilitated by the ability to directly simulate experimental images from the snapshots of the defect dynamics models [51].

Third, as both experiments and simulations are generating data at an unprecedented rate, tackling metal plasticity using the data science/ML approach appears highly promising. The adoption of data-driven methods is already happening in the broader field of computational materials science [52], e.g. the search for desirable alloy compositions based on first-principles datasets in the Materials Genome Initiative, as well as in the more specific context of dislocation simulations [53]. To take advantage of data analytics for understanding plasticity, it is necessary to develop platforms and protocols (e.g. through collaboration with computer/data scientists) that facilitate the exchange and mining of microstructural data, which are highly diverse and complex. High throughput on-the-fly computational techniques should be developed and combined with defect/microstructure simulations to efficiently span vast parameter spaces and to sample statistically representative ensembles. The goal is to identify key features and to test hypotheses generated by computation and experiments to aid the development of physics-based continuum models of plasticity.

Given the rapid progress in computational and experimental techniques, a new generation of multiscale models of CP is expected to emerge over the next 10-to-15 years that would connect defect physics with engineering-scale predictions, are validated by experiments, and offer valuable recommendations for materials processing and design. We note that multiscale models making precise microstructure-property predictions solely from first-principles without any 'tunable parameters' are perhaps unrealistic and likely unnecessary for engineering applications. Instead, the goal of the multiscale model should be to provide physics-based answers to questions (e.g. regarding qualitative or semi-quantitative trends) for specific material systems given the available experimental observations at various scales. Finally, we note that once multiscale models and coarse graining techniques are developed, UQ—in terms of errors introduced by specific models and by the coarse graining algorithms used to upscale information—will be an essential step to confidently providing such physics-based answers.

Simulating materials, especially fluids, using MD started soon after the birth of electronic computing. In time, this new tool allowed not only the calculation of very complex material properties, such as space and time correlation functions, but also to conduct 'computational experiments', in which no one is entirely sure what would happen before it is run [54].

In parallel, the ever more precise understanding of the interactions between large numbers of electrons and nuclei, as embodied in so-called first principles quantum mechanical methods and their computational implementation has revolutionized our ability to predict the properties of solid materials and individual molecules, and to use these predictions to understand phenomena at all length scales, including alloy design, corrosion, spectroscopy and transport properties, just to name a few [55].

Although such 'bottom-up' methods are gaining ground all the time, the length and/or time scales involved in MD simulations that are worth doing (hundreds to thousands of atoms for thousands to millions of individual time steps) impose such high computational costs that many such simulations are arguably out of reach for all but users of the largest supercomputers. In situations when the alternative, a length scale free description in terms of continuous fields is certain to fail to capture the correct mechanisms, we find the niche for using 'interatomic potentials', i.e. empirical, simplified models of forces acting on atomic nuclei, via a potential energy written as an explicit function of nuclear positions. The potential is supposed to include implicitly the energies of the electrons that are assumed to have relaxed into their ground state and follow adiabatically the slower evolution of nuclear coordinates.

Except in a few very simple situations, there are no theories on what functional forms such a potential should use. Practical models have been made using a combination of intuition, guess-work, and some trial and error. The assumed functional forms got more complicated over the decades, and nonlinear empirical parameters proliferated. Up until recently, it was widely felt that such potentials have reached a 'plateau', in terms of accuracy, reliability, and in general usefulness. Even for simple materials, although trends between them were captured, specific defect energies were too far off to be predictive, and the ability to draw valid conclusions relied in a heady mix of experience, artful use of transferability (fitting to one property and calculating another) and no doubt in some cases, just luck. More complex materials, such as oxides, interfaces, chemically modified surfaces, were largely out of bounds.

A new direction was taken starting about ten years ago, using the newly popularized tools of machine learning (ML): non-parametric function fitting in many dimensions (being conscious that when it comes to the number of dimensions, often one person's many is another's few, in this case I use 'many' to refer to tens to hundreds of dimensions) [56, 57]. Casting the problem of constructing interatomic potentials as a special kind of 'learning task' in which training data is generated using expensive first principles electronic structure calculations. This kind of fitting is sometime referred to as 'surrogate modeling'. There are some key differences with respect to the typical problems in ML. On the one hand, not only an arbitrary amount of essentially noise-free training data can be generated (at a fixed cost per item), but even the location of the training data can be chosen arbitrarily. On the other hand, accuracy demands are rather high: it turns out that '99%' accuracy in the potential energy of a 100-atom system is not particularly useful, and in order to predict properties better than existing models, ten or hundred times more accurate fits are often needed. Furthermore, such accuracy measures are not even that useful if understood in the statistical sense: while a model that makes a large error only very occasionally might be good for many independent tasks, this is not the case for much of materials modeling. For example, a large error in a MD simulation occurring just at a single time step could throw off the entire subsequent trajectory (for example by trapping it in an unphysical local minimum) and thus render the whole simulation useless. The key here is that these models are not evaluated on inputs from known and independently defined and generated probability distributions, but the models themselves are used to actually generate the distributions (or measures) on which the observables are evaluated. Transitions driven by the models (e.g. using Markov chain Monte Carlo (MC) or MD) are used to generate invariant measures, and errors in transitions that are small or rare in the statistical sense can lead to large errors in the invariant measures and even to broken ergodicity and thus to large errors in observables.

The current status in the capability of these new ML potentials is roughly as follows.

• •  
  Both shallow neural networks, using a wide variety of atomistic descriptors, and kernel learning (and also linear regression), using specially designed basis functions, have been successfully used to fit accurate potentials for a wide range of materials (see articles in [58] for recent examples). Impacts on real materials science problems are beginning to appear, notably for amorphous materials among others [59].
• •  
  The data are total energies and gradients, typically computed using DFT, or sometimes even more accurate wavefunction based quantum chemical methods.
• •  
  A large fraction of published works are still 'proof of fit' type, with little attempt at transferability: the accuracy of the potentials are tested on configurations very similar to the training set, often generated by the same protocol.

A critical ingredient of the non-parametric fits is the way in which the nuclear positions are represented: this needs to respect basic symmetries of the target potential energy function with respect to translation, rotation and permutation of like atoms. Almost all representations in current use start with a fixed radius neighborhood of an atom, represented as an 'atom-density' and project this onto a rotationally invariant finite dimensional basis, using e.g. spherical harmonics. (The one exception uses wavelet transforms to capture the global atom-density on multiple length scales [60].) This brings with it the first of a number of challenges.

• How can long range electrostatic (and dispersive) interactions be made part of the ML model? If the radius of the neighborhood in the representation is significantly enlarged, the dimensionality of the fit quickly becomes unmanageable, and the rotational invariants also lose their appeal. With a finite radius, not only is electrostatics not properly described, but long range charge transfer is also missed.A disconcerting feature of ML models is their 'fragility': predictions made even not very far outside the region of the training data are essentially random, that is the 'price' paid for accuracy within the region of training data using generic functional forms with a very large number of free parameters. The corresponding challenge is thus:
• How can we ensure that ML potentials, with their very narrow range of transferability, do not lead to non-sensical predictions, which would contaminate simulation results? Or, turning it around, can we ever make ML potentials that correctly describe a material in a very wide range of (perhaps all sensible) configurations? A first attempt of this is in [61], albeit for a single component material (silicon) and using a 'hand built' database. Can we design protocols for automatically generating training databases suitable for a given scientific problem? Can we quantify the extent to which a training database covers the relevant part of configuration space? The distance metric defined between configurations that is implicit in kernel-based fits would appear to be useful here. And finally, can we create a single ML model that covers a wide variety of materials?There are many problems in atomic scale materials modeling that cannot be tackled just by having a potential energy function of atomic positions. The challenge there is to extend the non-parametric high dimensional fitting approach to include:
• Spin degrees of freedom and magnetic interactions, electronic entropy, multiple oxidation states, excited electronic state potential energy surfaces, etc. These are situations where the electrons which were eliminated in defining interatomic potentials appear to make a comeback, but at the same time the reintroduction of fully explicit electronic degrees of freedom may not be necessary.

We now indicate some possible directions that could be taken in order to overcome the above challenges. For molecular systems, electrostatics has long been described by multipole expansions, and the corresponding response functions could be calculated and fitted using ML models. For solids, this is considerably more challenging, because there is no unique way to partition a strongly bound solid into a disjoint set of interacting electrostatic multipole sources. ML could be used to find the best such partitioning, and to fit its response functions, and such a model can then be added onto the current machinery for fitting the remaining short-range interactions. An example along these lines has been published originally for NaCl [62] and subsequently for other ionic materials.

The challenge of transferability is perhaps the thorniest. Its general solution, not only in materials modeling, but more widely in ML, might be in the form of creating models that operate in a hierarchy of spaces, starting with lower dimensional representations which afford less accurate but robust predictions, and sequentially refining this using richer representations and more accurate fits. For atoms, a good guess at some of these lower dimensional representations might be the interaction of pairs of atoms, then triplets, etc (a well worn idea in materials modeling, see also [63] for an example use with ML), but to go further we need to glean the right representations from the data itself to avoid the exponential blowup of the body-order expansion.

Finally, ML potentials may be able to link the worlds of reactive materials modeling, dominated by DFT that more or less correctly describes bond forming and bond breaking, with the world of wavefunction based quantum chemistry, which is the right approach when exquisite accuracy is required e.g. sufficient to obtain equations of state and dynamical properties of molecular liquids. Ultimately, along with the advances in modeling hard materials, this approach might also lead to the 'holy grail' that is a reactive organic molecular force field with the 'gold standard' accuracy of coupled cluster theory. See [64] for a first stab at this.

In recent years, a large number of atomistic-continuum coupling methods have been developed with the goal of reproducing the results of the fully-atomistic model at lower computational cost, which is particularly important in problems where a large size of computational domain is required, e.g. to deal with long-range stress fields around crystalline defects. Moreover, since the time scales of dynamics simulations of the atomistic system such as MD are limited to sub-microseconds due to the short vibration period of atoms (typically orders of picoseconds) and many systems exhibit rate-dependent behaviors (e.g. hardness in nanoindentation, friction in sliding test, etc), it has been a long-time dream to run dynamics simulations of atomistic models for a time length closer to macroscopic scales. While several noble schemes have been developed to extend the MD time scale, including graphics processing unit (GPU)-based algorithms [65–69], it was just recently that temporal acceleration methods began to be combined with coupled atomistic-continuum approaches. Below the scope of these temporal and spatial multiscale methods will be discussed by reviewing three outstanding examples.

The first example is the study of temperature-dependent dislocation nucleation at the crack tip of face-centered cubic metals by Warner and Curtin [70]. For spatial coarse-graining they employed the finite-temperature CADD (Coupled Atomistic Discrete Dislocation) method (figure 5(a)), where atoms in the atomistic domain dynamically evolve as in MD with some atoms near the atomistic/continuum interface thermostated to prevent unphysical wave reflection whereas the continuum fields are updated in a quasistatic way using the mean positions of the interface atoms. The acceleration in time was achieved using the parallel replica method [66] where statistically equivalent multiple systems are simultaneously monitored until one of them exhibits a transition. When it happens, the total simulation time is calculated as the sum of the times of all the replica systems. The acceleration factor scales, thus, linearly with the number of replicas, i.e. the computational resources that can be allocated to run these replicas at the same time.

Zoom In Zoom Out Reset image size

The second example, called hyper-QC, was constructed by combining the finite-temperature quasicontinuum (QC) method (hot-QC) for spatial extension with the hyperdynamics method for temporal acceleration [71, 72] (figure 5(b)). In hot-QC, an effective potential for representative atoms containing all atoms in the atomistic domain and a small subset of atoms in the continuum domain is defined based on the local-harmonic approximation of the free energy, which can reproduce the canonical ensemble equilibrium properties. These representative atoms, whether in the atomistic or continuum domains, dynamically evolve as in NVT MD simulations. Moreover, a bias potential is added to the original potential energy surface to reduce the energy barriers so that transitions are expedited. It was formally proved that under the assumptions of the transition state theory, hyperdynamics simulations [65] can preserve the original state-to-state dynamics, i.e. the biased system visits each metastable state with the same probabilities as in the original system. The acceleration factor depends on the quality of the bias potential.

The third example was based on the maximum entropy (max-ent) formalism and the mean field approximation, which provide the governing equations for the dynamic evolution of the mean positions of atoms [73, 74]. Since the short-time atomic-scale vibrational modes are already averaged out in this formalism, the resultant trajectories are smooth on microscopic time scales so that much larger time steps than used in conventional MD simulations can be used, leading to the extension in overall time length. The spatial coarse-graining was realized by adopting the cluster-QC method, which is a full non-local version of the spatial multiscale method (figure 5(c)).

In short, there exist dozens of different coupled atomistic-continuum methods, but only a few of them have been extended for accelerated dynamics simulations. A couple of these temporal and spatial multiscale methods employ accelerated schemes that are developed for fully-atomistic models so that atomic-scale thermal vibrational modes are still retained. In contrast, there also exists a method where the coarse-graining is applied to the time-domain such that the simulated dynamics is for the mean positions of atoms enabling longer-time evolutions.

Whereas a proper dynamic coupling of atomistic and continuum domains is an emerging area of research including several key challenges such as heat exchange between domains (from atomistic to continuum and vice versa), in this section we focus only on the time acceleration issue. An atomistic system in the solid state often evolves through 'infrequent' thermally-activated transitions from one potential energy basin (state) to another, i.e. the system spends most of time near basins before quickly transiting into other adjacent states. This is the case where all energy basins are separated by large energy barriers as seen in figure 6(a) and, for example, dislocation nucleation is understood as a thermally-activated process. Many acceleration schemes aim to make such transitions occur at an expedited pace while preserving the relative transition probabilities among the neighboring states. In light of the Arrhenius type dependence of the transition rate on energy barrier and temperature, exp(−ΔV/k_BT), two natural ideas for accelerating thermally-activated events are either lowering energy barrier (as in hyperdynamics [65]) or increasing temperature (as adopted in temperature-accelerated dynamics (TAD) [67]). As discussed above, hyperdynamics was coupled with a spatial multiscale method (hyper-QC). Since the boost factor in both hyperdynamics and hyper-QC depends on the bias potential, the key challenge in hyperdynamics is to develop bias potentials that are computationally inexpensive, but rigorous enough not to distort the original state-to-state dynamics. The original method proposed by Voter using the eigenvalue/eigenvector of the Hessian matrix is very versatile, but its computational cost is not trivial [65]. Even though several alternative cheaper schemes have been proposed, they are not as robust as the original method. Thus, it remains a major challenge to develop a versatile, but less expensive bias potential for more popular application of this method. On the contrary, the acceleration in TAD does not require any modification of the potential energy surface, but can be achieved by simply raising temperature, although running a TAD simulation is not trivial because of some computationally expensive tasks such as saddle point searching. Unlike hyperdynamics or TAD, the parallel replica method does not employ any assumption or demand high-cost computational tasks, but its acceleration factor is limited by the available computational resources.

Zoom In Zoom Out Reset image size

Another critical consideration in temporal acceleration is that not all dynamic evolutions of the atomistic system is through infrequent transitions between states surrounded by large energy barriers. Instead, the potential energy surface may have energy barriers of various scales as exemplified in figure 6(b), making the strategy of the aforementioned acceleration schemes very inefficient. Moreover, in many cases systems are 'driven' by time-changing external parameters (e.g. indentation, fracture test, uniaxial tensile test, etc) so that the potential energy landscape also changes in time. As illustrated in figure 6(c) the initial energy barriers may increase or decrease as time passes by and new wells may emerge with others disappearing. In all these cases multiple time-scales are so entangled that selectively accelerating only 'infrequent' events would be difficult or provide little gain in simulation time. Therefore, it is a major challenge to develop a more comprehensive theoretical framework in temporal acceleration that can deal with not only infrequent events, but also systems with multiple entangled time-scales.

With regard to the direction of temporal acceleration, the formulation based on the dynamical evolution of the mean positions (the max-ent formalism) appears very promising, but it is an open question whether the trajectory predicted by the max-ent formalism also agrees with the state-to-state dynamics from the direct atomic-time scale simulations. A more thorough unified theoretical formulation would able to answer these question in the future. Moreover, it should be noted that for the past several decades computational science and engineering research has been greatly benefited from the advancement of computer technology in both hardware and software. This progress in computing capability will remain as a key factor in the development of atomistic simulations.

To sum up, employing temporal acceleration in coupled atomistic and continuum methods is a highly promising strategy to broaden the range of potential applications of atomistic methods by providing a more efficient utilization of computational resources. So far only a few spatial multiscale methods have been coupled with the temporal acceleration methods. Coupling temporal acceleration and spatial coarse-graining strategies has synergistic effects to each other so that any computational gain obtained from the spatial coarse-graining scheme enables longer simulations and vice versa (see also section 9). In harmony with the continuous improvement of computing power, we envision that the theoretical and methodological development in atomistic simulation will continue to provide an opportunity for deeper understanding of our physical world.

Computational multiscale modeling techniques are essential for accurate predictions of how a material's macroscopic, effective properties emerge from microscale mechanisms and architectures in solids and structures (here and in the following, the terms 'macro' and 'micro' do not refer to particular length scales but to the, respectively, larger and smaller scales of interest to be linked by multiscale techniques). One of two strategies is commonly pursued: hierarchical scale-bridging assumes a separation of scales between micro- and macroscales and is hence based on the passing of effective information between scales. Classical examples are homogenization techniques which average over a representative volume element (RVE) at the microscale in order to link average kinematic and kinetic variables—either a priori or on the fly. Such techniques can be applied at several levels and have been successfully deployed not only at continuum scales but also down to the atomic scale where averaging in a statistical sense allows for the extraction of effective material behavior to be used as an input at larger scales [75, 76]. If a separation of scales cannot be assumed (e.g. when localization into macroscale features occurs on the microscale, or when material and structural feature sizes merge as in nanoscale metamaterials), then concurrent scale-bridging techniques are the method of choice. These decompose the simulation domain into two or more subdomains, each of which treated simultaneously by a different model. Microscale resolution is thereby confined to subdomains of interest, while—away from those regions—a less accurate description is introduced for computational efficiency. An example is the CADD method [77] which couples atomistic regions (treated by MD) with continuum regions (described by DDD). MAAD [78] couples a total of three constitutive frameworks: tight-binding, MD, and finite elements (FEs); see figure 7. Alternatively, the concept of coarse-graining is a hybrid technique that applies efficient continuum concepts to a lower-scale, discrete constitutive description. An example is the family of QC methods which retain full atomistic resolution in spatial subdomains of interest, while continuum-level interpolation schemes are used as kinematic constraints in the remaining subdomains in order to significantly reduce the number of degrees of freedom and thus gain efficiency [79]. The list of hierarchical and concurrent multiscale techniques is long and continuously growing, which demonstrates the demand for such methods in the computational design of materials systems—from optimizing process parameters to obtain favorable microstructures in conventional materials to engineering the microscale architecture of novel metamaterials with tailored macroscale properties. Each scale-bridging technique provides a compromise between accuracy and efficiency, and a myriad of open challenges invite further developments.

Zoom In Zoom Out Reset image size

The following provides a non-exhaustive list of current and future challenges in this field.

• Homogenization at continuum scales is probably the oldest multiscale technique and has matured considerably with the use of computational techniques such as the FE method. Yet, dynamic homogenization (i.e. the homogenization of time-dependent material behavior including inertial effects) still presents challenges, especially in the presence of nonlinearity and transient effects. With the emergence of metamaterials, this challenge has gained momentum (see also section 13), and only few recent approaches have demonstrated success [80, 81].
• Homogenization requires the identification of an RVE whose size is usually chosen for reasons of efficiency without meeting its statistical requirements. Statistically similar RVEs have been proposed for continua [82], while a proper definition of stochastically arranged matter (such as network models for soft tissues) is challenging and calls for non-trivial compromises between ensemble averaging and ensemble enlargement.
• If there is no or limited separation of scales (i.e. if micro- and macroscale characteristic feature sizes become comparable), higher-order homogenization and enriched continua provide powerful tools that have not been fully explored. Fracture (or, more generally, localization) on the microscale is a related problem where the separation of scales gradually fails as microscale cracks grow into macro-cracks. This calls for RVE adaptation methods and related techniques that apply in case of limited scale separation.
• Computational efficiency has been gained by techniques of model order reduction [83], which reduce the microscale number of degrees of freedom by introducing intelligent kinematic constraints specific to a given microstructure. Extensions to, e.g. finite deformations, dynamics, and FFT-based spectral representations are ongoing developments.
• Discrete-to-continuum coupling techniques are essential when descending to the smallest of all scales, since atomistic and quantum mechanical techniques (such as MD and DFT, e.g.) are inherently discrete in their description. While spatial coupling has been accomplished in a variety of ways (see the MADD, QC, CADD, BDM, BSM, CAC, etc, techniques), an ongoing challenge arises in the context of temporal scale-bridging (see e.g. section 9 and section 7). For example, atomic ensembles evolve at the level of femtoseconds, governed by the discrete equations of motion and with atomic vibrations indicating the system temperature as known from statistical mechanics. Continuum models operate at significantly larger time scales (as required by most applications), governed by deterministic systems of PDEs and with temperature being a continuous field. Linking those two distinct descriptions is challenging and has been addressed, e.g. by applying continuum thermodynamics at the level of atoms [84] or by phase-space formulations and statistical averaging [85]. While upscaling in space is a traditional concept in multiscale modeling, upscaling in time of the system kinetics is a major challenge across scientific disciplines. For example for the coarse-graining of atomistic systems, modeling heat and mass transport as well as defect interactions over relevant time scales still poses major challenges. First derived for complex fluids, the general equation for non-equilibrium reversible-irreversible coupling [86] has been one promising upscaling avenue in cases where a separation of scales applies (see also section 9).
• Any concurrent or hybrid scheme (especially those coupling discrete and continuous subdomains) introduces model interfaces that produce spurious physics in terms of, e.g. wave reflection and refraction, or spurious forces. These have often been dealt with in an ad hoc fashion. Therefore, rigorous approaches to reducing the impact of model interfaces are needed.

It is difficult to list specific required advances as the field of concurrent and hierarchical multiscale modeling is broad and, as described above, faces diverse challenges that are both theoretical and computational (see [87] for a review). Dynamic homogenization, intelligent order reduction techniques, as well as discrete-to-continuum coupling scenarios call for new theory and numerical schemes. A practical challenge is the availability of numerical codes, which are primarily developed by and restricted to individual research groups. Open-source multiscale codes are still the exception. In addition, available codes have often grown over decades with consequences for programming languages, specific computing architectures, etc so that their suitability for and performance on current computing architectures is limited (see section 15). High-performance codes for real-time 3D simulations are still a rare find. Combining several distinct codes is oftentimes problematic and requires application programming interfaces. The scalability of such heterogeneous numerical codes is a further computational challenge. A natural challenge in both hierarchical and concurrent scale-bridging schemes is the multi-physics nature of the problems at hand when going across wide ranges of length and time scales—oftentimes combining concepts from, among others, mechanics, chemistry, materials science, statistical physics, and computational science. Current and future solutions must effectively integrate multiple concepts and thus call for interdisciplinary developments, also leveraging advances in data-based approaches (big-data or data-driven modeling), where high-fidelity multiscale models can serve as data generators in lieu of costly experiments.

The recent advent of additive manufacturing techniques for multiscale and multi-material architectures has enabled the creation of mechanical metamaterials with features sizes across many length scales—ultimately dissolving the distinction between solids and structures (see section 13). Especially this growing area requires theoretical techniques and computational tools to predict, optimize, and reverse-engineer the effective metamaterial properties on demand. Especially the inverse problem of optimizing (meta-)materials architectures across scales towards a tailored, superior effective performance is an active area of research [88–90], and rigorous inverse design methods, especially for the complex nonlinear, dynamic, and inelastic macroscale material performance, are needed in conjunction with powerful multiscale techniques.

With only few exceptions, hierarchical and especially concurrent multiscale modeling techniques have not yet made their way into industrial settings and have remained of academic interest, primarily because of the associated computational expenses. Overcoming this disconnect is a further—not purely technical—challenge with great potential benefits for both sides. This is where highly-efficient computational techniques that make use of both intelligent multiscale modeling theory and fully exploiting available computing architectures (as well as suitable verification, validation and UQ protocols) will be game-changing.

Modeling the behavior of materials ideally encompasses several aspects at the same time: (1) It gives a description as compact as possible that represents a sufficiently large set of experimental observations. (2) The model is formulated in terms of physically intuitive quantities. (3) The model establishes links between the small-scale building blocks from which the material is made and the experimental observations, thereby potentially enabling tailored materials design. (4) The model has the potential for extrapolation, i.e. predictions.

Multiscale modeling has emerged as a route for addressing the above issues, because it takes into account different levels of detail, e.g. continuum as well as the discrete constituents of the material, and couples these levels either sequentially or concurrently (see section 8). Another approach for describing complex material behavior is coarse graining (see figure 8). Here, the goal is to start with a fine level of description where one has good physical understanding, and then to transition to a coarse level with less degrees of freedom. While the fine-level degrees of freedom are eliminated, the coarse-graining procedure aims at representing the most relevant physics appropriately on the coarse level. In general, coarse graining refers to both spatial and temporal scales. In this contribution, the temporal aspects are emphasized, to raise awareness and to underline its significance for bottom-up approaches for materials design.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The coarse-graining procedure addressed in this contribution has its roots in statistical mechanics and non-equilibrium thermodynamics [91]. Equilibrium statistical mechanics is a beautiful and well-known realization of a coarse-graining technique, and it is the basis of many powerful computational tools. However, if one shifts the focus to beyond equilibrium, time scales and lack of full equilibration begin to play a crucial role. A particularly interesting question in this respect is how irreversibility on the coarse level (e.g. in the form of a viscosity or diffusion coefficient) emerges from the reversible Hamiltonian mechanics of the atoms. This emergence of irreversibility has been established for close-to-equilibrium situations in terms of so-called fluctuation-dissipation relations [91–93]: fluctuations on the fine level (e.g. in the particle-based shear stress) gives rise to the coarse-level transport coefficients (e.g. shear viscosity). In the last two decades, these relations have been generalized to hold not only close to equilibrium [91, 94, 95], by a combination of statistical mechanics and projection-operator techniques. This serves as a basis for powerful simulations that are suitable for beyond-equilibrium simulations [96]. However, in order to take most benefit of these advances for the computer-aided molecular design of materials, some points need further consideration.

For coarse graining, the choice of variables on the coarse level is of great importance for two reasons. First, as can be seen in classical statistical mechanics, it affects the static properties, i.e. the free energy. While the energy of the system is identical on the fine and coarse levels if the coarse-level variables are sufficiently rich, different sets of coarse-level variables will give rise to different entropies; the less detailed the coarse level, the larger the entropy. This is because the entropy is about the comparison between two levels by, roughly speaking, counting the number of fine-level states that are compatible with the given coarse-level state. However, at least equally important for the non-equilibrium coarse-graining is the fact that by lumping fine-level detail into coarse-level quantities, the latter usually have slower time-evolution. The cases usually dealt with are such that there is a so-called separation of time-scales between the two levels. This implies that the rapid motions on the fine-level are noticeable on the coarse level as stochastic only, which enters the coarse-level model as transport coefficients, as a realization of the fluctuation-dissipation theorem.

The current challenge consists in choosing the coarse-level variables appropriately, implying that the neglected fine-level dynamics indeed occurs on much shorter times, i.e. that it can be considered as stochastic on the coarse level. Current practice is to choose coarse-level variables, execute the coarse-graining procedure, and check in hindsight whether the separation of time scales is indeed respected. This procedure has been followed on the example of dislocation-based plasticity, the fine and coarse levels being those of many discrete dislocations and of average densities of dislocations [97], respectively. The results indicate that there is no separation of time-scales. In another case, namely polymer melts, it was found that time-scale separation depends subtly on the choice of the coarse-level variable [96].

The following questions in relation to the time-scale separation should be addressed: Can for any system a set of coarse-level variables be found such that time-scale separation holds? On the one hand, can one formulate guidelines, or even a constructive (possibly iterative) procedure to work towards such a set of variables? On the other hand, if such a set of coarse-level variables that respects the time-scale separation does not exist or is not needed, what is a proper systematic coarse-graining procedure?

Coarse graining is goal-oriented and often subjective, in the following sense. Goal-oriented in the sense that, along the lines of the first paragraph, the model should be as simple as possible for describing the situations of interest. However, this requirement does not determine the level and variables of the model uniquely, but rather these may be chosen by the modeler, which adds the subjective point. While I believe that insight and intuition are indispensable in coarse graining, the result might still not be successful, and it could be useful to get support from other methods. It is not completely clear what those could be, but some possible directions for further development are given in the following.

One should carefully re-examine the essential steps in the derivation of the current coarse-graining procedure. For example, recently, a generalization has been proposed from diffusion processes to Markov processes that satisfy detailed balance and a large-deviation principle [98]. Further generalization is sought for cases where time-scale separation is lacking.

Concepts developed for pattern recognition could be useful for coarse graining. While humans are perceiving their environment by recognizing patterns continuously, it remains to be examined how pattern-recognition techniques can be employed to improve, in practice, the above mentioned coarse-graining procedure. Computer algorithms will naturally be part of this endeavor, in terms of machine-learning algorithms for pattern recognition [99]. It should be emphasized, however, that the purpose of computer simulations in this context is to support advancement of the theory; one often learns and gets insight from examples and applications. The idea cannot and should not be that large-scale computer simulations on the fine level combined with ML entirely replace physics-based coarse-graining. The latter remains an indispensable conceptual step, not only to foster fundamental understanding, but also to help tailor the design of materials with improved properties.

This contribution is concerned with coarse graining, in which a coarse-level model is derived based on fine-level information, by way of statistical mechanics and projection-operator techniques. In this procedure, the separation of time-scales is essential. It is advocated that the interplay between the choice of coarse-level variables on the one hand and the time-scale separation on the other hand be examined further. Specifically, guidelines or even constructive (possibly iterative) procedures are sought to work towards making good choices of variables in order to alleviate the violation of time-scale separation as much as possible. Furthermore, strategies for dealing with systems that lack time-scale separation are desirable.

The mechanical properties of amorphous metals, polymers, and oxide glasses exhibit a remarkable degree of universality. Through extensive particle scale modeling with either atomistic force fields or generic models for glass formers (Lennard-Jones mixtures and bead-spring polymers), it is now widely accepted that these materials exhibit heterogeneous elasticity at the nanoscale and share a common elementary mechanism of plastic deformation: swift localized irreversible rearrangements involving small (10–100) groups of particles occur in the background of a homogeneously deforming medium (see figure 9) [100]. These 'shear transformations' (ST) occur preferentially at predisposed 'soft' locations that can be associated with local structural features of the glass, such as the coordination of particles in the first neighbor shell or low shear moduli. The elastic response of the surrounding medium to such an ST exhibits anisotropic (quadrupolar) long-range behavior that is well described by the Eshelby inclusion model. Quantitative differences between materials with different chemical interactions are reflected in different types of disorder, size and energy scales of STs, and relaxation rates [101].

In a constant strain rate deformation, amorphous metals exhibit a maximum stress at strains of a few percent. This yield stress is therefore much larger than in crystalline materials. However, the stress maximum is frequently followed by a discontinuous transition in which shear localizes along shear bands, leading to brittle behavior. The tendency towards localization is controlled by the initial degree of annealing. Glasses that avoid fracture and flow, do this in a very intermittent fashion. Their stress-strain curves consist of elastic branches punctuated by sharp stress drops that signal individual or collective plastic events. All glasses undergo structural recovery (physical aging) at elevated temperatures, during which structural relaxation times increase. Plastic deformation rejuvenates the glass, and the molecular mobility increases by several orders of magnitude. This effect is best documented in polymer glasses and can be understood from a combination of repositioning and tilting of the underlying potential energy landscape. Post-yield deformation of polymer glasses differs from metallic and oxide glasses in one important aspect as polymers can exhibit strain hardening due to their macromolecular character. Glassy strain hardening originates not from entropic elasticity as in rubbers and elastomers, but rather from dissipative work (i.e. breaking of intermolecular van der Waals bonds) performed by polymer chains forced to deform affinely in the glassy matrix.

Mesoscale elastoplastic models (EPM) were introduced about 25 years ago [102], and coarse-grain amorphous solids into a lattice of blocks whose size corresponds to a typical ST. The blocks are endowed with (visco)elastic properties and a local yield threshold: upon yielding, they redistribute their load non-locally and anisotropically throughout the medium with an elastic propagator of the Eshelby form, G(r, θ) ∼ cos(4θ)/r^d (in the plane of the ST), with d the dimension of the system [103]. Yielding is either treated as an activated event (thermal glasses) or occurs upon crossing a threshold (athermal dynamics). Such cellular automaton models reproduce not only the bulk rheology of amorphous solids, but also enable the study of mechanisms of strain localization and the critical behavior at the yielding transition, namely the statistics of avalanches (distribution of stress drops) that dissipate energy collectively. Both reciprocal space and finite element methods (FEM) have been constructed to solve numerically for the long ranged elastic interactions in these models. This level of modeling is akin to a treatment of CP at the level of explicit dislocation dynamics.

For simulations at the atomistic scale, a persistent challenge (common for all simulations of condensed matter) consists of reaching realistic quench and deformation rates. On the timescales accessible with conventional MD, the prepared glasses are positioned too high on the energy landscape, and some thermal relaxation processes will therefore be pre-empted at high strain rates. Novel swap MC methods alleviate the former issue, but they are not yet generically applicable to all glass formers. While much has been learned about the nucleation, form, elastic consequences and correlations of STs with structural features of the amorphous packing, a point that still requires clarification is their relationship with the related concept of shear transformation zones, which considers STs as pre-existing defects in the structure that can be polarized, i.e. they are endowed with an orientational degree of freedom that reflects the local orientation of a group of atoms [104]. Gaps also remain in our understanding of the statistical properties of the mechanical noise produced by plastic activity and its role in activating plastic events. While particle simulations can capture the formation of transient slip bands through accumulations of multiple STs, it remains challenging to push them into a regime where permanent shear bands, their mechanisms of formation and dependence on parameters such as temperature or the degree of initial disorder can be studied.

EPMs formulated to date are believed to capture all generic qualitative features of amorphous plasticity, but they are not yet able to make quantitative predictions for specific materials. To this end, they must be systematically informed from atomistic simulations. This requires measurements of statistical properties such as the disorder distribution of local dynamical moduli, and also nonlinear quantities such as local yield thresholds and their spatial correlations. The phenomenological rules governing transitions between elastic and plastic states must also be better validated, possibly by connecting them to the energy landscape of a mesoscopic region. Lastly, EPMs must realistically model the dynamical propagation of yield events instead of imposing instantaneous elastic interactions and address the time scales associated with thermal activation. In the athermal limit, considerable debate still revolves around the universality class of the yielding transition, the nature of the associated control and order parameters, and different mechanisms contributing to transient or permanent shear localization and their dependency on control and material parameters. EPMs are also needed when spatially heterogeneous flow patterns arise or when spatial correlations are of interest.

The MSMSE community has expended considerable effort on accelerated MD techniques (parallel replicas, temperature acceleration, metadynamics; see also section 7) [105] and to efficiently explore the energy landscape of solids at zero temperature (activation relaxation technique). Unfortunately, most of these methods do not work well in the extremely rugged and highly variable energy landscape found in glasses. A deeper understanding is needed on the conditions and regimes of applicability of such methods.

Progress could be facilitated by collecting information about various models and analysis tools in a common repository that can be shared between different research groups. This could include characterization of the differences between various interaction potentials (flavors of Lennard-Jones mixtures, EAM potentials for metallic glasses and force fields for oxide glasses), diagnostics for detecting and characterizing plastic activity, and the extraction of locally heterogeneous properties like elastic moduli and yield thresholds.

EPMs similarly have been constructed in slightly different implementations by different research groups. It would be helpful if a common standard were developed and easily shared. While physicists have focused mostly on aspects of critical behavior in the athermal limit, materials scientists have emphasized viscoplastic effects at an operating temperature much closer to the glass transition. The next generation of EPMs should unify these regimes through the use of efficient dynamical FE codes to solve the Cauchy momentum equations in two and three-dimensions. This could benefit from synergy with the well-established field of dynamical FEM modeling in solid mechanics and by establishing standard simulation codes similar to LAMMPS for MD.

The ability to relate the rheological or plastic flow properties of amorphous materials to a microscopic description that can be treated with statistical methods would constitute major progress in materials physics. Theory and modeling of amorphous plasticity has been hampered by the absence of easily identifiable structural signatures of the elementary carriers of plasticity. On the other hand, amorphous materials do not exhibit texture that introduces an additional length scale in CP. It should therefore be easier to directly link mesoscopic and constitutive continuum models that operate on a mean field level.

Electronic, electrochemical and photovoltaic devices are systems of interfaces between different materials. Transistors, memory cells, capacitors, photovoltaic and fuel cells produce and conduct electrical current in heterogeneous environments of several materials. Due to constant scaling, the component materials of these devices are often nm thick films and their properties are strongly affected by interfaces with other materials. A typical example is shown in figure 10 [106]. The electronic properties of thin oxide films, such as one in figure 10, vary widely with their morphology and with the morphology of the interface. The film morphology, in turn, is very sensitive to small variations in chemical composition and to deposition and post deposition processing [106]. Most as-deposited thin oxide, nitride or chalcogenide films are amorphous and annealing at high temperature leads to their partial or full crystallization and formation of polycrystalline films. Structural disorder in amorphous films and grain boundaries in polycrystalline films help to reduce strain developing at the interface due to misfit of structural parameters of interfacing materials (see figure 10).

Zoom In Zoom Out Reset image size

Materials modeling has been instrumental in understanding the microstructure and defect properties of heterogeneous microstructures in ceramics, Si/oxide and metal/insulator interfaces [107]. The earlier works on metal/oxide interfaces emphasized the role of image interactions, which has later been included in atomistic simulations. More recent work is using DFT and MD to predict atomistic interface structures. There are broadly three approaches employed for constructing semiconductor/oxide and metal/insulator interfaces. All of them employ a periodic model along the interface plane and therefore are subject to constraints when trying to accommodate two or more materials in a periodic cell of certain dimension. Some simulations combine lattice matching with annealing using MD. This mimics a growth process and gives additional freedom for the interface atoms to relax. The most popular approach to constructing metal/oxide interfaces is using static DFT calculations for systems constructed by lattice matching [108] for different values of oxygen pressure [109]. The main criterion of the interface stability is the interface free energy as a function of oxygen pressure [110]. Yet another approach is based on simulating the film growth via atomic layer deposition or chemical vapor deposition using kinetic MC and other techniques [111].

Although atomistic modeling to date had established several important structure-property dependencies, the main challenge remains in identification and characterization of the realistic interface models as opposed to simplified (and possibly even unstable) structures and characterization of interface defects. Electrically active interface defect states may lead to degradation of device characteristics, such as the threshold voltage, the on-current, or the surface carrier mobility. Creation of an interface naturally leads to breaking some of the bonds in the interfacing materials. To further improve the interface quality, the number of dangling bonds can be further reduced by annealing the system in forming gas and passivation with nitrogen or hydrogen atoms.

Different communities use different methods of preparation and propose different ideas regarding the structure of these films. For example, the thermodynamic approach [112] suggests that an amorphous structure of oxide overgrowths on metals may be stabilized with respect to the corresponding crystalline oxide, up to a certain critical oxide-film thickness. This results from the reduction of mismatch strain for the amorphous oxide film (in contrast to an epitaxial or semi-coherent crystalline oxide film). The value of critical oxide-film thickness when initial amorphous oxide film transforms into a crystalline oxide film depends on the substrate orientation, temperature and metal–oxygen system under study [112]. However, density, thickness and morphology of thin films strongly depend on the substrate and method of deposition, anneal temperature, as well as other factors, such as oxygen and water pressure, which makes simulation of these structures challenging.

The biggest challenges for modeling defects in non-crystalline films are related to the variations in their local environment due to disorder and position with respect to interfaces as well as accounting for effects of interface strain and space charge. Using periodic cells is preferred in order to avoid border effects which may affect defect characteristics. However, this makes the amorphous structures quasi-periodic and induces constraints on the structural relaxation accompanying defect creation. To credibly predict distributions of properties, calculations at many sites in amorphous structure and of many models are required. If modeling amorphous films and interfaces is difficult, partial or complete crystallization of as-grown amorphous films and formation of grain boundaries and dislocations provides further challenges. Yet, there are no robust, comprehensive computational tools to predict and design how physical and functional properties of heterogeneous systems can be controlled by their microstructure.

Advances in applications of high resolution TEM, scanning probes and spectroscopies to studying the atomistic structures of hetero-structures, grain boundaries in ceramics as well as in determination of the 3D shape of nanoscale crystals with atomic resolution, make observations of atomistic structures and defect chemistry increasingly possible. Despite this, atomistic modeling of nanostructures is still dominated by the use of periodic models to study bi-crystalline arrangements of grains and adhesion between two nanocrystals. There are attempts to bridge nano- and microscale simulations through coarse graining methods [113]. However, much microstructure modeling is still based on Voronoi polyhedra, although some 3D FE modeling of poly-crystalline ceramics has been done using more realistic microstructures (e.g. [114, 115]). Some of the existing methods provide deep insight into the microstructure of materials but are rarely linked with their function. To improve design and performance of existing devices and develop future devices requires developing novel multiscale modeling tools combined with detailed experimental studies that can establish how materials properties depend on microstructure and chemical composition and how they change under electrical stress, heat and humidity. This requires the following.

• (1)  
  Deriving atomic-scale and coarse-grained models and databases of structures and properties of nanostructures, grain boundaries and metal/insulator interfaces, characterizing them by stability, structural topology, composition and electronic structure. This can be achieved using statistical approaches (informed by both simulation and experiment) to limit the search space to those interfaces which occur most frequently in polycrystalline materials and at contacts with electrodes.
• (2)  
  Developing new methodologies to simulate at the atomic level the response of individual structural features (interfaces and grain boundaries) and of complex systems to external fields, thermal gradients and mechanical constraints, and provide predictions and interpretation of experimental data.
• (3)  
  Bringing together atomic level detail for materials properties and processes in multilayer systems and validating multi-physics methods using 3D microstructural simulations (using, e.g. FE methods) to predict the response of full microstructures to electric currents, heat flow and mechanical constraints for electronic, ferro-, thermo- and piezo-electric ceramic materials and multifunctional structures.
• (4)  
  Investigating the role of dopants, impurities and vacancies in determining the material's function using e.g. grand canonical MC simulations at various fixed and uniform chemical potentials of O₂ (oxygen pressure), H₂0 (humidity), impurity concentrations to sample the structures and compensation mechanisms resulting from growth in variable oxygen/water pressure.

Most materials in electronic devices are either amorphous or made up of a complex assembly of small crystals (grains) which, together with the boundaries joining them, constitute the microstructure. Furthermore, to perform in electronic devices materials need electrical contacts; those are metallic or semiconducting and have their own microstructure. Material function is therefore multi-scale and is intimately linked to microstructural processes (diffusion, electron transfer, chemical reactions) at the internal grain boundaries and interfaces to electrodes or contacts. The roles of microstructure (particularly texture) are well established in metallurgy, but much less understood in functional ceramics, particularly when it comes to linking the atomic structure of polycrystalline materials to their function. The challenge is to merge this information at various length scales and to develop a fundamental understanding of the evolution of the microstructure and its effect on properties, such as thermal-mechanical response, electrical response, degradation, and failure. A step-change can be achieved only when it is possible to predict the effect of changing the microstructure on system function and so guide the materials processing route.

The application of fiber-reinforced composites in structural applications has grown continuously during the last forty years owing to their outstanding stiffness and strength as well as low density in addition to a reasonable damage tolerance. They have become the standard materials for high add value structural components in aerospace, sports, automotive, petrochemical, naval architecture, etc. From the viewpoint of manufacturing and design, composite materials present two singular features. Firstly, manufacturing of the material and of the component is carried out simultaneously, opening the possibility to optimize the properties of the material for each specific application. Secondly, fiber-reinforced composites are highly anisotropic and present a hierarchical structure with different failure mechanisms which depend on the length scale and the loading mode. Thus, accurate prediction of the strength of composite structures using standard simulation tools was not possible, leading to a very costly pyramid of experiments (beginning with composite plies and ending with full components) to qualify new materials and certify components.

These limitations were overcome in recent years with the development of novel multiscale modeling strategies [116, 117]. They are based on the consideration of three different length scales to account for the deformation and damage mechanisms in the composite material. The first one is the ply level, at which the control length scale is the fiber diameter (≈10 μm). The homogenized properties of the ply (stiffness, strength and toughness) can be obtained by means of computational micromechanics from the properties, volume fraction and spatial distribution of the matrix and fibers in the ply. This information is passed to the second level, which considers a multidirectional laminate in which the dominant length scale is the ply thickness (100–300 μm). The mechanical properties of the laminate are obtained by means of mesoscale simulations of the laminate based on the homogenized properties of the plies and on the interply behavior. This information is used as input to simulate the properties of the component made up of laminates by means of computational mechanics simulations using shell elements for each laminate. The dominant length scale in these simulations is the laminate thickness (>2 mm) and the homogenized properties of the laminates are provided by the mesoscale simulations.

These multiscale strategies can accurately determine the strength of components manufactured with composite materials and provide guidelines to optimize the design including the effect of matrix and fiber properties as well as of fiber architecture. Moreover, they lead to a significant reduction in the number of costly mechanical tests for certification and qualification of new materials.

Further expansion of fiber-reinforced composites to other industrial sectors and applications is limited by two factors. The first one is the high and recurrent manufacturing cost and, hence, methods to reduce processing costs (while mechanical properties are maintained) are mandatory. The second one is the limited functional properties, as compared with metallic materials. The use of fiber-reinforced composites in structural applications will be fostered if the thermal and electrical conductivity, fire resistance, etc, are equivalent to those of metals and new functionalities (energy harvesting and storage) are incorporated [118].

The first challenge is being addressed through the development of virtual processing tools. They are the perfect complement to the virtual testing tools presented above and can avoid or reduce expensive trial campaigns enabling right-first-time concepts in manufacturing (figure 11). Simulation of fiber-reinforced composites processing is a complex multiphysics problem that involves fluid mechanics (infiltration of the resin into the porous fiber preform), chemistry (resin curing), heat conduction and solid mechanics (fabric forming). Moreover, these problems have to be solved at different length scales (resin infiltration at micro and meso levels) and the current models and strategies available are still far from optimum [118]. There is a lack of robust computational methods for transferring information between the different length scales (micro-to-macro) as well as of fast and efficient strategies for coupling the different physical processes. Moreover, there is also a lack of reliable experimental methods to quantify uncertainties in processing parameters, and to a larger extent, to track uncontrolled manufacturing disturbances. The complexity of the simulation techniques impose large demands in terms of computer time and prevent their direct application to perform on-the-fly simulations during the manufacturing process. In this regard, surrogate models and strategies based on proper generalized decomposition and dimensionality reduction using principal component analysis become very relevant to speed-up manufacturing simulations [119–121].

Zoom In Zoom Out Reset image size

Current trends in manufacturing include the use of sensors to capture processing disturbances. They generate huge amounts of data that can be used to trigger control and corrective actions through the application of data analytics strategies based on experience and knowledge obtained during years in factories. Thus, deep learning and artificial intelligence algorithms can be used to solve the inverse problem to determine on-the-fly the actions to minimize defects during processing [122, 123]. Then, virtual translation/copies of a manufacturing process will become possible enabling the implementation of smart manufacturing concepts. As a result, processing will be controlled by expert systems that take decentralized and automated decisions leading to significant reductions in manufacturing costs by decreasing rejection, inspection, and/or repairing issues.

Extension of the virtual testing techniques to include functional properties requires coupling the different physical problems into a multiphysics framework. For instance, inclusion of carbon nanotubes or grapheme sheets into the polymer matrix improves the electrical conductivity. Typical electro-mechanical models use a network of standard resistive element (following Ohm's law) to account for the electron flow through the carbon fillers while the network ends are connected by special elements simulating electron tunneling between neighbor carbon nanotubes [124]. Such network can also be connected to mechanical models to simulate the effect of carbon fillers on the mechanical polymer matrix, which in turn modifies the properties of the composite [125].

Despite the progress achieved in recent years, substantial research is needed to expand the current simulation strategies and to provide more efficient tools that can be used by industry. These efforts should be centered on improvements in the multiscale and multiphysics approaches, as described below.

Extension of bottom-up multiscale methods becomes of primary importance when physically-based models are being demanded. They include atomistic simulations based on first-principles or molecular mechanics to determine the properties of the polymer matrix, fibers and interfaces. This information, combined with robust and reliable methods for mimicking real micro and meso-structures using statistical information obtained from 3D imaging methods (x-ray microtomography), can be used within the framework of computational homogenization to obtain the multifunctional properties of fiber-reinforced composite plies.

Multiscale simulation of composites is based on the analytical or computational homogenization of the behavior of a RVE of the microstructure and it is supported by the clear separation of length scales between micro, meso and macro levels (see section 8). This strategy is very efficient and accurate to determine average properties (elastic constants, thermal conductivity, etc) but it is not so efficient to predict behaviors that are controlled by the local details of the microstructure (fracture localization, electrical percolation). Thus, extension of homogenization theories that can account for localization effects are needed [126].

Multifunctional composites with high thermal and electrical conductivity are starting to be used in multifunctional applications (such as protection against lightning impact in aircrafts). The macroscopic simulation of these phenomena requires efficient multiphysics codes that couple electrical and heat transport equations with plasma physics and mechanical and acoustic interactions [127, 128] (figure 12). The integration of any computational tool into non-academic codes is necessary in order to consolidate and expand this knowledge to design and predict the properties of multifunctional composites in industrial environments.

Zoom In Zoom Out Reset image size

Multiscale modeling of the mechanical properties fiber-reinforced composites in the last decade has led to the development of accurate virtual testing tools. They were based on the clear separation of length scales—at the micro, meso and macro level—and multiscale modeling was carried out by successive homogenization and transfer of information between these length scales. Current challenges in the multiscale modeling of multifunctional composites are focused in the extension of the strategies to simulate processing and to account for multifunctional properties (thermal and electrical conductivity). Both scenarios require the use of efficient and robust multiphysics modeling tools to couple the different phenomena involved in the processing and the performance of multifunctional composites. Last, but not least, transfer of these strategies and tools to industry is necessary to support the expansion of composites in engineering applications.

Metamaterials are man-made materials, revealing remarkable effective properties that are considered as unconventional or counter-intuitive, i.e. going beyond the properties of their constituents. Whereas metamaterials have a long-standing history in, e.g. electro-magnetics and optics, their impact in mechanics emerged more recently. In the field of mechanical engineering, metamaterials are typically multi-scale architectured materials, whereby the structure and geometry across the scales govern the unconventional effective macro-scale engineering properties. They typically reveal (correlated) fine-scale fluctuation fields, that govern the macroscopic response beyond their mean value. The field is growing rapidly and has expanded its scope towards a wide range of static and dynamical properties. Metamaterials exploited in a quasi-static or transient regime (not focusing on wave propagation related properties) are here called mechanical metamaterials, whereas those affecting mechanical wave propagation are termed dynamical metamaterials.

The field of mechanical metamaterials gradually evolved from simple (auxetic) counter-intuitive mechanisms (30 years ago) to the considerably more complex mechanistic responses achieved today. Recent examples show how to convert compression into twist [129], or how to convert progressive instabilities into motion [130]. The modus operandi of these materials is to switch between different microstructural states under the influence of internal constraints or an externally applied load. Accordingly, the activation of these switching mechanisms can be passive or active. The mechanisms are either built-in or induced by material or geometrical instabilities in the microstructures [131, 132], leading to the formation of microstructural patterns, see figure 13 for a simple example. The complexity pops up when such microstructures pattern differently under the influence of internal frustration [134], different external loads or boundary conditions. Integrating advanced actuation mechanisms to externally control the desired mechanistic response, gives rise to a smart mechanical metamaterial serving various future applications, e.g. for soft robotics.

Zoom In Zoom Out Reset image size

The field of dynamical metamaterials focuses on the mitigation or manipulation of mechanical waves in substructured solid or fluid-solid systems [135–138]. The modus operandi is here twofold: Bragg scattering (i.e. phononic crystals) and local resonance (i.e. locally resonant metamaterials). The first mechanism typically operates at length scales of the same order as the size of the underlying periodic unit cells (meta-atoms), leading to the formation of band gaps. The local resonance mechanism on the other hand, operates in the deep-subwave length regime, interacting with microstructures that are much smaller than the wave length of the excitation wave. It reveals sharp localized band gaps, driven by negative (effective) mass and/or stiffness effects. Research in this field has given clear directions towards applications in sound attenuation, high-resolution acoustic imaging, transformation acoustics, cloaking, acoustic lenses, non-reciprocal transmission, vibration control, instabilities-induced large-amplitude wave propagation in dissipative media, mechanical analogs of topological insulators, etc [132, 138].

Even though the physical modi operandi of mechanical and dynamical metamaterials are well understood, the design space—spanning microstructural to engineering scales—is extremely wide, and various fine scale material properties influence the effective macro-scale response. The challenges reside both in the design, synthesis, manufacturing, testing and modeling of these materials. Here, focus is given on the challenges requiring advanced modeling approaches.

For mechanical metamaterials, the main research challenges are:

• Multi-stable patterning metamaterials: multi-stability, entailing non-convexity, should be further explored for numerous microstructural configurations. This leads to the controlled appearance of multiple complex mechanistic fine-scale patterns, enabling mode conversion (e.g. from tension to torsion, from compression to shear, etc).
• Effective continuum models: for engineering applications, it is computationally impossible to incorporate the microstructure in full-scale macro-scale computations. Predictive engineering analyses can only be achieved through appropriate homogenized constitutive models that take the effect of the microstructure into account. Scale separation becomes non-trivial in most cases, with non-trivial boundary conditions, affecting the emerging patterns.
• Actuation/activation: the fine scale exploitation of mechanical, electric, magnetic, thermal, or chemical fields to actively switch between different patterned microstructural states has to be scaled up to the effective emerging properties at the coarse scale.
• Load and power transmission: conventional smart materials are attractive in their kinematical actions, but are very limited in terms of the power they can exert. Metamaterials have the ability to deliver a relatively high output power, which makes them challenging and rewarding.
• Nonlinear and time-dependent phenomena: nonlinearities have a pronounced influence on the macro-scale response. Viscous damping and rate effects may strongly affect the mechanistic patterns, which needs an in-depth modeling analysis.

For dynamical metamaterials, important research challenges are:

• Effective continuum models: homogenized effective models are needed to enable real engineering applications, covering both Bragg scattering and local resonance [139]. Developing these in nonlinear regimes presents an even larger challenge.
• Modeling approaches for non-periodic and dynamic metamaterials with distributed heterogeneous resonators: to achieve wider frequency bands for attenuation or to exploit the potential of local resonance phenomena in an aperiodic microstructure.
• Multi-scale design and optimization: multi-scale optimization and reversed engineering (from macro to micro) offer great potential in identifying novel promising multi-material designs.
• Nonlinear material effects: nonlinearity makes the microstructural response amplitude-dependent and the overall response intrinsically size-dependent. Nonlinear aspects also impact mode-conversion and wave-wave interactions.
• Geometrical nonlinearities and instabilities: topological microstructural changes induced by large deformations or instabilities strongly alter the response of dynamical metamaterials [132].
• Dissipation: material and viscous damping are known to affect dynamical metamaterials to a large extent. The combination with thermo-viscous dissipation in entrapped fluids, gives rise to metafoams, holding great potential. Non-Hermitian acoustics aims to either dissipate or add energy to the system [136].
• Activation: locally deforming microstructures at the small scale in a controlled manner enables reconfigurability, adaptivity and tunability of the effective behavior at the coarse scale, e.g. magnetoactive acoustic metamaterials [140].

To address the challenges sketched above, a multi-scale approach is needed with expertise from different fields. Common to both mechanical and dynamical metamaterials, the following required advances in modeling can be identified:

• Multi-scale computational design and optimization: the microstructural degrees of freedom in these metamaterials span both the 3D micro-scale geometry and the integration of (different) materials. Dealing with nonlinearities, instabilities, dissipation, actuation, requires a systematic modeling-driven approach. Targeting particular macro-scale kinematics or dynamics necessitates a modeling and (topology and materials) optimization approach [141] that links the scales towards the desired engineering properties, ultimately leading to the reversed engineering of such metamaterials.
• Advanced homogenization approaches: the homogenization of metamaterials does not fit in classical first-order homogenization schemes, since scale separation does not truly apply anymore. First attempts to homogenize the response of these materials clearly reveal the emergence of a so-called micromorphic continuum with additional physical fine scale degrees of freedom. Moreover, the identification of proper boundary conditions at the macro-scale still remains difficult and needs a dedicated multi-scale approach.
• Small scale defects and heterogeneities: mechanical metamaterials are generally modeled in a strongly idealized setting. Processing induced defects, both in the materials and the geometry, and spatial heterogeneities need to be integrated in multi-scale modeling approaches.
• Actuation and activation: physical actuation mechanisms, their integration in the microstructure, and the control thereof, leading to different micromorphical states and accompanying macro-scale response, still require extensive research.
• Nonlinearities, anisotropy and time-dependence: the next generation mechanical/dynamical metamaterials is expected to systematically exploit nonlinearities, anisotropy and/or time-dependent effects. Nonlinear effects are hardly explored and may have a strong influence to the benefit of the overall response.

A major recent advancement in the field of mechanical and dynamical metamaterials relates to their manufacturability. Indeed, the progress made in additive manufacturing constitutes a real enabler for this class of materials, but it also requires extensive parallel modeling to achieve real predictive models.

The field of this class of metamaterials is expanding tremendously. The future of acoustic metamaterials looks sound [138], and mechanical metamaterials are expected to make a break-through, among others, in soft robotics applications. Numerous exciting effects have been demonstrated in the literature and the prospects are great. Yet, the path towards their real engineering exploitation is long and difficult. Metamaterials are intrinsically multi-scale in nature. The effective properties are due to the underlying physical, geometrical and topological details at the fine scale. Advanced modeling is required at both scales, and reliable predictive scale bridging methods are indispensable.

The market for mechanical or dynamical metamaterials is still in a premature state. In order to enable the many envisioned applications, computational efforts will need to be complemented by parallel progress in manufacturing (bulk processing and 3D additive manufacturing), multi-scale experiments and model validation.

In this section we will explore where we stand with regards to multiscale modeling of any material by discussing the immensely important material of steel. A simplified, partial, phase diagram of this iron/carbon alloy is shown in figure 14. We will focus on the task of describing the process of cooling a piece of steel in the austenite phase (carbon dissolved in γ-iron) to the room temperature steel phase. This is one of the most important processing steps to form steels of varying properties.

Zoom In Zoom Out Reset image size

The phase diagram in figure 14 contains two solid phases of iron, the body-centered cubic (bcc, α, ferrite), and the face-centered cubic (fcc, γ, austenite). For simplicity, high temperature solid phases and liquid will not be discussed here. Carbon exists in two phases, atomistic when dissolved in fcc-iron, and combined with three iron per carbon in the Fe₃C compound (cementite). When γ-iron with carbon in solid solution is cooled from the austenite phase, all carbon will seek to form Fe₃C and the remaining Fe will form a bcc-lattice, that is, α-iron. This is because α-iron cannot contain more than a maximum of around 0.1 atom% carbon. Depending on cooling rates a wide variety of microstructures can form in this process. Pearlite, a eutectoid (two-component) phase, at around 3.5 atom% carbon, with Fe₃C and α-iron forming a laminar structure, is the thermodynamically stable phase. Iron with higher than 9 atom% carbon is called cast iron. In addition to carbon, a variety of low atom% substitutional elements (replacing iron ions in the lattice), such as, manganese (Mn) and chromium (Cr), are used to fine tune the properties of the steel.

For proper modeling of the subset of processes important for steel production described above, we will need to investigate at least two things: (1) How the phase transition from fcc- to bcc-iron works depending on substitutional elements, carbon content and cooling rate, and (2) how substitutional elements affect the diffusion of interstitial carbon in fcc-iron.

When a grain of iron transforms from an fcc to a bcc lattice structure on cooling, the atoms can rearrange in several different ways and the resulting bcc lattice has a specific orientation vis-a-vis the original fcc lattice. Fast cooling follows the Bain path and the Bain orientation relation (OR) between fcc- (blue/white) and bcc- (red) iron is depicted top, right in figure 15. This fast cooling can trap carbon inside the octahedral space in the fcc lattice, resulting in martensite, a distorted bcc-lattice with an occasional carbon in one of the bonds (red). Slower cooling give preference to several other transformation paths resulting in bcc phases with, among others, Pitsch, NW, and KS ORs, shown in figure 15. Depending on cooling rate, a range of microstructures can be formed between martensite, with bcc-iron in the Bain OR, and the above-mentioned pearlite phase which has bcc-iron mostly in the KS OR [142]. An accurate description of these transition paths is crucial.

Zoom In Zoom Out Reset image size

These transformation paths in pure iron can be examined using DFT [143, 144] based methods ([145], see [146] for start and end points of such transformation paths), at least within the accuracy of the exchange-correlation functionals available [147]. However, how these paths change when substitutional elements are present is not possible to accurately examine with DFT due to the low concentration of these elements (often less than 0.1%). Similarly, diffusion of carbon in fcc-iron can possibly be done with DFT based methods, but how the diffusion is influenced by substitutional elements cannot. The key problem is the unfeasibly large supercell that would be needed to accurately explore low atom% substitutional elements in iron with DFT.

To examine the effects of substitutional elements on phase transition paths and carbon diffusion, the determination of ion–ion interactions by explicit treatment of electrons within DFT will need to be replaced by some computationally simpler methods, using DFT only to verify the accuracy of these larger scale methods. Most common today is the use of parameter fitted formulas to determine forces, so called interatomic potentials or force fields. Force fields can be used both in classical MD and MC methods. However, we are only trading the large supercell problem in DFT, for how to treat many different types of ions influencing each other, and how to treat temperature and magnetism [148].

The simplest interatomic forces are pair potentials and we need one specific potential for each type of pair. For iron/carbon, C–C, C–Fe, and Fe–Fe potentials are needed, and adding one substitutional element, S1, adds S1–S1, S1–C, and S1–Fe pair potentials. The number of pair potentials grows fast with the number of different elements present. In addition, in solid state settings pair potentials are never sufficient, and more complicated force fields, taking many- (more than two) body interactions into account, are needed. So far there are only preliminary works on including temperature [149] and spin dependence [150] into calculations using force fields.

On the next scale up from MD and MC methods, one common technique for investigating the microstructure of steel and other materials is phase field (PF) modeling. In figure 16, a very simple, hypothetical, MD investigation of a phase transition, of use as input for PF modeling, is depicted: What is the velocity of the phase transition between fcc- and bcc-iron, as a function of carbon content, substitutional elements, and cooling rate? With current methods [151] this question cannot be accurately addressed and our inability to computationally answer questions relevant for PF and other meso-scale methods has broken the multiscale chain. In order to be able to take this next step we will need to be able to model at the atomistic level with improved or new methods. A minimum effort would be to develop temperature dependent force fields that can deal with magnetism.

Zoom In Zoom Out Reset image size

The difficulty to accurately construct the plethora of force fields needed in a realistic calculation also compels us to ask, is there a better way? The Hohenberg-Kohn theorem of DFT [143] says that all information about the system is embedded in its electron density. Considering Kohn's near-sightedness principle [152] as well, maybe it is possible to construct a fully local model for how ions are moving within a lattice. Such a model, only dependent on electron density, would be independent of the exact type of neighboring atoms present, eliminating the need for specialized force fields for every type of environment, and enabling needed classical MD/MC simulations.

Even if we simplify the steel problem dramatically we still do not have the methods available for making a full investigation from quantum to continuum with only calculational input. For the foreseeable future we will need to continue doing smaller modeling efforts based on interactions with experiments. This situation is not a prohibiting issue for steel but for investigations of materials that are expensive, dangerous, and/or impossible to synthesize in enough quantities to perform experiments on, it is urgent to advance the science to the point where experiments are only needed for validation of the results of our purely calculational multiscale methods.

Modern multiscale materials modeling increasingly requires (i) powerful simulation engines that take advantage of modern, often heterogeneous, computing platforms, (ii) distributed, collaborative efforts to develop and sustain these simulation tools, (iii) complex workflows that integrate multiple models and data sources, including UQ, to tackle complex problems. In addition, as this Roadmap article has illustrated, there is a growing need to integrate tools from data science and ML to bridge between scales and to orchestrate simulations and workflows in materials design efforts. The complexity of the resulting workflows and the need to make these tools and the data they generate accessible and useful to the broad material community calls for the use of modern cyberinfrastructure to tackle increasingly complex multiscale challenges.

The US National Science Foundation is credited with having coined the term cyberinfrastructure as the collection of high-performance computing resources, networks, simulation software and data, connected experimental facilities and communication software designed to accelerate progress in science and engineering [153]. Around the world, significant investments in cyberinfrastructure, also called e-science [154], are being made seeking to accelerate the advance of science and engineering as well as enhance education and learning. The field of materials science and engineering is capitalizing from this digital revolution as it seeks to synergistically combine experiments, multiscale modeling, and data science to advance the discovery and optimization of new materials, understanding and prediction of materials behavior, and their deployment into new technologies. To significantly accelerate multiscale modeling of materials, new methods must leverage these collective capabilities and capitalize on growing cyber-resources. Here I discuss the available resources, opportunities, and the potential of cyberinfrastructure to transform materials modeling and single and multiple scales.

Recent efforts in cyberinfrastructure build, in part, on advances in materials modeling over many decades and the development of powerful community codes. For example, the development of the local density approximations to the exchange and correlation functional starting in the 1980s, followed by gradient corrections, provided a balance of accuracy and computational efficiency that transformed DFT in the workhorse electronic structure method for materials simulations. Today DFT powers some of the largest materials data repositories in the world such as the Materials Project [155]. Similarly, the development of force fields capable of describing a wide range of materials since the 1980s, see for example [156–161], added MD to the toolset available to materials researchers. Similar progress can be cited at mesoscales, with tools that are being adopted by the wider community [162, 163]. The impact of these models has been significantly enhanced by widely available, powerful and generally applicable codes that enabled researchers to use them as tools without the need to write software from the ground up. Examples of tools widely used across scales in materials modeling include Quantum Espresso [164], LAMMPS [7], ParaDIS [165], OOF2 [166], and MOOSE [167].

More recently, significant efforts have been devoted to making these tools and associated data universally available and useful. The US National Science Foundation's nanoHUB has delivered online simulations since 2002 and currently provides end-to-end services for tool developers to make their codes accessible via cloud computing with powerful and easy to develop graphical user interfaces and access to high-performance computing hardware. Users access these tools from a standard web-browser, without the need to download or install any software. Automatic UQ is available for a wide range of tools [168]. Researchers, educators and students can focus on their domain science and not worry about the technical challenges associated with setting up simulations or accessing leadership-class computing. In addition, incentivized greatly by efforts like the Materials Genome Initiative [169] and the European Materials Modeling Council [170], the current decade witnessed an explosive growth in the deployment of web-accessible data repositories, registries and analysis tools. Prime examples of these resources include the Materials Project [155], OpenKIM [9], NOMAD [171], AFLOW [172] and many others. This growing cyberinfrastructure enables the use of data science tools to uncover patterns, identify materials with desired properties, and discover new ones, see for example [173–175]. In the field of multiscale materials modeling, these simulation and data resources contributed significantly to reproducibility of published results and, more generally, shortened the time required to set up simulations. For example, using openKIM, researchers can utilize published interatomic potentials for MD simulations by simply linking openKIM models to supported MD codes, as opposed to copying parameters from a paper into an input file or writing scripts to create tables as used to be case until recently. Similarly, crystal structures from the Materials Project, pseudopotentials from a library like Quantum Espresso's PSLIBRARY, and online simulation tools [176] make it easier to set up and perform DFT calculations. Last, but certainly not least, AiiDA, a powerful infrastructure to develop, manage and share computational workflows and data [5], is contributing to reproducibility, data reuse and sharing, and insulating the researcher from the details of the underlying compute and storage systems. The field of mesoscale modeling would benefit greatly from the development of such an infrastructure of codes and data in areas like dislocation dynamics and PFs. The creation of the Phase Field Community Hub [177] and online simulation tools [178, 179] are auspicious initial steps. Challenges and opportunities in data management and utilization are not restricted to computational efforts; experimental work has identical needs. Often, a combination of experimental and computational data is needed to solve challenging problems. Examples of repositories of experimental materials data include the Materials Data Facility, developed by the Center for Hierarchical Materials Design (CHiMaD), currently hosting over 35TB of data [180], the NIMS Materials Database from Japan [181], and the Materials Data Repository from the US National Institute for Standards and Technology [182].

The growing availability of web-accessible data and simulation tools could launch a new era in multiscale modeling, one where researchers could create multiscale modeling workflows with distributed data and codes. Simulation codes could automatically request input data from remote servers via web-services. Smart caching of simulation results and surrogate models could be used to provide results in a timely fashion and avoid duplicate calculations. While the technology to perform these tasks exists today, and many cyber-resources are web-accessible [183], multiscale modeling workflows still require time consuming human intervention. This is exemplified graphically for the case of a nanoHUB simulation in the left panel of figure 17. To run and store a DFT calculation, a user may find a set of crystal structures in the Materials Project and pseudopotentials in the PSLIBRARY. She would have to download them to her desktop and then upload them in, for example, the Quantum Espresso tool in nanoHUB. After running the simulations, she would likely use ad hoc scripts to extract the quantities of interest and then, in the best-case scenario, upload the simulations in a repository such as NOMAD. A multiscale modeling workflow will involve the use of multiple tools and would necessarily become more complex and include a larger number of manual steps. A better solution would be creating a cyber-ecosystem where different cyberinfrastructures and tools can connect to each other without intervention of the user except in the selection of inputs and outputs. This is exemplified on the right panel of figure 17. Here the researcher goes to nanoHUB and interacts with a simulation tool capable of identifying possible data streams for inputs and outputs, and directly connects to the various resources as needed. Such ecosystem would have multiple entry points, for example, simulations could be launched automatically from a data exploration tool.

Zoom In Zoom Out Reset image size

Partially motivated by these needs, the nanoHUB team has included online support for Jupyter notebooks that researchers can use to create, document and publish scientific workflows, including multiscale simulations. These workflows can use web-services to request input data, setup a series of simulations, performing the necessary analysis and conversions between them and display the final results. The glass transition notebook in nanoHUB is one such example [184]. Using the tool, researchers can: (i) select a monomer to download from the ChemSpider database [185], (ii) build an amorphous polymer sample using MC using the Polymer Modeler tool [186], (iii) perform a LAMMPS calculation where the polymer is subject to a thermal cycle, and (iv) postprocess and plot the results to extract the glass transition temperature. Another recent effort involved connecting simulation tools with input data streams. The KimExplorer tool [187] and the Nanomaterials Mechanics Explorer tool [188] query the OpenKIM repository [13] and present the researcher with possible interatomic potentials for their system of interest. Once selected, these potentials are automatically downloaded and installed to perform the simulations.

In summary, modern cyberinfrastructure has the potential to fundamentally change how multiscale modeling of materials is practiced and accelerate progress in the field. Modern cyberinfrastructure can simplify the connection between tools across scales and with sources of data, enable publishing and sharing multiscale workflows, and can put these tools in the hands of domain experts who may not have computational expertize to access them otherwise. In addition, the ability to publish end-to-end workflows reduces human error, contributes to reproducibility, and accelerates progress leading to more reliable predictive simulations of materials behavior.