Local-buckling-induced elastic interaction between inclusions in a free-standing film

International Journal of Solids and Structures 50 (2013) 3742–3747

Contents lists available at ScienceDirect

International Journal of Solids and Structures journal homepage: www.elsevier.com/locate/ijsolstr

Local-buckling-induced elastic interaction between inclusions in a free-standing film Ling Qiao a,b, Linghui He b, Yong Ni b,⇑ a b

Department of Engineering Mechanics, College of Civil Engineering, Southeast University, Nanjing, Jiangsu 210096, PR China CAS Key Laboratory of Mechanical Behavior and Design of Materials, Department of Modern Mechanics, University of Science and Technology of China, Hefei, Anhui 230026, PR China

a r t i c l e

i n f o

Article history: Received 20 October 2012 Received in revised form 5 June 2013 Available online 31 July 2013 Keywords: Thin film Buckling Elastic interaction Inclusion

a b s t r a c t The local-buckling-induced elastic interaction between two circular inclusions in a free-standing film is reported using numerical simulation. The simulation relies on a continuum model based on the modified Föppl-von Kármán plate theory for a film with arbitrarily distributed eigenstrain and eigencurvature. It is shown that due to the overlapping of the nonlinear local buckling the elastic interaction between the two inclusions with the same eigencurvature is repulsive, while the interaction between them with the opposite eigencurvature is attractive. The interaction strength in both cases decays with their mutual distance. In addition, the inclusion with positive/negative eigenstrain above critical values can trigger an axisymmetric/non-axisymmetric buckling, respectively, and the buckling induced elastic interaction between the two inclusions with eigenstrain shows a nonmonotonic behavior. Ó 2013 Elsevier Ltd. All rights reserved.

1. Introduction Defects such as inclusions with eigenstrain contained in the bulk tend to distort its surrounding environment and produce an elastic strain field (Eshelby, 1957). The overlapping effect of the elastic field produced by two inclusions usually results an effective elastic interaction between them, dependent on their mutual distance d. There are only few special cases where the elastic interaction vanishes, i.e. two isotropic defects in an infinite isotropic medium (Bitter, 1931; Eshelby, 1955; Hardy and Bullough, 1967). The elastic interaction between inclusions in a bulk or in a film, or on a substrate surface can be repulsive or attractive, exhibiting different power laws, 1/dn (Bartolo and Fournier, 2003; Goulian et al., 1993; He, 2008; Kralchevsky and Nagayama, 2000; Lau and Kohn, 1977; Marchenko and Misbah, 2002; Peyla and Misbah, 2003; Peyla et al., 1999). Such kind of the interaction is found to play an important role in mechanical properties of materials (Mura, 1987), functional properties of cell membrane (Dan et al., 1993; Kim et al., 1998; Park and Lubensky, 1996), quantum dot self-assembly (Ni and He, 2010) and physisorption (He, 2010). When two inclusions are mounted in a thin film, the elastic strain energy can be further released by buckling the film into the third dimension, with formation of regions of Gauss curvature (Seung and Nelson, 1988). Each inclusion can trigger its own Gauss curvature, and the like-charge interactions between the positive or negative Gauss curvatures appear and thus produce an indirect ⇑ Corresponding author. Tel.: +86 551 63603172. E-mail address: [email protected] (Y. Ni). 0020-7683/$ - see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijsolstr.2013.07.017

