Stability of 2D triangular lattice under finite biaxial strain Текст научной статьи по специальности «Физика»

UDC 539.3

STABILITY OF 2D TRIANGULAR LATTICE UNDER FINITE BIAXIAL STRAIN

E. A. Podolskaya1, A. Yu. Panchenko1, A. M. Krivtsov1

institute for Problems in Mechanical Engineering (IPME RAS), Saint-Petersburg, Russia

katepodolskaya@gmail.com

PACS 46.25.Cc, 68.35.Rh

Stability of 2D triangular lattice under finite biaxial strain is investigated. In this work only diagonal strain tensor is regarded. The lattice is considered infinite and consisting of particles which interact by pair force central potential. Dynamic stability criterion is used: frequency of elastic waves is required to be real for any real wave vector. Two stability regions corresponding to horizontal and vertical orientations of the lattice are obtained. It means that a structural transition, which is equal to the change of lattice orientation, is possible. The regions' boundaries are explained: wave equation coefficients change their signs at the border, as well as Young modulae and shear modulae. The results are proved by direct numerical simulation.

Keywords: stability, triangular lattice, finite strain, biaxial strain, pair potential, elastic wave, structural transition.

1. Introduction

In this work stability of plane triangular lattice under finite biaxial strain is investigated. This problem and its possible solutions have been considered earlier in [1]. Then in [2] a qualitatively new result, that is proved and discussed in detail in this work, was obtained. In spite of the fact that real objects with such lattice have not been discovered, triangular lattice is convenient for both constructing the theory and carrying out numerical simulations. Finite strain may lead to discontinuities, so in this work analysis and modeling based on discrete atomistic methods [3] is proposed. The medium is represented by a set of particles interacting by a pair force central potential, in particular Lennard-Jones and Morse. Direct tensor calculus [4] is used.

Let us introduce the following notation to describe the geometry:

ak = rk-r0, (1)

where rk is radius vector of a particle k, r0 is radius vector of reference particle. If a lattice is simple, then any particle can be named "reference", each particle k has a pair —k and a_k = — ak. Triangular lattice is simple and close-packed, which means that it coincides with its Bravais lattice and possesses maximum concentration of nods in elementary volume V0 with the given minimum distance between the nods. Let us refer to the geometry which is described by rk and ak as reference configuration.

o

Let V and V be Hamilton's operators in reference and current configurations [4]:

o d d

V = e —, V = e —. (2)

dxi ~ldX,

Vectors ej form an orthonormal basis. If vector r has projections .r,: in reference configuration, then in current configuration r will turn into R with projections Xi in the same basis.

Suppose that the lattice is subject to strain characterized by VR. According to long-wave approximation [3,5]

Âk = E (r - Qk) — £ (r) ~ «fc-Vñ.

(3)

Long-wave approximation takes into account those wave lengths that are much greater than the interatomic distance. The thermal motion is neglected.

Morse and Lennard-Jones potentials are used in this work to describe the interaction between particles

n(r) = D

r2^"1) -

nLJ (r) = D

12

- 2

^ 6 r

(4)

Here a is equilibrium distance in the system of two particles, D is the depth of potential well, B characterizes the width of the well. If B = 6, these potentials coincide in the elastic zone. If r ^ 0, then nLJ ^ to, but Morse potential remains finite. This is very important for numerical simulations of strong compression. If r ^ to, then Morse potential decreases faster, so less particles may be taken into consideration. That is why Morse potential is preferable in this work.

Let Fk = F(Ak) = -n'(Ak) be interaction force and Ck = C(Ak) = n"(Ak) be the bond stiffness, both calculated in current configuration. Then we can introduce

4

ÂkÂk>

%

ÂkÂkÂkÂk->

$ = -

Fk Ak

i.k

$

Qk—k ">

2

*fe

lB

Qk—k—k—k

I>44-

(5)

2. Stability criterion

Following [1,2], let us write the equation of motion in Piola form [4]

Pou = V-E, (6)

where p0 is density in the reference configuration, u = R — r is displacement vector, P is Piola stress tensor. According to [3]

1

(7)

To investigate the stability of current configuration the first variation of (6) is found [1]

p06u = V-6E. (8)

o

If initial deformation characterized by VR is uniform, the following equation is obtained from (8) using (5)

mv=4Q---VVv, (9)

where m = p0V0, v = 5u, 4C = E$+4B.

The solution of (9) can be represented in the following form

v = v0eluJtelK'E, (10)

where uj is frequency, K is wave vector. Substitution of solution (10) into equation (9) leads to a system of homogeneous linear equations for amplitude v0. This system has non-trivial solution, if

