On a -convergence analysis of a quasicontinuum method††thanks: This work is supported by a grant of the Deutsche Forschungsgemeinschaft (DFG) SCHL 1706/2-1.
In this article, we investigate a quasicontinuum method by means of analytical tools. More precisely, we compare a discrete-to-continuum analysis of an atomistic one-dimensional model problem with a corresponding quasicontinuum model. We consider next and next-to-nearest neighbour interactions of Lennard-Jones type and focus on the so-called quasinonlocal quasicontinuum approximation. Our analysis, which applies -convergence techniques, shows that, in an elastic setting, minimizers and the minimal energies of the fully atomistic problem and its related quasicontinuum approximation have the same limiting behaviour as the number of atoms tends to infinity. In case of fracture this is in general not true. It turns out that the choice of representative atoms in the quasicontinuum approximation has an impact on the fracture energy and on the location of fracture. We give sufficient conditions for the choice of representative atoms such that, also in case of fracture, the minimal energies of the fully atomistic energy and its quasicontinuum approximation coincide in the limit and such that the crack is located in the atomistic region of the quasicontinuum model as desired.
Key words. Quasicontinuum method; atomistic-to-continuum; –convergence; fracture.
AMS subject classifications. 49J45, 74R10, 74G15, 74G10, 74G65, 70C20
The quasicontinuum (QC) method was introduced by Tadmor, Ortiz and Phillips  as a computational tool for atomistic simulations of crystalline solids at zero temperature. The key idea is to split the computational domain into regions where a very detailed (atomistic, nonlocal) description is needed and regions where a coarser (continuum, local) description is sufficient. The QC-method and improvements of it are successfully used to study crystal defects such as dislocations, nanoindentations or cracks and their impact on the overall behaviour of the material, see e.g. .
There are various types of QC-methods: Some are formulated in an energy based framework, some in a force based framework; further, different couplings between the atomistic and continuum parts and different models in the continuum region are considered. In the previous decade, many articles related to the numerical analysis of such coupling methods were published. We refer to [15, 23] for recent overviews, in particular on the large literature including work on error analysis.
In this article, we consider a one-dimensional problem and focus on the so-called quasinonlocal quasicontinuum (QNL) method, first proposed in . The QNL-method and further generalizations of it (see e.g. [16, 30]) are energy-based QC-methods and are constructed to overcome asymmetries (so called ghost-forces) at the atomistic/continuum interface which arise in the classical energy based QC-method.
We are interested in an analytical approach in order to verify the QNL-method as an appropriate mechanical model by means of a discrete-to-continuum limit. This is embedded into the general aim of deriving continuum theories from atomistic models, see e.g. [3, Section 4.1], where also the need of a rigorous justification of QC-methods is addressed.
Our approach, announced in , is based on -convergence, which is a notion for the convergence of variational problems, see e.g. . We start with a one-dimensional fully atomistic model problem which takes nearest and next-to-nearest neighbour interactions into account. The limiting behaviour of the corresponding discrete model was analyzed by means of -convergence techniques in  for a large number of atoms. In particular the -limit and the first order -limit are derived there, which take into account boundary layer effects.
From the fully atomistic model problem, we construct an approximation based on the QNL-method. In particular, we keep the nearest and next-to-nearest neighbour interactions in the atomistic (nonlocal) region and approximate the next-to-nearest neighbour interactions in the continuum (local) region by certain nearest neighbour interactions as outlined below. Furthermore, we reduce the degree of freedom of the energy by fixing certain representative atoms and let the deformation of all atoms depend only on the deformation of these representative atoms.
It turns out that the choice of the representative atoms has a considerable impact on the validity of the QC-method, see Theorem LABEL:5:cor, which is the main result of this work. This theorem asserts that the QC-method is valid if the representative atoms are chosen in such a way that there is at least one non-representative atom between two neighbouring representative atoms in the local region and in particular at the interface between the local and nonlocal regions. In Proposition LABEL:lemma:fracture, we prove that the mentioned sufficient condition on the choice of the representative atoms is indeed sharp by showing that in cases where the condition is not satisfied the limiting energy functional of the QC-method does not have the same minima as the limiting energy of the fully atomistic model and thus should not be considered an appropriate approximation. This implies by means of analytical tools that in numerical simulations of fracture one has to make sure to pick a sufficiently large mesh in the continuum region and at the interface.
The outline of this article is as follows. In Section 2 we present the two discrete models, namely the fully atomistic and the quasicontinuum model, in detail.
In Sections 3 and 4 we investigate the limiting behaviour of the quasicontinuum energy functional by deriving the -limits of zeroth and first order. It turns out that the -limit of zeroth order of the fully atomistic and the quasicontinuum model coincide (Theorem LABEL:theorem:zero). If the boundary conditions are such that the specimen behaves elastically, we prove that both models also have the same -limit of first order (Theorem LABEL:theorem:elastic).
If the boundary conditions are such that fracture occurs, the quasicontinuum approximation leads to a -limit of first order (Theorem LABEL:theorem:fracture) that is in general different from the one obtained earlier for the fully atomistic model (, cf. Theorem LABEL:th:fracturennn). To compare the fully atomistic and the quasicontinuum model also in this regime, we analyze the -limits of first order further in Section 5. As mentioned above, it turns out that if we use a sufficiently coarse mesh in the continuum region, the minimal energies of the two first order -limits coincide (Theorem LABEL:5:cor). In fact we are able to show that in our particular model problem it is sufficient that the mesh size in the continuum region is at least twice the atomistic lattice distance. With this choice, fracture occurs always in the atomistic region as desired.
Furthermore, the -convergence results imply, under suitable assumptions, a rate of convergence of the minimal energy of the quasicontinuum model to the minimal energy of the fully atomistic model (Theorem LABEL:th:limmin). Finally, we show that the condition on the mesh size is sharp. In Proposition LABEL:lemma:fracture, we provide examples where the corresponding -limit has a different minimal energy and different minimizers than the fully atomistic system, which is due to poorly chosen meshes. This yields an analytical understanding of why meshes have to be chosen coarse enough in the continuum region.
Similar models as the one we consider here, were investigated previously in terms of numerical analysis. We refer especially to [14, 21, 26, 28, 29] where the QNL method is studied in one dimension. By proving notions of consistency and stability, those authors perform an error analysis in terms of the lattice spacing. To our knowledge, most of the results do not hold for “fractured” deformations. However, in  a Galerkin approximation of a discrete system is considered and error bounds are proven also for states with a single crack of which the position is prescribed. Recently, a different approach based on bifurcation theory is used in  to study the QC-approximation in the context of crack growth.
In , a different one-dimensional atomistic-continuum coupling method is investigated. Similar as in the QC-method the domain is splitted in a discrete and a continuum region. In the discrete part the energy is given by nearest neighbour Lennard-Jones interaction and in the continuum part by an integral functional with Lennard-Jones energy density. It is shown that fracture is more favourable in the continuum than in the discrete region. To overcome this, the energy density of the continuum model is modified by introducing an additional term which depends on the lattice distance in the discrete region. Furthermore, in [5, p. 420] it is remarked that if the continuum model is replaced by a typical discretized version, the fracture is favourable in the discrete region. As mentioned above, we here treat a similar issue in the QNL-method, see in particular Theorem LABEL:5:cor, Proposition LABEL:lemma:fracture.
The techniques of our analysis of the QNL method are related to earlier approaches based on -convergence for the passage from discrete to continuum models in one dimension, see [8, 9, 10, 11, 12, 31, 32]; see also [18, 19] for a treatment of two dimensional models. Recently, -convergence was used in  to study a QC approximation. In  a different atomistic model, namely a harmonic and defect-free crystal, is considered. Under general conditions it is shown that a quasicontinuum approximation based on summation rules has the same continuum limit as the fully atomistic system.
Common in all those works based on -convergence is that primarily information about the global minimum and minimizers are obtained. Since atomistic solutions are not necessary global minimizers, it would be of interest to obtain also results for local minimizers, for instance in the lines of [7, 9]. In this article, we treat systems with nearest and next-to-nearest neighbour interaction. A natural question is how the sufficient conditions on the choice of representative atoms change if we consider also interacting neighbours, . Therefore the corresponding fully atomistic model has first to be studied, which is part of ongoing research.
2 Setting of the Problem
First we describe our atomistic model problem which is the same as in . We consider a one-dimensional lattice given by with and interpret this as a chain of atoms. We denote by the deformation of the atoms from the reference configuration and write as shorthand. We identify such functions with their piecewise affine interpolations and define
The energy of a deformation is given by
where and are potentials of Lennard-Jones type which will be specified in [LJ1]–[LJ4] below. Moreover, we impose boundary conditions on the first and last two atoms. For given we set
To consider only deformations which satisfy (LABEL:def:boundarycond), we define the functional
The goal is to solve the minimization problem
which we consider as our atomistic problem.
The idea of energy based quasicontinuum approximations is to replace the above minimization problem by a simpler one of which minimizers and minimal energies are good approximations of the ones for . Typically this new problem is obtained in two steps:
Define an energy where the long range (in our case next-to-nearest neighbour) interactions are replaced by certain nearest neighbour interactions in some regions.
Reduce the degree of freedom by choosing a smaller set of admissible functions.
To obtain (a), the next-to-nearest neighbour interactions are approximated as
see e.g. . While this approximation turns out to be appropriate in the bulk, this is not the case close to surfaces, where the second neighbour interactions create boundary layers. This motivates to construct a quasinonlocal quasicontinuum model accordingly: For given let with . For we define the energy by using the above approximation for , cf. Fig. LABEL:fig:qcchain and keeping the atomistic descriptions elsewhere
Analogously to we define the functional
A crucial step for the following analysis is to rewrite the energy in a proper way. By defining
and , sometimes called Cauchy-Born energy density (see ), we can write
for satisfying (LABEL:def:boundarycond). To emphasize the local structure of the continuum approximation, we rewrite the summation over the terms with in (LABEL:hnlqc) as an integral. To this end we use the fact that is constant on for and thus
for satisfying (LABEL:def:boundarycond).
To obtain (b) we consider instead of the deformation of all atoms just the deformation of a possibly much smaller set of so called representative atoms (repatoms). We denote the set of repatoms by with and define
Since we are interested in the energy for deformations , we define
In the following sections we study as tends to infinity. Therefore, we will assume that is such that
Hence, in particular . The above assumption corresponds to the case that the size of the atomistic region becomes unbounded on a microscopic scale (i), but shrinks to a point on a macroscopic scale (ii). While assumption (i) is crucial, see also Remark LABEL:rem:elastic, the assumption (ii) can be easily replaced by , and . In this case the analysis is essentially the same, but in the case of fracture, see Theorem LABEL:theorem:fracture, one has to distinguish more cases. We assume (LABEL:ass:kn) (ii) here because it is the canonical case from a conceptual point of view. Otherwise the atomistic region and continuum region would be on the same macroscopic scale.
3 Zero-Order -Limit
In this section we derive the -limit of the discrete energy (LABEL:def:hnlqct), which we refer to as zero-order -limit. This limit involves the convex and lower semicontinuous envelope of the effective potential energy which is already introduced in  defined by
We state the assumptions on the functions , and under which the following results are obtained.
] (strict convexity)
] (uniqueness of minimal energy configurations) For every such that we have where is defined as
] (regularity and behaviour at , ). be in on their domains such that on its domain. Let and . Moreover, we assume the following limiting behaviour
] (structure of , and ). , are such that there exists a convex function
and there exist constants such that
Further, there exist such that
and is strictly convex in on its domain for . Moreover, it holds and for all .
(a) The main examples we think of are Lennard-Jones interactions, defined classically as
and . The calculations in [31, Remark 4.1] show that defined as above satisfy [LJ1]–[LJ4]. Another example of interatomic potentials which satisfy the above assumptions, see [31, Remark 4.1], are Morse-potentials, defined for as
(b) The assumptions [LJ1]–[LJ4] imply that . In particular, we have
(c) Note that [LJ4] and (LABEL:lim:j012) imply that either or that there exists such that or for . In [LJ3], we assume for simplicity. However, this could be dropped making suitable assumptions on in the following statements.
To define appropriate function spaces, we use a similar notation as in  and . Let be a function with bounded variation. Then we say that if satisfies the Dirichlet boundary conditions and . To allow jumps in respectively , the boundary conditions are replaced by respectively in this case. Analogously, we define for special functions with bounded variations and the above boundary conditions. Let (or in ), then we denote by the jump set of in , and for we set . Moreover we denote by the singular part of the measure with respect to the Lebesgue measure.
Let us now state and prove the zeroth-order -limit of the functional . It turns out that the limiting functional is equal to the -limit of the functional , cf. .
Suppose [LJ1]–[LJ4] are satisfied and let . Let satisfy (LABEL:ass:kn) and let with be such that
Then the -limit of defined in (LABEL:def:hnl) and of defined in (LABEL:def:hnlqct) with respect to the –topology is the functional defined by
Proof. The result for follows from [31, Theorem 3.1]. Thus we prove the result for . The following compactness property and lower bound follow from [10, Theorem 3.7] and [11, Theorem 3.1]. For the readers convenience, we present direct proofs here.
Compactness. Let be a sequence with equibounded energy . The definition of and the properties of , imply that . Define the set . Next, we make use of the fact that are bounded from below and that the energy is equibounded. Moreover, we apply (LABEL:ass:psi2) and Jensen’s inequality to obtain
for some independent of . By (LABEL:ass:psi1), we have that for some constant independent of . Moreover, by using the boundary conditions, we obtain
Since , we obtain by the Poincaré-inequality that is equibounded. Thus, we can extract a subsequence of which converges weakly to some , see [2, Theorem 3.23]. As argued in [31, Theorem 3.1], we have .
Liminf inequality. Let and be a sequence with equibounded energy which converges to in . The above compactness property and [2, Proposition 3.13] imply that converges to weakly in . By using [LJ3], [LJ4], we obtain for the recession function
with arbitrary. For every there exists such that for every . For large enough, we deduce from (LABEL:hnlqc:integral) by the definition of and [LJ4]
Note that by it follows for all , thus there exists such that
The last inequality is a direct implication of [2, Theorem 2.34], using that weakly converges to . By using that the right-hand side above is finite only if , we obtain the liminf inequality from the arbitrariness of .
Limsup inequality. To show the existence of a recovery sequence, we first do not take the boundary conditions into account. Therefore, we define the functional by
For every we show existence of a sequence converging to in such that
As outlined in the proof of [10, Theorem 3.5] it is enough to show the above inequality for linear and for with a single jump: by density, this proves the statement for and the general estimate follows by relaxation arguments. Firstly, we consider functions with a single jump. Let with , and . By (LABEL:T_n) there exists with and such that for . We define now a sequence by
Obviously we have in . The functions are defined such that for and for all with . Using , (LABEL:T_n), [LJ3] and [LJ4] this implies
Now let for some . For every sequence satisfying (LABEL:T_n) we find a sequence of natural numbers such that
We define for every a set with , where such that there exist which satisfy
From we deduce and thus . Let us define such that and
By using and , we obtain
and thus in . Indeed, by , and , the last term tends to zero as . For the limsup inequality we argue similarly as in the case of a jump before. By definition, we have