nonlinear elastic interaction between the two inclusions (Vitelli and Turner, 2004). Recently, understanding the buckling-induced elastic interaction between inclusions in thin films attract intense interest in functional biomembrane and graphene nanotechnology (Lusk and Carr, 2008; Ursell et al., 2009). In biomembranes, this elastic interaction was found to significantly regulate domain phase separation (Minami and Yamada, 2007), lateral domain organization (Ursell et al., 2009) and anisotropic aggregation of spherical nanoparticles (Saric and Cacciuto, 2012). In graphene sheets, it was reported that the formation energies of the dislocation and disclination defects, as well as the interaction energies between them depend on the film buckling modes, which severely modifies its mechanical and electrical property (Chen and Chrzan, 2011; Yazyev and Louie, 2010; Zhang et al., 2012). When the thin film is a fluid membrane or a lipid bilayer (biomembrane), the membrane induced elastic interaction between inclusions has been extensively addressed, wherein the elastic energy of the film is written as just the Helfrich’s curvature energy and the inclusions are assumed to be rigid and to create local force multipoles (Bartolo and Fournier, 2003; Brannigan and Brown, 2007; Dan et al., 1993; Dommersnes and Fournier,1999; Goulian et al., 1993; Kim et al., 1998; Marchenko and Misbah, 2002; Park and Lubensky, 1996; Weikl et al., 1998). However when the thin film is an elastic plate, it can endure in-plane shear stress and the elastic free energy of the thin film should include the in-plane strain energy plus the Helfrich energy. Furthermore, the inclusions in the film induce not only the out-of-plane fluctuation but also the in-plane fluctuation of the film. The equilibrium configuration of the film perturbed by the inclusions is governed by the total elastic

3743

L. Qiao et al. / International Journal of Solids and Structures 50 (2013) 3742–3747

free energy minimization and we need to simultaneously resolve the in-plane and out-of-plane equilibrium equations for the film. The only solution of the Kirchhoff plate containing an elliptic inclusion or a rotational symmetrical inclusion with eigencurvature was obtained analytically (Beom, 1998; Xu and Wang, 2005). An infinite elastic film containing a circular inclusion with positive and negative eigenstrain always adopts a flat morphology if the film is treated as a Kirchhoff plate, while if the film is undergone moderate out-of-plane deformation, it has to be modeled as a Föppl–von Kármán (FvK) plate, there is an axisymmetric or non-axisymmetric buckling instability driven by the supercritical positive or negative eigenstrain, respectively. The critical conditions for these buckling instabilities driven by the eigenstrain were derived (Li et al., 2010; Qiao et al., 2011), however, the post-buckling profile of the film, as well as the case of a film with multiple inclusions remains elusive. Quantitative analysis of the bucklingmediated elastic interaction requires taking into account the buckling deformation of the film into a fully nonlinear region. In this paper we propose a computational methodology to study the thin film with multiple inclusions based on a modified FvK plate theory with arbitrarily distributed eigenstrain and eigencurvature. The inplane equilibrium of the FvK plate is solved using an Airy stress function and the out-of-plane equilibrium is solved via a Ginzburg–Landau kinetic equation, which steady solution provides the equilibrium buckling profile of the film. Numerical solutions based on the current model indicate that there is a localbuckling-induced elastic interaction between two circular inclusions with eigenstrain and/or eigencurvature in the film.

As shown in Fig. 1, we consider a free-standing isotropically elastic film of thickness h containing multiple inclusions. A shape function h(x1, x2) is defined to characterize the inclusion configuration, i.e. h(x1, x2) = 1 in the inclusions, h(x1, x2) = 0 outside of them. The inclusions are prescribed a uniform eigenstrain, e0ab and eigen0 curvature, kab with a, b = 1, 2. The total strain in the film consists of 0 the transformation strain induced by e0ab and kab , and the strain due el to elastic deformation, eab . Under the assumption of small strain, the total strain is simply the sum of the two parts:



0

0



eab ¼ e þ eab þ x3 kab hðx1 ; x2 Þ:

ð1Þ

According to the FvK plate theory (Mansfield, 1989), the total strain can be expressed by the middle-plane displacement, u = (ua, f)

1 2

1 2

eab ¼ ðua;b þ ub;a Þ þ f;a f;b  x3 f;ab :

ð2Þ

1 2

Z

h=2

h=2

Z

1 1





eelab ¼ ðua;b þ ub;a Þ þ f;a f;b  e0ab h  x3 f;ab þ k0ab h :

1

1 C abdc eelab eeldc dx1 dx2 dx3 : 2

1

ð4Þ

Substituting Eq. (3) into Eq. (4), and integrating with respect to x3 leads to

Etot ¼ Estretch þ Ebend ;

ð5Þ

where

Estretch ¼

Ebend ¼

1 2

D 2

Z

1

Z

1

Z

1

1

1

Nab eab dx1 dx2 ;

ð6Þ

1

Z

1

1