det [22 —

where Q = mui2, D=4C--K, K = KK.

0,

(11)

a

r

1

Dynamic stability criterion is: frequency of elastic waves is required to be real for any real wave vector. Thus, the following has to be fulfilled for any real vector K, so that v will remain infinitesimal

Q > 0. (12)

3. Deformation of triangular lattice

Fig. 1. Reference and current configurations

Fig. 1 shows the typical part of triangular lattice before and after deformation. In reference configuration a = 60°.

Let e1 and e2 be the unit vectors of directions 1 and 2 respectively. In 2D case (11) takes the form

Q2 - ntrU + detQ = 0. (13)

According to (12) roots of equation (13) are positive for stable current configurations. Thus, stability criterion is

tr22 > 0, det|2>0, 2trg2 - (tr|2)2 ^ 0.

(14)

Inequality 2tr D2 — (tr D)2 ^ 0 is always true in 2D. Let G = E--4C.~The equations (14) yield

trQ > 0 ^ GnK'i + G22K22 > 0, detD > 0 ^ AK\ + 2BK'2K'2 + CK\ > 0,

where

G ii = C11 + C21,

G

22

C12 + C'

22,

A = C11C21, 2B = C11C22 + ^12^21 _ 4G44, C = C12C22

(15)

(16)

Here Cij are the components of tensor4^1 from (9).

The stability regions for Morse potential with 6 = 6 are shown in Fig. 2. Here e1 and e2 are linear parts of Cauchy-Green tensor. Two coordination spheres are considered, because taking only the first sphere into account is not sufficient in the case of finite strain [2]. Fig. 3 shows the transition from the vertical orientation to the horizontal orientation of the lattice. Particles from first coordination sphere in vertical lattice are marked by circles, stars correspond to second sphere. After deformation that is equal to rotation all particles change their position, and some "stars" occupy the place of "circles" in close vicinity to reference particle. If the second sphere is neglected, there will be no equilibrium at horizontal lattice [1,2].

Let us draw a series of stress-strain diagrams, e.g. Fig. 4. According to [3] Cauchy stress tensor has the form

i = (17)

Fig. 2. Stability regions

Fig. 3. Coordination spheres

where V = y/?>/2{l + £i)(l + £2)- Grey zones in Fig. 4 correspond to stability regions. In Fig. 4 we can see, that the loss of stability is strongly connected with the sign of da fde.

Let us evaluate the signs of Young modulae and shear modulae (the plural is due to anisotropy of the deformed lattice).

In Fig. 5 we can see, that dynamic and "static" stability (positivity of local elastic modu-lae) coincide in the case of diagonal stress tensor.

The number and the shape of stability regions strongly depend on the potential (see Fig. 6). For example, intermediate small stability area vanishes with the decrease of d (increase of width of potential well).

4. MD simulation

In order to verify the obtained results, numerical simulation is carried out. The simulation method is based on molecular dynamics (MD) technique and it is described in [3,6]. A set of current configurations with periodic boundary conditions to model an infinite lattice is exposed to additional minor strains caused by small chaotic velocities, given to each particle. At each

Fig. 4. Uniaxial loading (u2 = 0)

Fig. 5. Young modulae and shear modulae

step the following system is integrated numerically

N

mrt

J>(lrfc-rj)

n=1

— k —n

I Tfc ~~ ij

k = N,

(18)

where N is the total number of particles in the simulation. The total kinetic energy of the system is calculated. If a sudden growth of kinetic energy is observed, the configuration is considered unstable. An example of unstable configuration evolution is presented in Fig. 7. Black ovals mark "crack" initiation zones.

In Fig. 8 the results for diagonal strain tensor are presented. We can see a good agreement between two approaches. Discrepancies at high compression and at the borders are caused by

Fig. 6. Stability regions for different potentials

Fig. 7. Configuration e1 = 0.506, e2 = —0.108 after 105 and 3• 105 integration steps 0.2-,

MD

Theory .:»;»::::!:!:

0.0-

-0.2-1

'l

-0.4 A

-0.6-

-0.8

—i-1-1-1-1-1-1-1-1-1-1-1-1-1-1-1-1-.-1

-0.8 -0.6 -0.4 -0.2 0.0 o 0.2 0.4 0.6 0.8 1.0

Fig. 8. Results of numerical simulation

computational errors, perhaps by insufficient number of steps and also by the difference between the number of configurations regarded (106 in theoretical analysis and 104 in numerical simulation due to big time consumption of the latter).

5. Concluding remarks

Stability of 2D triangular lattice under finite biaxial strain was investigated. Two stability regions, which correspond to vertical and horizontal orientations of the lattice, were obtained both analytically and using MD simulation. It was shown that taking more than one coordination sphere into account leads to a new effect: possibility of structural transition, which is equal to the change of lattice orientation. It was noticed that dynamic and static stability coincide in the case of diagonal strain tensor. It was also shown that stability loss during hydrostatic compression is connected with shear modulus, which coincides with results of [7] obtained for 3D case. To sum up, the analytical approach was proved to be consistent and thus it can be applied to more complex cases, e.g. arbitrary strain of triangular lattice and FCC lattice.

References

[1] Tkachev P.V., Krivtsov A.M. Stability criterion of internal structure of the material // XXXIII Science week, St. Petersburg State Polytechnical University, Part V, 2004. P. 4-6.

[2] Podolskaya E.A., Tkachev P.V., Krivtsov A.M. FCC latice stability analysis // XXXIX Science week, St. Petersburg State Polytechnical University. Part V, 2010. P. 111-112.

[3] Krivtsov A.M. Deformation and fracture of solids with microstructure [in Russian]. — M.: Fizmatlit, 2007. — 304 p.

[4] Lurie A.I. Nonlinear theory of elasticity. — Amsterdam: North-Holland, 1990. — 617 p.

[5] Born M., Huang K. Dynamical theory of crystal lattices. — Oxford: Clarendon Press, 1954. — 420 p.

[6] Krivtsov A.M. Molecular dynamics simulation of plastic effects upon spalling // Physics of the Solid State, 2004. V. 46(6). P. 1055-1060.

[7] Milstein F., Rasky D. Theoretical study of shear-modulus instabilities in the alkali metals under hydrostatic pressure // Phys. Rev. B., 1996. V. 54, No. 10. P. 7016-7025.