h i 0 H2  2ð1  v f Þkg dx1 dx2 :

ð7Þ

Estretch is the stretching energy of the film with the membrane strain of the film at the middle plane (x3 = 0) and the membrane stress defined as, respectively,

1 1 ðua;b þ ub;a Þ þ f;a f;b  e0ab hðx1 ; x2 Þ; 2 2   ¼ S ð1  v f Þeab þ v f ecc dab ;

eab ¼ Nab

where S ¼

2hlf 1v f

ð8Þ ð9Þ

lf h3

and D ¼ 6ð1v Þ are the 2D Young’s modulus and the f

bending stiffness, respectively, with lf and

vf being the shear mod0

ulus and Poisson’s ratio of the film. H ¼ ðDf þ kcc Þ is the total mean 0

0

0

0

2

curvature and kg ¼ ðf;11 þ k11 Þðf;22 þ k22 Þ  ðf;12 þ k12 Þ is the total Gauss curvature. By using an Airy stress function, v(x1, x2), the membrane stress can be expressed as

ð3Þ

The total elastic strain energy in the film is obtained by integrating through its volume

ð10Þ

where eam is an antisymmetric unit rank-2 tensor. The variation of Etot in Eq. (5) with respect to ua, f yields the in-plane and out-ofplane equilibrium equations for the film, respectively

r4 v ¼ Sð2am 2bn e0ab @ m @ n h þ kg Þ;

ð11Þ

Dr4 f ¼ Nab f;ab h    i 0 0 0 0 0  D k11 þ v f k22 h;11 þ k22 þ v f k11 h;22 þ 2ð1  v f Þk12 h;12 ; ð12Þ f2;12 .

0

0

where kg ¼ f;11 f;22  If eab ¼ eT dab and kab ¼ kT dab , substituting Eq. (11) into Eq. (10), one has the form in the Fourier space

~ ab ¼ Sðna nb  dab Þ N





eT ~h þ k~g =n2 ;

ð13Þ

where

 1=2 n ¼ n21 þ n22 ;

The elastic strain thus has the form

1 2

Z

Nab ¼ 2am 2bn @ m @ n vðx1 ; x2 Þ;

2. Modeling of a FvK plate containing multiple inclusions

el ab

Etot ¼

n1 ¼ n1 =jnj;

n2 ¼ n2 =jnj;

ð14Þ

~ are the Fourier transform of N , h, k , respectively. ~ ab , ~ h, k and N ab g Throughout the paper, Greek indices are used to indicate 1 or 2. The usual summation convention applies for repeated indices, and a comma stands for differentiation with respect to the suffix index. The membrane stress can be solved given the quantities of h and f

Nab ðx1 ; x2 Þ ¼

1 8p 2

ZZ

S





eT ~h þ k~g =n2 ðna nb  dab Þeinx dn1 dn2 :

ð15Þ

Usually it is very difficult to directly solve the coupled integral equations (11) and (12). We can firstly get the solution of the equilibrium in-plane deformation in terms of f(x1, x2), and h(x1, x2) through Eq. (15). After that we can find the solutions of Eq. (12) by solving the following Ginzburg–Landau kinetic equation (Ni et al., 2011) Fig. 1. Sketch of an infinite plate containing two circular inclusions with their radius R and mutual distance d, in which uniform eigenstrain and eigencurvature are prescribed.

@f dEtot ¼ C ; @t df

ð16Þ

L. Qiao et al. / International Journal of Solids and Structures 50 (2013) 3742–3747

where C is a kinetic coefficient which characterizes the relaxation rate of the buckling process in the over-damped dynamics. The steady-state solution of Eq. (16) recovers Eq. (12) and provides the equilibrium distributions of f. In the numerical simulations, the dimensional quantities are scaled by the film thickness (r0 = r/ h), and the time t is scaled by s ¼ h=Clf e2T (i.e. t0 = t/s). Square computational domains of size 1024  1024 are used with periodical boundary conditions and the time step Dt 0 ¼ 0:2. A semi-implicit spectrum method is adopted and the iterative scheme of Eq. (16) is written as

~fðnþ1Þ ¼

~fðnÞ  Dt0

hðN

Dð1þv Þ

ab f;ab Þn

þ l e2 f kT n2 ~h f T

lf e2T

1 þ n4 Dt 0 D=lf e2T

i

FvK plate at kTh=0.002 Kirchoff plate at kTh=0.002

0.010

FvK plate at kTh=0.004 Kirchoff plate at kTh=0.004 FvK plate at kTh=0.02

0.005

Kirchoff plate at kTh=0.02

0.000 -0.005 -0.010

ð17Þ

:

-0.015 0

We numerically solved the coupled Eqs. (15) and (17) with a small random fluctuation of f prescribed, it is clear that the solutions of membrane stress and out-of-plane displacement for the in-plane and out-of-plane elastic equilibrium in the film can be obtained after the steady solution of Eq. (17) (~fðnþ1Þ ¼ ~fðnÞ ) is satisfied. By choosing the film system with different inclusion configuration, the total elastic energy for each case can be calculated. A buckling-induced elastic interaction between inclusions appears when the total elastic energy depends on the mutual distance of the inclusions. 3. Results and discussion 3.1. Buckling of the FvK plate containing a circular inclusion We have noted that the iterative scheme in Eq. (17) for the solution of f is quickly convergent driven by the monotonous minimization of the total elastic energy. Figs. 2 and 3 present the efficiency and accuracy of the numerical solution. Fig. 2 shows the simulated energy changes in the plate as a function of the iteration number during the development of the buckling instability induced by a circular inclusion of the radius R = 30h with a supercritical positive eigenstrain, eT = 8ec,

2

ec ¼ 1am2 ðRhÞ with a is a conf

stant on the order 1, dependent on the boundary condition (Qiao et al., 2011). At the beginning the total elastic energy 1þv

E0 ¼ 2lf e2T 1v f hpR2 is fully stored as the stretching energy in the f

film, the subsequent sharp drop in the stretching energy curve is corresponding to the buckling process. It typically takes about several hundred iterations to reach convergences and thus the stable post-buckling profile can be obtained. Fig. 3 plots the equilibrium radial curvature distribution in the plate containing a circular inclusion of R = 20h with different values of positive eigencurva-

1.0

25

50

75

r/h Fig. 3. The equilibrium radial curvature distribution in the plate containing a circular inclusion with a positive eigencurvature in the case of the FvK plate in comparison with the case of the Kirchhoff plate.

ture, in comparison with the analytical result in the case of Kirchhoff plate. It shows that when the eigencurvature is kT = 0.002/h the numerical solution is almost reduced to the analytical result of the Kirchhoff plate (Beom, 1998), where the radial curvature distribution has a positive constant value in the inclusion and a negative value decaying away from the inclusion. For a large value of the eigencurvature kT = 0.02/h, the numerical solution shows a non-uniform radial curvature in the inclusion and there is a significant deviation from the result of a Kirchhoff plate which is only valid in the case of linear buckling. As we have checked, the convergence and accuracy of the numerical solution are not sensitive on the inclusion configuration. This could allow us to study the case with multiple inclusions. 3.2. Interaction between two circular inclusions with eigencurvature Fig. 4 plots the reduced total elastic energy as a function of d/R in the plate containing two circular inclusions of the same radius R = 20h with different values of eigencurvature. Fig. 4(a) and (b) shows that the total elastic energy has a minimum at d/R = 0 when the eigencurvatures in the two inclusions are opposite while it has a maximum at d/R = 0 when the eigencurvatures in the two inclusions are equal, respectively. This indicates that there is a repulsive elastic interaction between the two inclusions with the same eigencurvature and an attractive elastic interaction between them with the opposite eigencurvature due to the overlapping of the nonlinear local buckling. When the value of the eigencurvature is very small, the buckling is linear as is the case of Kirchhoff plates and the elastic interaction between the two circular inclusions vanishes. 3.3. Buckling of the FvK plate containing inclusions with eigenstrain

bending energy stretching energy total elastic energy

0.8 0.6 Eelas/E0

0.015

krh

3744

0.4 0.2 0.0 2

10

3

10

4

t/ t

10

Fig. 2. The simulated energy changes in the plate with the iteration number during the development of the buckling instability induced by a circular inclusion with a positive eigenstrain.

For an elastic plate containing an inclusion with eigenstrain, the plate may adopt a flat configuration and there is no buckling-induced elastic interaction between inclusions. When the compressive stress due to the eigenstrain is above a critical value, the buckled state is energetically favorable. There are two typical post-buckling shapes of the film containing a circular inclusion with eigenstrain, as shown in the simulated contour maps of f in Fig. 5(a) and (b). Positive eigenstrain in the circular inclusion triggers an axisymmetric buckling, while negative eigenstrain in the inclusion tends to produce a circumferential buckling. Fig. 6 shows that the reduced critical eigenstrain for the circumferential buckling in the film containing a circular inclusion of radius R = 30h is proportional to the second power of circumferential buckling number, n2. We have confirmed that the critical eigenstrain leading

3745

L. Qiao et al. / International Journal of Solids and Structures 50 (2013) 3742–3747

(a)

1.00

25

0.98

20

2

0.96

κT/h=0.005 FvK plate

15

2

2

εT(1-vf )R /h

Etot/Etot

max

Kirchhoff plate κT/h=0.003 FvK plate

0.94

10 5

0.92 0

1

2 d/R

3

4 0 10

(b)

1.00

30 n

Kirchhoff plate κT/h=0.003 FvK plate

0.98

40

50

60

70

2

Fig. 6. The reduced critical eigenstrain for the circumferential buckling in the film containing a circular inclusion of radius R = 30h as a function of the second power of circumferential buckling number, n2. The insertions are the corresponding postbuckling profiles.

κT/h=0.005 FvK plate max

0.96 Etot/Etot

20

0.94 0.92 0.90 0

1

2

3

4

d/R Fig. 4. The reduced total elastic energy as a function of d/R in the plate containing two circular inclusions of the same radius R = 20h with different values of eigencurvature for (a) two inclusions have the opposite eigencurvature; (b) two inclusions have the equal eigencurvature.

If the film is under an additional tension, the compressive region, as well as the buckling area shrinks. Fig. 7 shows the range and amplitude of the buckling in the film containing a circular inclusion with positive eigenstrain decrease significantly as the biaxial isotropic external tension exerted on the film increases. The contour maps for f in Fig. 8 indicate that a uniaxial tension can abruptly change the post-buckling shape of the film containing a circular inclusion with negative eigenstrain. This indicates that the external tension can regulate the local buckling profile and as well as the local-buckling-induced elastic interaction between the inclusions.

to the onset of axisymmetric/non-axisymmetric buckling in the film containing an isolated circular inclusion is 2 n

2

ec ¼ 1Am2 ðRhÞ with f

A = a > 0 for the positive eigenstrain and A  a < 0 for the negative eigenstrain, where an is a coefficient dependent on the circumferential buckling number, n. The critical condition is consistent with our previous result based on the buckling stability analysis and scaling concept (Li et al., 2010; Qiao et al., 2011). The axisymmetric/non-axisymmetric buckling amplitudes are found to decrease to zero up to the infinite distance. We noted that the membrane stress in the pre-stressed film due to the eigenstrain and the external tension before the buckling is

lf h eT for r < R; ð1  v f Þ lf heT R 2 for r > R: Nr ¼ Nh ¼ N0h  ð1  v f Þ r

Nr ¼ Nh ¼ N0r 

(a)

ð18Þ

3.4. Interaction between two circular inclusions with eigenstrain and eigencurvature Fig. 9 plots the reduced total elastic energy as a function of d/R in the plate containing two circular inclusions of the same radius R = 20h with eigenstrain, eT = 0.005 > ec. The buckling can be in the same or opposite direction, the total elastic energy shows a nonmonotonic behavior and has a minimum at d/R = 0 in both cases. However, the case with the buckling in the opposite direction has a lower energy at most values of the mutual distance. This is consistent with the previous report where the domain coalescence driven by minimization of domain boundary energy and the domain-induced interaction is mainly observed between domains with opposite buckling directions (Minami and Yamada, 2007).

(b)

Fig. 5. The simulated contour map of the buckling profile, f(x1, x2), in the plate containing a single isolated circular inclusion of the radius R = 30h with (a) eT = 0.002, (b) eT = 0.004.

3746

L. Qiao et al. / International Journal of Solids and Structures 50 (2013) 3742–3747

5 εΤ=0.005

4

⎯ε=0 ⎯ε=0.0005 ⎯ε=0.001

0.40 Eelas/E0

3

/h

opposite budding same budding

0.41

2

0.39

1 0.38

0 0

100

200

300

400

0

500

1

2

3

4

5

6

d/R

r Fig. 7. The simulated buckling profile, f(r, h), along the radial direction in the plate containing a single isolated circular inclusion of the radius R = 30h with a positive eigenstrain, eT = 0.005 under different equibiaxial external strain boundary conditions.

Fig. 9. The reduced total elastic energy as a function of d/R in the plate containing two circular inclusions of the same radius R = 20h with the same positive eigenstrain, eT = 0.005.

Eq. (19) can be solved by using Green’s function When the inclusions in the film contain both eigenstrain and eigencurvature the buckling profile as well as the buckling-induced elastic interaction can be regulated. For example, the inclusion with eT > 0 and kT > 0 favors buckling into a budding profile, while the inclusion with eT > 0 and kT < 0 favors buckling into a dimple profile. Fig. 10 plots the reduced total elastic energy as a function of d/R in the plate containing two circular inclusions of the same radius R = 20h with the opposite eigencurvature and the same positive eigenstrain of different values. The buckling-induced elastic interaction is found to be the sum contributions due to the eigenstrain and eigencurvature, but is nonadditive.

fðxÞ ¼ 

ZZ

0

0

Gðx  x0 Þð1 þ mf ÞkT r2 hðx0 Þdx1 dx2 ;

ð20Þ

where the Green’s function G is given by

G¼

1 2 ðx  x0 Þ lnðx  x0 Þ: 8p

ð21Þ

Substituting Eq. (20) into Eq. (12) and subtracting the elastic energy due to each inclusion, we obtain the elastic interaction energy between the two inclusions

3.5. Discussion The buckling-induced elastic interaction between two circular inclusions in the elastic plate is different from that in a fluid film or lipid bilayer. One difference is that the inclusion in the former is elastically deformable, adaptive with the in-plane and outof-plane deformation of the film, while the inclusion in the latter is imposed to be rigid. The other difference is that the former involves in-plane strain energy relaxation of the film and the latter cannot. If we neglected the contribution of the in-plane strain energy in the current model, for the two inclusions with isotropic 0 eigencurvature kab ¼ kT dab the equilibrium configuration of the film dominated by minimization of the curvature elastic energy shown in Eq. (7) is determined by the simplified Eq. (12)

r4 f ¼ ð1 þ v f ÞkT r2 h:

(a)

ð19Þ

y

=0

(b) y

1.00

Etot/Etot

max

0.98

κT/h=0.005, εT=0.001 κT/h=0.005, εT=0.002 κT/h=0.005, εT=0.003

0.96

κT/h=0.005, εT=0.004

0.94

0.92 0

1

2

3

d/R 4

Fig. 10. The reduced total elastic energy as a function of d/R in the plate containing two circular inclusions of the same radius R = 20h with the opposite positive eigencurvature and the same eigenstrain of different values.

= −0.025

T

(c)

y

= −0.05

T

Fig. 8. The buckling profile, f(x1, x2), in the plate containing a single isolated circular inclusion of the radius R = 30h with a negative eigenstrain, eT = 0.01 under different levels of uniaxial stretching, (a) ey ¼ 0, (b) ey ¼ 0:025eT , (c)  ey ¼ 0:05eT .

L. Qiao et al. / International Journal of Solids and Structures 50 (2013) 3742–3747

Eint ¼ 

ZZ

2

0

0

Dð1 þ mf Þ2 kT hð1Þ hð2Þ dx1 dx2 ;

ð22Þ

where h(1) and h(2) are the shape functions of the two inclusions respectively. It is clear that the interaction energy of Eq. (22) is zero. The result is consistent with the prediction that the buckling-induced elastic interaction vanishes for two isotropic inclusions in the membrane (Dommersnes and Fournier,1999; Marchenko and Misbah, 2002; Park and Lubensky, 1996). If the contribution of inplane strain energy of the film is noticeable, from Eq. (15) it is clear that the membrane force Nab is strongly coupled with the Gauss curvature kg. The equilibrium equation for f becomes highly nonlinear

Dr4 f ¼ Nab f;ab  Dð1 þ v f ÞkT r2 h:

ð23Þ

Taking the solution of Eq. (19) as the zeroth-order solution of Eq. (23), the next order solution f(n+1) of Eq. (23) can be obtained from a lower-order solution f(n). fðnþ1Þ ðxÞ ¼

ZZ

 1 0 0 ðnÞ Gðx  x0 Þ N ab ðfðnÞ ðx0 ÞÞf;ab ðx0 Þ  ð1 þ mf ÞkT r2 hðx0 Þ dx1 dx2 : D

ð24Þ

ðnþ1Þ

From Eq. (24), we can see that the higher order solution fð1Þ for the equilibrium out-of-plane displacement induced by inclusion 1 is also influenced by inclusion 2. Thus the interaction energy between the two inclusions is not zero anymore and depends on the mutual distance even for two isotropic inclusions, as shown in our numerical simulation. There is an induced interaction consistent to the report by Marchenko and Misbah (2002), wherein the term for the induced effect is different. The induced repulsive interaction between two isotropic inclusions with the equal eigencurvature is similar to the interaction of two conical membrane inclusions of equal orientation, also the interaction of two inclusions with the opposite eigencurvature is similar to the interaction of two conical membrane inclusions of opposite orientation (Weikl et al., 1998). 4. Conclusions In summary, we have developed a continuum model to study the post-buckling of the film with arbitrarily distributed eigenstrain and eigencurvature. Based on the modeling and simulation, the critical buckling condition and post-buckling profile for an infinite thin film containing a circular inclusion with eigenstrain and eigencurvature have been identified. The buckling-induced elastic interaction between two circular inclusions with eigenstrain and eigencurvature containing in the film is reported. The elastic interaction can be attractive or repulsive, dependent on the distribution of the eigenstrain and eigencurvature in the inclusions. Finally under given the distribution of multiple inclusions with the eigenstrain and eigencurvature in the film, the resultant post-buckling profile of the film, as well as the buckling-induced domain interaction can also be simulated. Acknowledgements This work was supported by the Basic Research Program of China (Grant No. 2011CB302101), the Chinese Natural Science Foundation (Grant Nos. 11072232, 11132009), the ‘‘Hundred of Talents Project’’ of the Chinese Academy of Sciences, L. Qiao acknowledges the supports of the Priority Academic Program Development of Jiangsu Higher Education Institutions, and the open research fund of Key Laboratory of MEMS of Ministry of Education, Southeast University.

3747

References Bartolo, D., Fournier, J.B., 2003. Elastic interaction between ‘‘hard’’ or ‘‘soft’’ pointwise inclusions on biological membranes. Eur. Phys. J. B 11, 141–146. Beom, H.G., 1998. Analysis of a plate containing an elliptic inclusion with eigencurvatures. Arch. Appl. Mech. 68, 422–432. Bitter, F., 1931. On impurities in metals. Phys. Rev. 37, 1527–1547. Brannigan, G., Brown, F.L.H., 2007. Contributions of Gaussian curvature and nonconstant lipid volume to protein deformation of lipid bilayers. Biophys. J. 92, 864–876. Chen, S., Chrzan, D.C., 2011. Continuum theory of dislocations and buckling in graphene. Phys. Rev. B 84, 214103. Dan, N., Pincus, P., Safran, S.A., 1993. Membrane-induced interactions between inclusions. Langmuir 9, 2768–2771. Dommersnes, P.G., Fournier, J.B., 1999. N-body study of anisotropic membrane inclusions: membrane mediated interactions and ordered aggregation. Eur. Phys. J. B 12, 9–12. Eshelby, J.D., 1955. The elastic interaction of point defects. Acta Metall. 3, 487–490. Eshelby, J.D., 1957. The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proc. R. Soc. London Ser. A 241, 376–396. Goulian, M., Bruinsma, R., Pincus, P., 1993. Long-range forces in heterogeneous fluid membranes. Europhys. Lett. 22, 145–150. Hardy, J.R., Bullough, R., 1967. Point defect interactions in harmonic cubic lattices. Philos. Mag. 15, 237–246. He, L.H., 2008. Elastic interaction between force dipoles on a stretchable substrate. J. Mech. Phys. Solids 56, 2957–2971. He, L.H., 2010. Mechanics of physisorption on elastomer surface. J. Mech. Phys. Solids 58, 1195–1211. Kim, K.S., Neu, J., Oster, G., 1998. Curvature-mediated interactions between membrane proteins. Biophys. J. 75, 2274–2291. Kralchevsky, P.A., Nagayama, K., 2000. Capillary interactions between particles bound to interfaces, liquid films and biomembranes. Adv. Colloid Interface Sci. 85, 145–192. Lau, K.H., Kohn, W., 1977. Elastic interaction of two atoms adsorbed on a solidsurface. Surf. Sci. 65, 607–618. Li, K., Yan, S.P., Ni, Y., Liang, H.Y., He, L.H., 2010. Controllable buckling of an elastic disc with actuation strain. Europhys. Lett. 92, 16003. Lusk, M.T., Carr, L.D., 2008. Nanoengineering defect structures on graphene. Phys. Rev. Lett. 100, 175503. Mansfield, E.H., 1989. The Bending and Stretching of Plates. Cambridge University Press, Cambridge. Marchenko, V.I., Misbah, C., 2002. Elastic interaction of point defects on biological membranes. Eur. Phys. J. E 8, 477–484. Minami, A., Yamada, K., 2007. Domain-induced budding in buckling membranes. Eur. Phys. J. E 23, 367–374. Mura, T., 1987. Micromechanics of Defects in Solids. Kluwer, Dordrecht. Ni, Y., He, L.H., 2010. Spontaneous formation of vertically anticorrelated epitaxial islands on ultrathin substrates. Appl. Phys. Lett. 97, 261911. Ni, Y., He, L.H., Liu, Q.H., 2011. Modeling kinetics of diffusion-controlled surface wrinkles. Phys. Rev. E 84, 051604. Park, J.M., Lubensky, T.C., 1996. Interactions between membrane inclusions on fluctuating membranes. J. Phys. I 6, 1217–1235. Peyla, P., Misbah, C., 2003. Elastic interaction between defects in thin and 2D films. Eur. Phys. J. B 33, 233–247. Peyla, P., Vallat, A., Misbah, C., Muller-Krumbhaar, H., 1999. Elastic interaction between surface defects in thin layers. Phys. Rev. Lett. 82, 787–790. Qiao, L., He, L.H., Ding, K.W., 2011. Axisymmetric buckling of an elastic plate containing a circular inclusion with dilative eigenstrain. J. Appl. Mech. 78, 014503. Saric, A., Cacciuto, A., 2012. Fluid membranes can drive linear aggregation of adsorbed spherical nanoparticles. Phys. Rev. Lett. 108, 118101. Seung, H.S., Nelson, D.R., 1988. Defects in flexible membranes with crystalline order. Phys. Rev. A 38, 1005–1018. Ursell, T.S., Klug, W.S., Phillips, R., 2009. Morphology and interaction between lipid domains. PNAS 106, 13301–13306. Vitelli, V., Turner, A.M., 2004. Anomalous coupling between topological defects and curvature. Phys. Rev. Lett. 93, 215301. Weikl, T.R., Kozlov, M.M., Helfrich, W., 1998. Interaction of conical membrane inclusions: effect of lateral tension. Phys. Rev. E 57, 6988–6995. Xu, B.X., Wang, M.Z., 2005. The quasi Eshelby property for rotational symmetrical inclusions of uniform eigencurvatures within an infinite plate. Proc. R. Soc. A 461, 2899–2910. Yazyev, O.V., Louie, S.G., 2010. Topological defects in graphene: dislocations and grain boundaries. Phys. Rev. B 81, 195420. Zhang, J.F., Zhao, J.J., Lu, J.P., 2012. Intrinsic strength and failure behaviors of graphene grain boundaries. ACS Nano 6, 2704–2711.