Critical buckling load of honeycomb sandwich panels reinforced by three-phase orthotropic skins enhanced by carbon nanotubes

Journal Pre-proofs Critical buckling load of honeycomb sandwich panels reinforced by three-phase orthotropic skins enhanced by carbon nanotubes Michele Bacciocchi, Angelo Marcello Tarantino PII: DOI: Reference:

S0263-8223(19)34650-1 https://doi.org/10.1016/j.compstruct.2020.111904 COST 111904

To appear in:

Composite Structures

Received Date: Revised Date: Accepted Date:

7 December 2019 7 January 2020 8 January 2020

Please cite this article as: Bacciocchi, M., Marcello Tarantino, A., Critical buckling load of honeycomb sandwich panels reinforced by three-phase orthotropic skins enhanced by carbon nanotubes, Composite Structures (2020), doi: https://doi.org/10.1016/j.compstruct.2020.111904

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2020 Published by Elsevier Ltd.

Critical buckling load of honeycomb sandwich panels reinforced by three-phase orthotropic skins enhanced by carbon nanotubes Michele Bacciocchi1,2*, Angelo Marcello Tarantino3

Abstract The buckling analyses illustrated in this research aim to provide useful results for the design and application of sandwich plates with a honeycomb core and three-phase orthotropic skins, which are reinforced by both carbon nanotubes (CNTs) and straight oriented fibers. A multiscale approach is developed to this aim, which is based on the Eshelby-Mori-Tanaka scheme and the Hahn homogenization technique. The outcomes are presented in terms of critical buckling loads for different boundary conditions, lamination schemes, fiber orientation and mass fraction of CNTs, in order to prove that all these elements represent fundamental design parameters in the analysis, manufacturing and behavior of these sandwich plates. The theoretical framework is based on the Reissner-Mindlin theory for laminated plates and on the von Kármán hypothesis as far as the nonlinear terms are concerned. The kinematic model includes the Murakami’s function to deal with such peculiar mechanical configurations, which is required to capture the zig-zag effect due to the different mechanical properties of the core and the external skins. The numerical approach and the theoretical methodology are validated by means of the comparison with the experimental and the theoretical results available in the literature.

Keywords Critical buckling load; Multiscale model; Honeycomb sandwich panels; Numerical analysis; Carbon nanotubes reinforced structures.

DICAM Department, University of Bologna, Bologna, Italy DESD Department, University of San Marino, San Marino 3 DIEF Department, University of Modena and Reggio Emilia, Modena, Italy * Corresponding author. Email address: [email protected] (Michele Bacciocchi) 1 2

1

1. Introduction The critical buckling load of a structural system represents typically the value of applied force that causes its loss of stability. In general, when the critical load is reached, an evolution from one equilibrium state to another occurs, such as a transition from a stable equilibrium path to an unstable one [1-7]. The first studies of buckling phenomena, which were focused on the stability of slender beams subjected to compression loads, were attributable to Euler [8]. Then, important contributions to this topic were provided by Poincaré [9], Schmidt [10], Liapunov [11], von Kármán [12], and Koiter [13], as stated in the book by Bloom and Coffin [14]. The structural response of flat panels or plates is affected by the application of in-plane compression forces. In these circumstances, if the critical buckling load is exceeded, a sudden transverse deflection may occur due to the buckling instability [15-18]. Such phenomenon is influential also in laminated composite plates made of a sequence of thin layers [19, 20]. Due to the typical slenderness that can be reached because of the higher level of stiffness which require less material than a conventional medium, the buckling analysis should be conveniently performed to avoid the buckling instability [21-24]. Recent advances in this regard have been presented by Vescovini and Bisagni [25], Vescovini and Dozio [26], Vescovini et al. [27], and Ruocco and Reddy [28]. A peculiar class of laminated plates is the one of sandwich panels, as illustrated in the technical report by Hemp [29] and in the books by Plantema [30], Allen [31] and Vinson [32]. These structures are made of a quite thick internal core, usually characterized by weaker mechanical properties, and two external reinforcing plies or skins, which provide additional strength to the element. This aspect has been clearly described in the papers [33-38]. Such structures are widely used in many engineering fields, such as aerospace, mechanical and automotive engineering, in which the buckling resistance is crucial [39-42]. This kind of analysis is even more interesting if the internal core of the sandwich panels has a honeycomb structure [43-47]. Therefore, the current paper is focused on the Buckling response of sandwich plates with a honeycomb core. A brief and concise literature review of such topic is presented. Apologies are given in advances for any important contribution that the authors may have omitted. In particular, the critical buckling loads of clamped and simply supported compressed sandwich plates with a honeycomb core were evaluated by Boller through experimental tests [48]. Then, the same results were obtained and validated analytically by Thurston [49]. An analytical study of the buckling characteristic of an orthotropic simply supported sandwich panel subjected to compressive loads and a uniform lateral pressure was presented by Robinson [50]. On the other hand, the paper by Chang et al. [51] illustrated a general stability analysis of orthotropic sandwich panels by means of buckling coefficient curves based on an energy method for several

2

boundary conditions. Benson and Mayers developed the governing equations by means of a variational approach based on potential energy and Lagrange multipliers to study the stability of sandwich plates with honeycomb core [52]. Pearce and Webber investigated the buckling response of simply supported sandwich panels made of fiber reinforced skins and honeycomb cores. In their papers [53, 54], they presented the results of both numerical and experimental tests. The Rayleigh-Ritz method was used instead by Rao to study the effect of the skin anisotropy on the critical buckling loads of sandwich plates [55]. A wide experimental campaign was carried out by Minguet et al. [56] to study the compressive failure of sandwich plates. Various boundary conditions, mechanical configurations, and geometries were taken into account to this aim. Then, an experimental investigation focused on the buckling characteristics of sandwich structures made of a honeycomb core reinforced by fibrous orthotropic skins was described in the research carried out by Yeh and Wu [57]. An analytical model was developed and presented in the paper by Kwon et al. [58] to investigate the critical loads of sandwich composites with a honeycomb core and glass-reinforced plastic skins. Their results were validated also by means of experimental tests carried out by the same authors. A higher-order theory was developed instead by Frostig to analyze the buckling behavior of sandwich panels with a soft-core [59], and a closed form solution was presented for simply supported plates loaded by in-plane forces. In the papers by Dawe and Yuan [60] and by Yuan and Dawe [61], the B-spline Finite Strip Method was employed to compute the buckling loads of rectangular sandwich plates with fiber-reinforced composite laminated skins. Analytical and numerical considerations were described and presented in detail. Another experimental campaign that should be mentioned was presented in the paper by Roberts et al. [62], which was focused on the buckling analysis of rectangular orthotropic sandwich panels. The buckling behavior of sandwich panels with sinusoidal honeycomb core was investigated by Davalos and Chen [63] taking into account partially restrained boundary conditions. Their results were obtained both analytically and experimentally. A similar approach was followed by Kaman et al. [64] for the evaluation of the critical buckling load of honeycomb sandwich panels. A refined two-dimensional plate model was developed by means of a layerwise kinematics within a partially mixed variational statement for the buckling analysis of sandwich panels by D’Ottavio and Polit [65]. In particular, a closed form solution based on the Navier approach was used to obtain the buckling loads. Analogously, the Sublaminate Generalized Unified Formulation was employed in the paper by D’Ottavio et al. [66] to investigate the buckling behavior of the same structures. Experimental tests were carried out and presented in the papers by Muc et al. [67] and by Nazari et al. [68] for the buckling behavior of sandwich cylindrical shells and panels, respectively. Finally, the global buckling behavior, as well as the wrinkling

3

phenomenon, of sandwich plates made of honeycomb and soft cores reinforced by anisotropic skins was studied by means of a higher-order plate theory based on a sublaminate variable-kinematic approach by Vescovini et al. [69]. Several configurations were analyzed and discussed, and the results of many numerical applications obtained through the Ritz method were presented. In particular, the aforementioned paper should be considered for the thorough and systematic illustration of the wrinkling behavior of such composites. For completeness purposes, it should be also recalled that some buckling phenomena can be related to the viscoelastic features of the constituents [70]. Further details concerning the theory of viscoelasticity can be found in the papers [71-75]. The approach presented in this paper aims to provide a fast, immediate and easily implementable tool for the evaluation of the global buckling response of sandwich structures. The critical loads are computed numerically by means of a Finite Element (FE) code [76, 77] developed by the authors, for different boundary conditions, multiaxial load cases, and fiber orientations. The proposed methodology does not require complex kinematic models such as the ones that characterize Higher-order Shear Deformation Theories, but it is based on the wellknown Reissner-Mindlin (RM) theory for laminated plates [78, 79] enriched only by the Murakami’s function [80, 81]. This approach is termed Reissner-Mindlin Zig-zag (RMZ) model in the current paper for conciseness purposes. As illustrated by Carrera [82-85] and shown in the papers [86-89], this zig-zag model represents an accurate and reliable approach to deal with sandwich structures characterized by layers with noticeably differences in terms of mechanical properties, without resorting to higher-order approaches, such as Layer-Wise (LW) models, or more demanding techniques, such as three-dimensional FE methodologies or local strip models [69]. The key aspects of this approach are presented in Section 2. With respect to conventional configurations, in the current paper the external skins are modeled as a multiphase material [90-93]. In particular, these layers are made of a polymer matrix enriched by both oriented straight fibers and randomly oriented carbon nanotubes (CNTs) [94, 95]. A highlighted in the papers [96-99], this reinforcing phase allows to improve noticeably the mechanical behavior of the structures in which they are applied. At the nano-scale level, the Eshelby-Mori-Tanaka approach is employed to characterize the engineering constants of the isotropic composite made of a polymer matrix enriched by CNTs [100, 101], whereas the overall properties taking into account also the effect of the straight fibers are given by the Hahn methodology [102]. The introduction of a multiphase approach allows to consider two different reinforcing phases. Differently from the usual configurations described and analyzed in the literature, the mechanical features of fiber-reinforced layers are further enhanced by

4

means of the nanofibers scattered in matrix. This inclusion can change completely the global buckling response of the structure under consideration. The mechanical characterization of these composites is illustrated in Section 3. The presence of CNTs has a great influence in the buckling response of these structures and allows to change completely the solutions, given in terms of critical buckling loads. The numerical applications, which are presented in Section 4, prove that a small amount of these reinforcing nanofibers is able to affect noticeably the buckling behavior of these structures. Therefore, the results illustrated in this paper should be considered as design criteria by engineers and manufacturers. Finally, an Appendix is included in the paper to illustrate briefly the FE formulation of the current approach.

2. Reissner-Mindlin zig-zag model for laminated plates The theoretical model considered in the research is based on the Reissner-Mindlin theory for laminated plates enriched by the Murakami’s function. This function is required to capture the zig-zag effect that characterizes the mechanical behavior of multilayered sandwich plates [82]. Consequently, the acronym RMZ is introduced to denote the zig-zag Reissner-Mindlin approach. According to the current methodology, the displacement field of a laminated plate is given by U x  x, y, z   u x  x, y   zx  x, y   Fz  z  x  x, y  U y  x, y, z   u y  x, y   z y  x, y   Fz  z  y  x, y  ,

(1)

U z  x, y , z   u z  x, y 

where U x , U y , U z are the three-dimensional displacements along the principal coordinates x, y, z . The Cartesian reference systems is shown in Figure 1, in which a generic sandwich plate is represented. The seven degrees of freedom of the model are defined on the plate middle surface, which is the reference domain of the problems under consideration, and are given by three displacements u x , u y , u z , two rotations  x ,  y , and the extent of the zigzag effect  x ,  y . The Murakami’s function Fz  z  assumes the following aspect [80]:

Fz  z    1

k

z z 2z  k 1 k , zk 1  zk zk 1  zk

(2)

in which zk 1 , zk stand for the upper and lower thickness coordinates of the k -th layer, and the thickness hk of the k -th ply is hk  zk 1  zk . Consequently, the total thickness of a laminated structure made of N L layers is given

5

by h   k L1 hk . If a sandwich structure made of an internal core and two reinforcing skins is considered, one gets N

NL  3 . The strain-displacement non-linear relations associated with the kinematic field (1) can be computed assuming moderate strains and rotations according to the von Kármán hypothesis [12, 20]. As far as the membrane strains are concerned, the following definitions are needed:

 xx

2

2

2

u y 1  u z   y  y 1  U z     Fz     z y y 2  y  y y 2  y  , U y U x U z U z u y u x u z u z        x y x y x y x y

 yy   xy

2

U x 1  U z  u x 1  u z  x  x       Fz    z x x 2  x  x x 2  x  U y



  y x    y  x  z     Fz   x y  y   x

(3)

  

whereas the shear strains are given by

U z U x u z F    x  z  x z x z x .  U U z u F y    z  y  z  y y z y z

 xz   yz

(4)

As highlighted in the book by Reddy [20], it should be specified that the strains in (3) do not denote a complete non-linear approach [103, 104], but include only the nonlinear terms related to U z according to the aforementioned von Kármán hypothesis. The following compact notation could be used for conciseness purposes:

 xx   xx 0  z xx1  Fz  xx 2  yy   yy 0  z yy1  Fz  yy 2  xy   xy 0  z xy1  Fz  xy 2 ,

(5)

Fz  2  xz z F  z  yz 2 z

 xz   xz 0   yz   yz 0

             in which  xx ,  yy ,  xy are the in-plane strains,  xx ,  yy ,  xy are the bending strains (curvatures),  xz ,  yz the 0

0

1

0

1

1

0

0

          transverse shear strains, whereas  xx ,  yy ,  xy ,  xz ,  yz are the strain components associated with the zig-zag 2

2

2

2

2

effect. Their meaning can be easily deduced by comparing relations (5) and (4). It can be noted that the strains  xx ,  yy ,  xy are linear along the thickness direction, while a constant behavior of the shear strains  xz ,  yz is

6

assumed. The Murakami’s function, instead, allows a piecewise continuous displacement field along the plate thickness [80, 81]. Once quantities in (5) are defined, the constitutive equations for the k -th orthotropic layer can be introduced as well, in order to relate the stresses to the corresponding strain components

 xx k   Q11 k  xx  Q12 k  yy  Q16 k  xy  yy k   Q12 k  xx  Q22 k  yy  Q26 k  xy  xy k   Q16 k  xx  Q26 k  yy  Q66 k  xy , k 

k 

(6)

k 

 xz  Q44  xz  Q45  yz  yz k   Q45 k  xz  Q55 k  yz where  xx k  ,  yy k  ,  xy k  are the membrane stresses, whereas  xz k  ,  yz k  denote the shear stresses. The mechanical properties of the materials are included within the terms Qij k  , which represent the plane stress-reduced stiffness coefficients in the Cartesian reference system. It should be recalled that a generic orientation 

k

can be assigned

to the layer and the proper coordinate changes are required to compute the stiffnesses Qij k  , as shown in the book by Reddy [20]. These quantities are strictly related to engineering constants of the layer, which are the Young’s k  k  moduli E11 k  , E 22 , the Poisson’s ratio  12k  , and the shear moduli G12 k  , G13 k  , G 23 , as illustrated in [20, 79].

Their definition will be presented in the following Section, since a multiphase approach is developed to this aim. The governing equations are deduced by means of the principle of virtual work [20], assuming that the only applied forces are the in-plane loads Nˆ xx , Nˆ yy . Seven partial differential equations are obtained in terms of stress resultants. In particular, the translational equilibrium along the principal coordinates x, y, z provides   0 N xx  N xy  0 x y 0

N xy  0

x

  N yy 0



y

,

0

(7)

  0 Tx  Ty u    0 u z       Nˆ xx z    Nˆ yy 0 x y x  x  y  y  0

where Nxx 0 , N yy 0 , Nxy 0 are the in-plane force resultants, Tx 0 , Ty 0 the transverse force resultants. On the other hand, the rotational equilibrium is satisfied by the following relations:

M xx M xy   Tx 0  0 x y M xy x



M yy y

7

 0

 Ty  0

,

(8)

in which Mxx , M yy , Mxy represent the moment resultants. Finally, the equilibrium linked to the zig-zag additional degrees of freedom is accomplished as follows  1 N xx  N xy 1   Tx   0 x y 1

N xy 

 N yy

1

x

1



,

(9)

1

 Ty  0

y

where N xx1 , N yy1 , N xy1 and Tx1 , Ty1 are higher-order stress resultants associated with the zig-zag effect. The stress resultants just introduced can be conveniently written as a function of the aforementioned strain components. In particular, the in-plane force resultants are given by

N xx 0  A11 xx 0  A12 yy 0  A16 xy 0  B11 xx1  B12 yy1  B16 xy1  F11 xx 2  F12 yy 2  F16 xy 2 N yy 0  A12 xx 0  A22 yy 0  A26 xy 0  B12 xx1  B22 yy1  B26 xy1  F12 xx 2  F22 yy 2  F26 xy 2 .  0

 0

 2

1

1

1

 0

 0

 2

(10)

 2

N xy  A16 xx  A26 yy  A66 xy  B16 xx  B26 yy  B66 xy  F16 xx  F26 yy  F66 xy On the other hand, the moment resultants can be written as follows

     M xx  B11 xx   B12 yy  B16 xy   D11 xx   D12 yy  D16 xy   G11 xx   G12 yy  G16 xy  0

0

0

1

1

1

2

2

2

     M yy  B12 xx   B22 yy  B26 xy   D12 xx   D22 yy  D26 xy   G12 xx   G22 yy  G26 xy  . 0

0

0

1

1

1

2

2

2

(11)

     M xy  B16 xx   B26 yy  B66 xy   D16 xx   D26 yy  D66 xy   G16 xx   G26 yy  G66 xy  0

0

0

1

1

1

2

2

2

Finally, the transverse force resultants assume the following aspect:

 T    A

Tx    A44  xz   A45 yz   H 44  xz   H 45 yz  0

0

0

2

2

0

0 xz

0 yz

2 xz

2 yz

y

45



 A55

 H 45

 H 55

, 

(12)

where,  is the shear correction factor assumed equal to 5 6 . As far as the higher-order terms related to the zigzag effect are concerned, the following definitions should be taken into account for the corresponding stress resultants. The forces linked to the in-plane stresses are given by  0 1  2  F16 xy 0  G11 xx1  G12 yy  G16 xy1  L11 xx 2  L12 yy  L16 xy 2 N xx1  F11 xx 0  F12 yy  0 1  2  F26 xy 0  G12 xx1  G22 yy  G26 xy1  L12 xx 2  L22 yy  L26 xy 2 , N yy1  F12 xx 0  F22 yy 1

 0

 0

 0

1

1

1

 2

 2

(13)

 2

N xy  F16 xx  F26 yy  F66 xy  G16 xx  G26 yy  G66 xy  L16 xx  L26 yy  L66 xy whereas the ones associated with the shear stresses are defined as

 T     H

.    P    P  

Tx1   H 44 xz 0   H 45 yz 0   P44 xz 2   P45 yz 2  1

y

0 45  xz  H 55

8

0 yz

45

2 xz

55

2 yz

(14)

Definitions (10)-(14) require clearly the following terms: the extensional, coupling and bending stiffnesses Aij , Bij , Dij can be computed as follows, respectively N L zk 1

 A , B , D     Q  1, z, z  dz ij

ij

ij

k

2

ij

k 1 zk

.

(15)

The same quantities Fij , Gij , Lij related to the zig-zag effect, instead, are given by N L zk 1

 F , G , L     Q   F , zF , F  dz . ij

ij

ij

k 1 zk

k

ij

z

z

2 z

(16)

Finally, the additional stiffnesses H ij , Pij assume the following definitions: N L z k 1

 H ij , Pij   



k 1 z k

 F k Qij   z   z

2  F   ,  z   dz .  z  

(17)

At this point, the governing equations (7)-(9) could be written in terms of displacements having in mind the constitutive relations (10)-(14) and the definitions of the strain components that can be deducted from relations (4)-(5). These equations are solved numerically by means of the Finite Element (FE) method, once the corresponding weak formulation is developed and the spatial approximation of the degrees of freedom is introduced [77]. These details are summarized in the Appendix. Mathematically speaking, the fully discretized FE model for the buckling analysis of the whole structure is given by Kd  Gd  0 ,

(18)

where K is the stiffness matrix, whereas G is the geometric stiffness matrix. On the other hand, the vector d collects the discrete nodal values of the displacements. The parameter  stands for the critical increment of the in-plane compressive loads Nˆ xx , Nˆ yy which identify the global buckling of the structure. Relation (18) represents a generalized eigenvalue problem that can be solved numerically once the proper boundary conditions are enforced 0 in terms of displacements only, since a C requirement is needed by the FE formulation. Further details can be

found in the recent papers [92, 98]. In particular, fully clamped and simply supported edges are considered. The former is characterized by null displacements along each boundary edge, whereas the latter can be modeled assuming u x  u z   x   x  0 , for y  0, L y and 0  x  Lx , and u y  u z   y   0  y  L y , where Lx , Ly denote the size of the plate along x , y , respectively.

3. Mechanical characterization of the composite plate

9

y

 0 , for x  0, Lx and

The composite plates under consideration are made of three layers, in which the central one is characterized by a honeycomb structure, as shown in Figure 1. Its mechanical characterization is provided in terms of engineering constants as mentioned in the previous Section. On the other hand, the two external skins are modeled through a multi-scale approach, since they are made of a polymer matrix reinforced at the nano-scale level by randomly oriented carbon nanotubes (CNTs), and at the micro-scale by oriented straight glass fibers. A general framework which describes in detail such approach can be found in [91] for completeness purposes. The Eshelby-Mori-Tanaka scheme is applied at the nanoscale level [100, 101], whereas the Hahn approach is employed as homogenization technique to obtain the overall mechanical properties of these layers in terms of engineering constants [102]. The complete procedure is shown in detail in the previous papers by the authors [92, 98], but the main aspects are summarized below for completeness purposes. Based on the results presented by Shi et al. [95], CNTs are randomly oriented in the polymer matrix (epoxy resin). Consequently, this enriched matrix has isotropic features and its mechanical properties are given by two independent parameters. Once the Young’s modulus EM and the Poisson’s ratio  M of the matrix are known (Table 1), it is possible to obtain the corresponding bulk modulus KM and shear modulus GM as follows

KM 

EM , 3 1  2 M 

GM 

EM . 2 1   M 

(19)

These quantities are affected by the presence of CNTs, whose mass fraction is denoted by wC . The reinforcing nano-fibers are characterized by the following volume fraction VC 1

    VC   C  C  1 ,  wC M M 

(20)

in which C ,  M stand for the densities of the CNTs and the matrix, respectively. In the current research, the values of C  1400 kg m3 and  M  1200 kg m3 are taken into account to this aim. The mechanical properties of the enriched matrix can be computed as shown below in terms of bulk KM* and shear moduli GM* K M*  K M 

VC  C  3K M  C  3 VM  VC C 

,

GM*  GM 

VC C  2GM  C  2 VM  VC  C 

,

(21)

recalling that VM  VC  1 . Once these quantities are evaluated, the Young’s modulus EM* and the Poisson’s ratio  M* of this composite are given by

10

EM* 

9 K M* GM* , 3K M*  GM*

 M* 

3K M*  2GM* , 6 K M*  2GM*

(22)

The mechanical properties of the CNTs are included in the following quantities, which are required to evaluate the parameters  C ,  C ,  C ,C in (21):

C 

3  K M  GM   k C  lC

,

3  GM  kC 

(23)

 2  GM  6 K M  8GM   4GM 1  4GM  2kC  lC  ,   GM  pC GM  3K M  GM   mC  3K M  7GM   5  3  GM  kC   

C  

1 3

 C   nC  2lC 

 2kC  lC  3K M  GM  lC  

 , 

GM  kC

(24)

(25)

 8mCGM  3KM  4GM  2  k  l  2G  l  8G p 1 2  nC  lC   M C  C C M C   , (26) GM  pC 3KM  mC  GM   GPM  7mC  GM   3 GM  kC  5 3

C   

where kC , lC , mC , nC , pC are the Hill’s elastic moduli of the single fiber of CNT modeled as a transversely isotropic material, as illustrated in the paper by Odegard et al. [94]. Their values are listed in Table 2 for a single-walled CNT, characterized by a chiral index of 10 and armchair structure [91]. Finally, the mass density of this isotropic composite can be evaluated as  M*    C   M  VC   M . At this point, the effect of the straight reinforcing fibers can be included at the composite by means of the homogenization procedure proposed by Hahn [102]. In the current paper, this reinforcing phase is made of glass fibers, whose mechanical properties are listed in Table 1, in terms of Young’s modulus EF , shear modulus GF , and Poisson’s ratio  F . Their bulk modulus KF is also required and can be evaluated using the first definition shown in (19). It can be observed that this constituent is modeled as isotropic, as well. The volume fraction of the fibers depends on their mass fraction wF , and assumes the following aspect:  F   VF    *F  1 * w   M  F M 

1

,

(27)

where  F stands for the mass per unit volume of the fibers (  F  2450 kg m 3 ). As mentioned above, the relation VM*  VF  1 is valid also at this level. The homogenization procedure can be now applied. First of all, this approach

requires the plane strain bulk modulus KT of the composite given by

11

KT 

VF   KVM* , K F   KVM* K M*

VF

(28)

in which the parameter  K is needed

K 

1  GM* K F



2 1   M*

.



(29)

The overall engineering constants of the k -th three-phase composite layer of the structure can be evaluated by applying the following definitions:

E11   E FVF  EM* VM* ,

(30)

12 k    FVF  M* VM* ,

(31)

k

G12   G13   k

k

  G23  k

VF  1VM* , VF GF  1VM* GM*

(32)

VF   2VM* , VF GF   2VM* GM*

(33)

where the quantities 1 ,  2 are required

1 

1  GM* GF , 2

3  4 M*  GM* GF

2 



4 1   M*

(34)



It can be easily observed that the composite layer is transversely isotropic, since five independent mechanical properties are introduced, which are the longitudinal Young’s modulus E11  , the shear moduli G12   G13  , G23  , the k

k

k

k

Poisson’s ratio 12  , and the plane strain bulk modulus KT  KT  . The transverse Young’s modulus can be evaluated k

k

as follows    E 22

  4 E11  K T G23 k   k  k  k   G23  E11  4 K T  12  k

k

E11  K T k

k

k

k

 

2

  

.

(35)

Further details concerning this approach can be found in [91]. It should be recalled that several alternatives are available in the literature and can be applied for the same purpose [91, 92, 98].

4. Numerical applications The results of the numerical tests obtained through the application of the current methodology are shown in this Section. Firstly, the proposed approach is validated by means of the comparison with the results available in the literature. Then, several analyses are carried out to investigate the critical loads of global buckling for various

12

sandwich plates with a honeycomb core. In all the following computations, the plates are discretized by using a regular mesh made of ten finite elements along each side.

4.1 Validation of the procedure The first test aims to validate the numerical approach and the theoretical model with respect to the buckling response of a square clamped plate with Lx  L y  946.66 mm . In this circumstance, the external skins are isotropic ( v    v    0.3 1

3

and E 1  E  3   68.2580972 GPa ) and characterized by

h1  h3  0.305 mm . The

2 honeycomb core, instead, is described only by its shear rigidities, since G13 2   G 23  24.132 M Pa , and its

thickness is given by h2  6.198 mm . The critical buckling load N x ,cr applied in the x direction obtained by means of the current approach (RMZ theory and FE method), is compared with the numerical solution by Yuan and Dawe [61], the theoretical solution by Thurston [49], and the experimental output provided by Boller [48]. As it can be observed from the values in Table 3, the present solution is in good agreement with all these references. The accuracy of the current methodology is proven especially by the agreement with the experimental results. Then, sandwich plates with orthotropic oriented skins are considered with the same purpose. In particular, the next application aims to verify the need of the Murakami’s function when these peculiar mechanical configurations are analyzed. A simply supported square plate is considered now with Lx  L y  225 mm . The mechanical properties of the layers are listed in Table 4 and are taken from the paper by Vescovini et al. [69]. All the engineering constants of the materials are specified even if they are not required by the current RMZ theory. In particular, the skins are made of “Material 1” with h1  h3  0.2 mm and characterized by a variable orientation          , whereas 1

3

the “Material 2” is used for the honeycomb core with h2  10 mm . The critical buckling load N x ,cr ( N mm ) is shown in Table 5 for different values of    0 , 90   . It can be observed that the results associated with the RMZ theory are in good agreement with the reference solution provided by the B-Spline Finite Strip Method by Yuan and Dawe [61] and the corresponding percentage difference ( % diff ) is lower than 3.7% (from 1.87% to 3.61% depending on the orthotropic angle). On the other hand, the percentage difference is noticeably higher if the well-known RM theory is employed (from 13.25% to 19.75%). Consequently, these results should suggest the use of the zig-zag theory to obtain an accurate evaluation of the critical buckling load for sandwich structures with

13

a honeycomb core. The same results are presented in graphical form in the paper by Vescovini et al. [69] obtained by the Ritz method, which are in agreement with the ones shown in Table 5. As far as the last comparison test is concerned, a simply supported square plate is taken into account with Lx  L y  225 mm . The core is made of “Material 3” (Table 4) and it is characterized by variable thickness values

h2 . Two different orientations are considered for the orthotropic skins made of “Material 1” (Table 4) and thickness h1  h3  0.2 mm , which are           0 and           40 . The buckling loads N x ,cr ( 1

3

1

3

N mm ) are shown in Figure 2 for different values of h2 . In the same graph, the results provided by Yuan and Dawe [61] and by Vescovini et al. [69] are also shown. A good agreement is observed between the present solution and the reference ones, as far as the global buckling response is concerned. It should be recalled that the reference solutions are obtained through more complex models, whereas the present results are carried out by means of the application of the RMZ theory. A more refined mesh is required to capture the wrinkling phenomenon of these anisotropic sandwich plates, which is represented by the descending curves in Figure 2 as clearly illustrated in the paper by Vescovini et al. [69]. Such phenomenon is activated only for higher ratios of h2 Lx . Therefore, as emphasized in the introduction, the current approach aims to provide an accurate and fast evaluation of the global buckling response by using a linear plate model. A complete discussion about the wrinkling phenomenon is presented in the paper by Vescovini et al. [69].

4.2 Buckling analysis of three-phase composite plates The multiscale approach previously presented is employed in this section to compute the mechanical properties of the external skins of the plates. These layers have the same thickness h1  h3  1mm and are made of an epoxyresin (Table 1) reinforced by randomly oriented CNTs (Table 2) at the nano-scale level, characterized by variable mass fraction wC . This enhanced matrix is also reinforced by straight glass fibers assuming wF  0.85 as mass fraction, whose properties are listed in Table 1. Their orientations  1  ,   3  are assumed as variable and depend on the parameter    0 , 90   . The honeycomb core, instead, is made of “Material 3” (Table 4) and is characterized by h2  45 mm . The plate is defined by Lx  L y  1m . Three stacking sequences are considered, which are a symmetric scheme  / core /   , an antisymmetric scheme  / core /   , and a generic one



/ core / 90     . As far as the external loads are concerned, both uniaxial ( Nˆ xx or Nˆ yy ) and biaxial ( Nˆ xx  Nˆ yy

14

applied simultaneously) configurations are investigated. Each analysis is performed for clamped and simply supported boundary conditions. The main aim of the following analyses is the evaluation of the critical buckling loads N x ,cr or N y ,cr in the case of uniaxial loads applied along x or y , respectively, or Ncr for the biaxial case. The following graphs, one for each value of wC   0, 0.20  , will investigate their variation as a function of the angular parameter  , which is representative of the orthotropic angle of the external skins. Figures 3-5 are associated with the clamped plates. In particular, Figure 3 is related to the critical buckling loads N x ,cr applied in the x direction. On the other hand, the variations of N y ,cr along the y coordinate are depicted

in Figure 4. Finally, Figure 5 shows the variations of the biaxial critical loads Ncr . The following observation can be deduced: -

The symmetric and antisymmetric lamination schemes are more affected by the variation of the parameter  , if compared to the corresponding results obtained for the third stacking sequence.

-

In general, the highest value of the critical buckling load is reached when the reinforcing straight fibers are in the same direction of the applied load, as expected. This circumstance is verified for the symmetric and antisymmetric lamination sequences, if the uniaxial load case is analyzed (Figures 3 and 4). Conversely, the lowest values can be obtained when the fibers are perpendicular to the external load.

-

The variation of the critical buckling loads as a function of  is more emphasized for lower values of

wC . This aspect is due to the higher level of orthotropy of the external skins. In other words, the difference between the longitudinal and transverse Young’s moduli is higher, and the variation of the orthotropic angle is more influential. If the value of wC is greater, instead, the polymer matrix has more enhanced and predominant features, and the effect of the fibers is reduced. This aspect is really evident for

wC  0.10  0.20 , where the graphs related to the symmetric and antisymmetric schemes tend to the one of the generic sequence and the variation of the critical loads are extremely contained. Therefore, the solutions slightly depend on the angle  . In addition, the three curves tend to overlap by increasing the CNT mass fraction. Nevertheless, higher values of wC allow to increase noticeably the critical buckling loads. -

The critical load variation associated with the third lamination scheme is accomplished in a limited interval, approximatively about the average value of the critical loads related to the other sequences.

15

-

The results related to the investigation of N x ,cr and N y ,cr are not symmetrical due to the mechanical properties of the honeycomb core, which are not equal in the corresponding directions. In particular, the highest critical load is reached if the in-plane load is applied in the x direction.

The same analyses are performed for simply supported plates, following the previous approach (Figures 6-8). The observations listed below can be gathered from the corresponding graphs: -

With respect to the clamped configurations, these ones are characterized by more symmetrical output independently from the lamination scheme. The axis of symmetry in the graphs is approximatively located at   45 . In addition, Figures 6-7, which are related to N x ,cr and N y ,cr respectively, are quite similar. Therefore, the simply supported boundary conditions reduce the effect of fiber orthogonality.

-

In general, the maximum critical load is obtained for an angular parameter about   45 . This situation is overturned for wC  0.20 , for which it defines the minimum critical load (however, it should be noted that for higher values of wC the range in which the critical loads are defined is extremely limited, and there is no such a big variation).

-

The generic lamination scheme is characterized by a wider range of variability, if compared to the other two sequences. This difference is gradually reduced by increasing the value of wC , since the orthotropic effect is minimized as explained before.

-

As illustrated for the previous case, the variation of the critical buckling loads as a function of  is more evident for lower values of wC . In addition, the three curves related to the different stacking sequences tend to move closer by increasing the CNT mass fraction. Finally, higher values of wC allow to increase the critical buckling load.

Final remarks The proposed approach based on the Reissner-Mindlin theory for laminated plates embedding the Murakami’s function for the zig-zag effect has proven to be an accurate and easily applicable tool for the analysis of the global buckling response of sandwich plates with a honeycomb core and three-phase orthotropic skins, as shown by the comparison with the experimental results available in the literature.

16

The results presented in this paper have proven that the boundary conditions, the lamination schemes, and the fiber orientation have a great influence with respect to the buckling response of the structures under consideration. In addition, the role of CNTs used as nano reinforcing phase is very influential since the proper choice of their mass fraction allows to increase the value of critical loads and reduce the orthotropy of the external skins. In fact, the orthotropic effect can be reduced for higher values of the CNT mass fraction and the external skins tend to show isotropic features. As a consequence, in addition to the possibility to reach higher critical buckling loads, the variation range in terms of critical load is noticeably limited since the directional effects tend to be negligible. The investigations carried out in the research could provide useful information for the design and application of such panels in many engineering fields, in order to optimize and maximize their critical buckling loads. In particular, the mass fraction of CNTs represents a fundamental design parameter in the manufacturing of these sandwich plates.

Appendix The displacements are approximated by means of the quadratic Lagrange interpolation functions for a nine-node rectangular element. The displacements are all approximated by using the same degrees of interpolation function as shown below m

m

u x   u x , j j ,

u y   u y , j j ,

j 1

m

 x    x , j j , j 1

j 1

m

u z   u z , j j , j 1

m

m

m

j 1

j 1

j 1

 y    y , j j ,  x   x , j j ,  y   y , j j

(36)

where  j stands for the j -th Lagrange polynomial. The discrete values of the displacements at the j -th node are indicated by u x , j , u y , j , u z , j ,  x , j ,  y , j ,  x , j ,  y , j . This approximation is performed in each finite element in which the structure is subdivided. Once the governing equations are obtained, the corresponding weak form can be easily developed as illustrated in the book by Reddy [20]. The discrete equation (18) should be first written for the generic element e as follows

K ede  G ede  0 ,

(37)

e in which the various matrices have the same meaning illustrated before. The vector d collects the nodal

displacements of the element as shown below d e  u x , j

u y, j

T

u z , j x , j  y , j  x , j  y , j  ,

17

(38)

for j  1,..., m , where m is the number of nodes of the e -th element. The element stiffness matrix, instead, is given by

 K11   K17        Ke   ,       71    K 77  K

(39)

where the generic element K ij , for i, j  1,..., m and  ,   1,...,7 , assume the following aspect:

   j   j   j   j  K ij11     A11 i  A16 i  A16 i  A66 i  dxdy x x x y y x y y  x y    j   j   j   j  K ij12     A12 i  A16 i  A26 i  A66 i  dxdy x y x x y y y x  x y K ij13  0    j    j   j   j  B66 i K ij14     B11 i  B16 i  B16 i  dxdy x  y y  x x x y y      x y

(40)

   j   j   j   j  K ij15     B12 i  B16 i  B26 i  B66 i  dxdy x y x x y y y x  x y    j   j   j   j  K ij16     F11 i  F16 i  F16 i  F66 i  dxdy x x x y y x y y  x y    j   j   j    j  F66 i  F26 i  F16 i K ij17     F12 i  dxdy x x y y y x  x y x y    j   j   j   j  K ij21     A16 i  A66 i  A12 i  A26 i  dxdy x x x y y x       y y  x y    j   j   j   j  K ij22     A26 i  A66 i  A22 i  A26 i  dxdy x y x x y y y x  x y K ij23  0    j    j   j   j  B26 i K ij24     B16 i  B66 i  B12 i  dxdy x  y y  x x x y y      x y    j   j   j   j  K ij25     B26 i  B66 i  B22 i  B26 i  dxdy x y x x y y y x  x y    j   j   j   j  K ij26     F16 i  F66 i  F12 i  F26 i  dxdy x x x y y x y y  x y    j   j   j    j  F66 i  F22 i  F26 i K ij27     F26 i  dxdy x x y y y x  x y x y

18

(41)

K ij31  0 K ij32  0    j   j   j   j  K ij33      A44 i  A45 i  A45 i  A55 i  dxdy x x x y y x y y   x y     K ij34      A44 i  j  A45 i  j  dxdy y x   x y

(42)

    K ij35      A45 i  j  A55 i  j  dxdy x y   x y     K ij36      H 44 i  j  H 45 i  j  dxdy x y   x y     K ij37      H 45 i  j  H 55 i  j  dxdy y x     x y    j   j   j   j   B16 i  B16 i  B66 i K ij41     B11 i  dxdy       y y  x x x y y x x y    j   j   j   j   B16 i  B26 i  B66 i K ij42     B12 i  dxdy x y x x y y y x  x y  j  j    A45i K ij43      A44i  dxdy x y   x y     j   j   j   j K ij44     D11 i  D16 i  D16 i  D66 i   A44i j  dxdy         x x x y y x y y  x y     j   j   j   j  D16 i  D26 i  D66 i   A45i j  dxdy K ij45     D12 i         x y x x y y y x  x y     j   j   j   j  G16 i  G16 i K ij46     G11 i  G66 i   H 44i j  dxdy x x x y y x y y  x y     j   j   j   j  D16 i  D26 i  D66 i   H 45i j  dxdy K ij47     D12 i x y x x y y y x  x y

19

(43)

   j   j   j     B66 i  B12 i  B26 i j  dxdy K ij51     B16 i       y y  x x x y y x x y    j   j   j   j   B66 i  B22 i  B26 i K ij52     B26 i  dxdy x y x x y y y x  x y  j  j    A55i K ij53      A45i  dxdy x y   x y     j   j   j   j K ij54     D16 i  D66 i  D12 i  D26 i   A45i j  dxdy x x x y y x y y  x y

(44)

    j   j   j   j  D66 i  D22 i  D26 i   A55i j  dxdy K ij55     D26 i         x y x x y y y x  x y     j   j   j   j  F66 i  F12 i K ij56     F16 i  F26 i   H 45i j  dxdy x x x y y x y y  x y     j   j   j   j  F66 i  F22 i  F26 i   H 55i j  dxdy K ij57     F26 i x y x x y y y x  x y    j   j   j   j   F16 i  F16 i  F66 i K ij61     F11 i  dxdy x x x y y x y y  x y    j   j   j   j   F16 i  F26 i  F66 i K ij62     F12 i  dxdy       y x  x y x x y y x y  j  j    H 45i K ij63      H 44i  dxdy x y   x y     j   j   j   j K ij64     G11 i  G16 i  G16 i  G66 i   H 44i j  dxdy x x x y y x y y  x y     j   j   j   j  G16 i  G26 i  G66 i   H 45i j  dxdy K ij65     G12 i x y x x y y y x  x y     j   j   j   j  L16 i  L16 i K ij66     L11 i  L66 i   P44i j  dxdy     y x y y x x x y      x y     j   j   j   j K ij67     L12 i  L16 i  L26 i  L66 i   P45i j  dxdy x y x x y y y x  x y

20

(45)

   j   j   j   j   F66 i  F12 i  F26 i K ij71     F16 i  dxdy       y y  x x x y y x x y    j   j   j   j   F66 i  F22 i  F26 i K ij72     F26 i  dxdy x y x x y y y x  x y  j  j    H 55i K ij73      H 45i  dxdy x y   x y     j   j   j   j K ij74     G16 i  G66 i  G12 i  G26 i   H 45i j  dxdy . x x x y y x y y  x y

(46)

    j   j   j   j  G66 i  G22 i  G26 i   H 55i j  dxdy K ij75     G26 i         x y x x y y y x  x y     j   j   j   j  L66 i  L12 i K ij76     L16 i  L26 i   P45i j  dxdy x x x y y x y y  x y     j   j   j   j  L66 i  L22 i  L26 i   P55i j  dxdy K ij77     L26 i x y x x y y y x  x y e On the other hand, the geometric stiffness matrix G is characterized only by the following non null term:

   j    j Gij33     Nˆ xx i  Nˆ yy i  dxdy . y y  x x x y

(47)

It should be recalled that the integrals which appear in the definitions (40)-(47) are solved numerically by means of the Gauss-Legendre quadrature rule, as specified in [77]. This procedure is accomplished by using the isoparametric formulation to facilitate the numerical computation of integrals [76, 77]. The elements of the stiffness matrix that involve the shear stresses are evaluated through the reduced integration scheme in order to avoid the shear locking issue. Analogously, the reduced integration is used also for the evaluation of the geometric stiffness matrix. Once the fundamental matrices are computed for the e -th element, the assembly procedure 0 provides the global operators of the discrete structure ensuring the C compatibility requirements among the

element interfaces. Equation (18) is obtained as a result.

21

References 1.

Brush D.O., Almroth B.O., Buckling of Bars, Plates and Shells, McGraw-HiII, New York, 1975.

2.

Åkesson B., Plate Buckling in Bridges and other Structures, Taylor & Francis Group, London, UK, 2007.

3.

Kubiak T., Static and Dynamic Buckling of Thin-Walled Plate Structures, Springer International Publishing, Switzerland, 2013.

4.

Tarantino A.M., Homogeneous equilibrium configurations of a hyperelastic compressible cube under equitriaxial dead-load tractions, J. Elasticity, 2008, 92, 227-254.

5.

Tarantino A.M., Equilibrium paths of a hyperelastic body under progressive damage, J. Elasticity,2014, 114, 225-250.

6.

Lanzoni L., Tarantino A.M., Equilibrium configurations and stability of a damaged body under uniaxial tractions, Z. Angew. Math. Phys., 2015, 66, 171-190.

7.

Lanzoni L., Tarantino A.M., A simple nonlinear model to simulate the localized necking and neck propagation, Int. J. Nonlin. Mech., 2016, 84, 94-104.

8.

Euler L., Additamentum I de curvis elasticis, methodus inveniendi lineas curvas maximi minimivi proprietate gaudentes, Opera omnia I, 1744, 24, 232-297.

9.

Poincaré H., Sur l’equilibre d’une masses fluide animés d’un mouvement de rotation, Acta Math., 1885, 7, 259- 380.

10. Schmidt E., Zur Theorieder Linearen und nichtlinearen Integralgleichungen, Math. Ann., 1908, 65, 370-379. 11. Liapunov A., Probleme Général de la Stabilité du Mouvement, Princeton University Press, 1947. 12. Von Kármán T., Festigkeitsprobleme im maschinenbau, Encyklopädie der Mathematischen Wissenschaften, 1910, 311-385. 13. Koiter W.T., The stability of elastic equilibrium, Stanford Univ. Ca Dept. of Aeronautics and Astronautics, 1970. 14. Bloom F., Coffin D., Handbook of Thin Plate Buckling and Postbuckling, Chapman & Hall/CRC, Boca Raton, Florida, USA, 2001. 15. Bert C.W., Francis P.H., Composite Material Mechanics: Structural Mechanics, AIAA J., 1974, 12, 11731186. 16. Librescu L., Elastostatics and Kinetics of Anisotropic and Heterogeneous Shell-Type Structures, Springer, 1975.

22

17. Reddy J.N., Energy and Variational Methods in Applied Mechanics, John Wiley, 1984. 18. Leissa A.W., Buckling of Laminated Composite Plates and Shell Panels, Air Force Wright Aeronautical Laboratories, AFWAL-TR- 85-3069, 1985. 19. Chou T.-W., Microstructural Design of Fiber Composites, Cambridge University Press, New York, New York, USA, 1992. 20. Reddy J.N., Mechanics of Laminated Composite Plates and Shells. Theory and Analysis. Second Edition, CRC Press, 2004. 21. Khdeir A.A., Free vibration and buckling of symmetric cross-ply aminated plates by an exact method, J. Sound Vib., 1988, 126, 447-461. 22. Khdeir A.A., Reddy J.N., Frederick D., A study of bending, vibration and buckling of cross-ply circular cylindrical shells with various shell theories, Int. J. Engng, Sci., 1989, 27, 1337-1351. 23. Touratier M., An Efficient Standard Plate Theory, Int. J. Engng. Sci., 1991, 29, 901-916. 24. Kahya V., Buckling analysis of laminated composite and sandwich beams by the finite element method, Compos. Part B.-Eng., 2016, 91, 126-134. 25. Vescovini R., Bisagni C., Buckling analysis and optimization of stiffened composite flat and curved panels, AIAA J., 2012, 50, 904-915. 26. Vescovini R., Dozio L., Exact refined buckling solutions for laminated plates under uniaxial and biaxial loads, Compos. Struct., 2015, 127, 356-368. 27. Vescovini R., Dozio L., D’Ottavio M., Polit O., On the application of the Ritz method to free vibration and buckling analysis of highly anisotropic plates, Compos. Struct., 2018, 192, 460-474. 28. Ruocco E., Reddy J.N., A closed-form solution for buckling analysis of orthotropic Reddy plates and prismatic plate structures, Compos. Part B-Eng., 2019, 169, 258-273. 29. Hemp W.S., On a theory of sandwich construction, The College of Aeronautics Cranfield, Report No. 15, 1948. 30. Plantema F.J., Sandwich Construction, Wiley, New York, 1966. 31. Allen H.G., Analysis and Design of Structural Sandwich Panels. Pergamon Press, 1969. 32. Vinson J.R., The Behavior of Sandwich Structures of Isotropic and Composite Material, Technomic Publishing Company, Inc., PA, USA, 1999.

23

33. Wempner G.A., Baylor J.L., General theory of sandwich plates with dissimilar facings, Int. J. Solids Struct., 1965, 1, 157-177. 34. Bazant Z.P., Beghini A., Stability and finite strain of homogenized structures soft in shear: Sandwich or fiber composites, and layered bodies, Int. J. Solids Struct., 2006, 43, 1571-1593. 35. Ji W., Waas A.M., Accurate buckling load calculations of a thick orthotropic sandwich panel, Compos. Sci. Technol., 2012, 72, 1134-1139. 36. Moleiro F., Mota Soares C.M., Carrera E., Three-dimensional exact hygro-thermo-elastic solutions for multilayered plates: Composite laminates, fibre metal laminates and sandwich plates, Compos. Struct., 2019, 216, 260-278. 37. Guo K., Zhu L., Li Y., Yu T.X., Numerical study on mechanical behavior of foam core sandwich plates under repeated impact loadings, Compos. Struct., 2019, 224, 111030. 38. Yuan L., Batra R.C., Optimum first failure load design of one/two-core sandwich plates under blast loads, and their ultimate loads, Compos. Struct., 2019, 224, 111022. 39. Zhang J., Ashby M.F., Buckling of Honeycombs Under In-Plane Biaxial Stresses, Int. J. Mech. Sci., 1992, 34, 491-509. 40. Zenkour A.M., A comprehensive analysis of functionally graded sandwich plates: Part 2-Buckling and free vibration, Int. J. Solids Struct., 2005, 42, 5243-5258. 41. Guo J., Sun T., Pan E., Three-dimensional nonlocal buckling of composite nanoplates with coated onedimensional quasicrystal in an elastic medium, Int. J. Solids Struct., 2019. In press. 42. Waddar S., Pitchaimani J., Doddamani M., Barbero E., Buckling and vibration behaviour of syntactic foam core sandwich beam with natural fiber composite facings under axial compressive loads, Compos. Part B.Eng., 2019, 175, 107133. 43. Fan H., Yang L., Sun F., Fang D., Compression and bending performances of carbon fiber reinforced latticecore sandwich composites, Compos. Part A, 2013, 52, 118-125. 44. Li Z., Wang B.L., Guo S.L., Li J.E., Thermal shock resistance of ceramic foam sandwich structures: Theoretical calculation and finite element simulation, Int. J. Solids Struct., 2019, 176-177, 108-120. 45. Catapano A., Montemurro M., A multi-scale approach for the optimum design of sandwich plates with honeycomb core. Part I: homogenisation of core properties, Compos. Struct., 2014, 118, 664-676.

24

46. Catapano A., Montemurro M., A multi-scale approach for the optimum design of sandwich plates with honeycomb core. Part II: the optimisation strategy, Compos. Struct., 2014, 118, 677-690. 47. Xue X., Zhang C., Chen W., Wu M., Zhao J., Study on the impact resistance of honeycomb sandwich structures under low-velocity/heavy mass, Compos. Struct., 2019, 226, 111223. 48. Boller K.H., Buckling loads of flat sandwich panels in compression, Forest Products Laboratory Report No. 1525-D, 1947. 49. Thurston G.A., Bending and Buckling of Clamped Sandwich Plates, J. Aero. Sci., 1957, 24, 407-412. 50. Robinson J.R., The Buckling and Bending of Orthotopic Sandwich Panels With All Edges Simply-Supported, Aeronaut. Quart., 1955, 6, 125-148. 51. Chang C.C., Ebcioglu I. K., Haight. C.H., General Stability Analysis of Orthotropic Sandwich Panels for four Different Boundary Conditions (Extension of the March‐Ericksen Approach). ZAMM-Z. Angew. Math. Me., 1962, 42, 373-389. 52. Benson A.S., Mayers J., General Instability and Face Wrinkling of Sandwich Plates: Unified Theory and Applications, AIAA J., 1967, 5, 729-739. 53. Pearce T.R.A., Webber J.P.H., Buckling of Sandwich Panels with Laminated Face Plates, Aeronaut. Quart., 1971, 23, 148-160. 54. Pearce T.R.A., Webber J.P.H., Experimental Buckling Loads of Sandwich Panels with Carbon Fibre Faceplates, Aeronaut. Quart., 1973, 24, 295-312. 55. Rao K.M., Buckling Analysis of Anisotropic Sandwich Plates Faced with Fiber-Reinforced Plastics, AIAA J., 1985, 23, 1247-1253. 56. Minguet P., Dugundji J., Lagace P.A., Buckling and Failure of Sandwich Plates with Graphite-Epoxy Faces and Various Cores, J. Aircraft, 1988, 25, 372-379. 57. Yeh W.-N., Wu Y.-E., Enhancement of buckling characteristics for sandwich structure with fiber reinforced composite skins and core made of aluminum honeycomb and polyurethane foam, Theor. Appl. Fract. Mech., 1991, 15, 63-74. 58. Kwon Y.W., Murphy M.C., Castelli V., Buckling of Unbalanced, Sandwich Panels With Titanium and GRP Skins, J Press. Vess.-T. ASME, 1995, 117, 40-44. 59. Frostig Y., Buckling of Sandwich Panels with A Flexible Core-High-Order Theory, Int. J. Solids Struct., 1998, 35, 183-204.

25

60. Dawe D.J., Yuan W.X., Overall and local buckling of sandwich plates with laminated faceplates, Part I: Analysis, Comput. Methods Appl. Mech. Engrg., 2001, 190, 5197-5213. 61. Yuan W.X., Dawe D.J., Overall and local buckling of sandwich plates with laminated faceplates, part II: Applications, Comput. Methods Appl. Mech. Engrg., 2001, 190, 5215-5231. 62. Roberts J.C., Boyle M.P., Wienhold P.D., White G.J., Buckling, collapse and failure analysis of FRP sandwich panels, Compos. Part B-Eng., 2002, 33, 315-324. 63. Davalos J.F., Chen A., Buckling Behavior of Honeycomb FRP Core with Partially Restrained Loaded Edges under Out-of-plane Compression, J. Compos. Mater., 2005, 39, 1465-1485. 64. Kaman M.O., Solmaz M.Y., Turan K., Experimental and Numerical Analysis of Critical Buckling Load of Honeycomb Sandwich Panels, J. Compos. Mater., 2010, 44, 2819-2831. 65. D’Ottavio M., Polit O., Linearized global and local buckling analysis of sandwich struts with a refined quasi3D model, Acta Mech., 2015, 226, 81-101. 66. D’Ottavio M., Polit O., Ji W., Waas A.M., Benchmark solutions and assessment of variable kinematics models for global and local buckling of sandwich struts, Compos. Struct., 2016, 156, 125-134. 67. Muc A., Stawiarski A., Romanowicz P., Experimental Investigations of Compressed Sandwich Composite/Honeycomb Cylindrical Shells, Appl. Compos. Mater., 2018, 25, 177-189. 68. Nazari A.R., Hosseini-Toudeshky H., Kabir M.Z., Experimental investigations on the sandwich composite beams and panels with elastomeric foam core, J. Sandw. Struct. Mater., 2019, 21, 865-94. 69. Vescovini R., D’Ottavio M., Dozio L., Polit O., Buckling and wrinkling of anisotropic sandwich plates, Int. J. Engng. Sci., 2018, 130, 136156. 70. Minahen T.M., Knauss W.G., Creep buckling of viscoelastic structures, Int. J. Solids Struct., 1993, 30, 10751092. 71. Dezi L., Menditto G., Tarantino A.M., Homogeneous structures subjected to successive structural system changes, J. Eng. Mech.-ASCE, 1990, 116, 1723-1732. 72. Dezi L., Tarantino A.M., Time dependent analysis of concrete structures with variable structural system, ACI Mater. J., 1991, 88, 320-324. 73. Dezi L., Menditto G., Tarantino, A.M., Viscoelastic heterogeneous structures with variable structural system, J. Eng. Mech.-ASCE, 1993, 119, 238-250.

26

74. Dezi L., Tarantino A.M., Creep in continuous composite beams. Part II: Parametric study, J. Struct. Eng.ASCE, 1993, 119, 2112-2133. 75. Forcellini D., Tarantino A.M., Assessment of stone columns as a mitigation technique of liquefactioninduced effects during Italian earthquakes (may 2012), Sci. World J., vol. 2014, Article ID 216278, 8 pages. 76. Ferreira A.J.M., MATLAB codes for finite element analysis: solids and structures, Springer Science & Business Media, 2008. 77. Reddy J.N., An Introduction to the Finite Element Method. Third Edition, McGraw Hill, New York, USA, 2006. 78. Fantuzzi N., Tornabene F., Bacciocchi M., Neves A.M.A., Ferreira A.J.M., Stability and accuracy of three Fourier expansion-based strong form finite elements for the free vibration analysis of laminated composite plates, Int. J. Numer. Meth. Eng., 2017, 111, 354-382. 79. Bacciocchi M., Tarantino A.M., Natural Frequency Analysis of Functionally Graded Orthotropic Cross-Ply Plates Based on the Finite Element Method, Math. Comput. Appl., 2019, 24, 52. 80. Murakami H., Laminated composite plate theory with improved in-plane responses, J. Appl. Mech., 1986, 53, 661-666. 81. Carrera E, On the use of the Murakami’s Zig-Zag function in the modeling of layered plates and shells, Comput. Struct., 2004, 82, 541-554. 82. Carrera E., C0 reissner-mindlin multilayered plate elements including Zig-Zag and interlaminar stress continuity, Int. J. Numer. Meth. Eng., 1996, 39, 1797-1820. 83. Carrera E., Developments, ideas and evaluations based upon the Reissner’s mixed theorem in the modeling of multilayered plates and shells, Appl. Mech. Rev., 2001, 54, 301-329. 84. Carrera E., Historical review of Zig-Zag theories for multilayered plates and shells, Appl. Mech. Rev., 2003, 56, 287-308. 85. Carrera E., Theories and finite elements for layered plates and shells: a unified compact formulation with numerical assessment and benchmarking, Arch. Comput. Meth. Eng., 2003, 10, 215-296. 86. Maturi D.A., Ferreira A.J.M., Zenkour A.M., Mashat D.S., Analysis of laminated shells by murakami’s Zig– Zag theory and radial basis functions collocation, J. Appl. Math., 2013, 2013, 14. 87. Filippi M., Carrera E., Valvano S., Analysis of multilayered structures embedding viscoelastic layers by higher-order, and zig-zag plate elements, Compos. Part B-Eng., 2018, 154, 77-89.

27

88. Hu H., Belouettar S., Daya E.M., Potier-Ferry M., Evaluation of kinematic formulations for viscoelastically damped sandwich beam modeling, J. Sandw. Struct. Mater., 2006, 8, 477-495. 89. Tornabene F., Fantuzzi N., Bacciocchi M., Foam core composite sandwich plates and shells with variable stiffness: Effect of the curvilinear fiber path on the modal response, J. Sandw. Struct. Mater., 2019, 21, 320365. 90. Falope F.O., Lanzoni L., Tarantino A.M., Modified hinged beam test on steel fabric reinforced cementitious matrix (SFRCM), Compos. Part B-Eng., 2018, 146, 232-243. 91. Tornabene F., Bacciocchi M., Fantuzzi N., Reddy J.N., Multiscale Approach for Three-Phase CNT/Polymer/Fiber Laminated Nanocomposite Structures, Polym. Compos., 2019, 40, E102-E126. 92. Bacciocchi M., Tarantino A.M., Time-dependent behavior of viscoelastic three-phase composite plates reinforced by Carbon nanotubes, Compos. Struct., 2019, 216, 20-31. 93. Capecchi D., Ruta G., Trovalusci P., Voigt and Poincaré’s mechanistic-energetic approaches to linear elasticity and suggestions for multiscale modelling, Arch. Appl. Mech., 2011, 81, 1573-1584. 94. Odegard G.M., Gates T.S., Wise K.E., Park C., Siochi E.J., Constitutive modeling of nanotube-reinforced polymer composites, Compos. Sci. Technol., 2003, 63, 1671-1687. 95. Shi D.L., Huang Y.Y., Hwang K.C., Gao H., The effect of nanotube waviness and agglomeration on the elastic property of carbon nanotube–reinforced composites, J. Eng. Mater. T. ASME, 2004, 126, 250-257. 96. Civalek Ӧ., Baltacıoğlu A.K., Vibration of carbon nanotube reinforced composite (CNTRC) annular sector plates by discrete singular convolution method, Compos. Struct., 2018, 203, 458-465. 97. Di Sciuva M., Sorrenti M., Bending, free vibration and buckling of functionally graded carbon nanotubereinforced sandwich plates, using the extended Refined Zigzag Theory, Compos. Struct., 2019, 227, 111324. 98. Bacciocchi M., Luciano R., Majorana C., Tarantino A.M., Free Vibrations of Sandwich Plates with Damaged Soft-Core and Non-Uniform Mechanical Properties: Modeling and Finite Element Analysis, Materials, 2019, 12, 2444. 99. Zakaria M.R., Akil H.M., Kudus M.H.A., Ullah F., Javed F., Nosbi N., Hybrid carbon fiber-carbon nanotubes reinforced polymer composites: A review, Compos. Part B-Eng., 2019, 176, 107313. 100. Eshelby J.D., The determination of the elastic field of an ellipsoidal inclusion, and related problems, P. Roy. Soc. Lond. A Mat., 1957, 241, 376-396.

28

101. Mori T., Tanaka K., Average Stress in Matrix and Average Elastic Energy of Materials with Misfitting Inclusions, Acta Metall., 1973, 21, 571-574. 102. Hahn H.T., Simplified Formulas for Elastic Moduli of Unidirectional Continuous Fiber Composites, J. Compos. Technol. Res., 1980, 2, 5-7. 103. Lanzoni L., Tarantino A.M., Finite Anticlastic Bending of Hyperelastic Solids and Beams, J. Elasticity, 2018, 131, 137-170. 104. Falope F.O., Lanzoni L., Tarantino A.M., Bending device and anticlastic surface measurement of solids under large deformations and displacements, Mech. Res. Commun., 2019, 97, 52-56.

29

Figure 1 – Sandwich plate with a honeycomb core and reinforcing skins.

Figure 2 – Variation of the critical buckling load N x ,cr [N/mm] for a simply-supported sandwich square plate with orthotropic skins (   0 and   40 ), as a function of the geometric parameter h2 Lx : comparison with the literature.

30

wC  0.00

wC  0.01

wC  0.02

wC  0.05

wC  0.10

wC  0.20

Figure 3 – Variation of the uniaxial critical buckling load N x ,cr [kN/m] of a fully clamped sandwich square plate as a function of the angular parameter  of the orthotropic skins, for several values of the CNT mass fraction wC .

31

wC  0.00

wC  0.01

wC  0.02

wC  0.05

wC  0.10

wC  0.20

Figure 4 – Variation of the uniaxial critical buckling load N y ,cr [kN/m] of a fully clamped sandwich square plate as a function of the angular parameter  of the orthotropic skins, for several values of the CNT mass fraction wC .

32

wC  0.00

wC  0.01

wC  0.02

wC  0.05

wC  0.10

wC  0.20

Figure 5 – Variation of the biaxial critical buckling load N cr [kN/m] of a fully clamped sandwich square plate as a function of the angular parameter  of the orthotropic skins, for several values of the CNT mass fraction wC .

33

wC  0.00

wC  0.01

wC  0.02

wC  0.05

wC  0.10

wC  0.20

Figure 6 – Variation of the uniaxial critical buckling load N x ,cr [kN/m] of a simply supported clamped sandwich square plate as a function of the angular parameter  of the orthotropic skins, for several values of the CNT mass fraction wC .

34

wC  0.00

wC  0.01

wC  0.02

wC  0.05

wC  0.10

wC  0.20

Figure 7 – Variation of the uniaxial critical buckling load N y ,cr [kN/m] of a simply supported sandwich square plate as a function of the angular parameter  of the orthotropic skins, for several values of the CNT mass fraction wC .

35

wC  0.00

wC  0.01

wC  0.02

wC  0.05

wC  0.10

wC  0.20

Figure 8 – Variation of the biaxial critical buckling load N cr [kN/m] of a fully clamped sandwich square plate as a function of the angular parameter  of the orthotropic skins, for several values of the CNT mass fraction wC .

36

Young’s modulus

Matrix (epoxy resin) EM  3.27GPa

Straight fibers (glass fibers) E F  71G Pa

Shear modulus

-

G F  30 G Pa

Poisson’s ratio

 M  0.38

 F  0.22

Property

Table 1 – Mechanical properties of the layers used in the comparison tests.

kC

271GPa

lC

88GPa

mC

17 GPa

nC

1089 GPa

pC

442 GPa

Table 2 – Hill’s elastic moduli for a single-walled CNT [91].

Methodology

Critical load N x ,cr [N/mm]

Numerical solution [61] Theoretical solution [49] Experimental solution [48] Present solution

43.19 – 43.22 42.4 39.05 – 42.38 41.478

Table 3 – Critical buckling load for a clamped sandwich square plate with isotropic skins: comparison with the literature.

Property

Material 1 [69] (skins, k  1, 3 )

Material 2 [69] (core, k  2 )

Material 3 [69] (core, k  2 )

E11 k  [MPa]

229000

0

0

k

13350

0

0

k

13350



828.0

k

G12 [MPa]

5249

0

0

G13 k  [MPa]

5249

146

146

G23 k  [MPa]

3000

90.40

90.40

 12 k 

0.3151

0.00

0.00

k 

0.3151

0.00

0.00

k 

0.3151

0.00

0.00

E22 [MPa] E33 [MPa]

 13

 23

Table 4 – Mechanical properties of the layers used in the comparison tests.

37

 0° 10° 20° 30° 40° 50° 60° 70° 80° 90°

Yuan and Dawe [61] 424.6 416.5 398.8 382.6 372.8 300.3 255.0 222.7 211.7

Present solution (RMZ) 409.284 402.513 388.161 375.245 365.813 342.532 292.867 248.829 217.449 206.78

Present solution (RM) 480.859 472.112 453.329 439.386 433.255 418.793 359.612 304.041 263.565 249.168

% diff

% diff

(RMZ) 3.61% 3.36% 2.67% 1.92% 1.87% 2.48% 2.42% 2.36% 2.32%

(RM) 13.25% 13.35% 13.67% 14.84% 16.22% 19.75% 19.23% 18.35% 17.70%

Table 5 – Critical buckling loads for a simply-supported sandwich square plate with orthotropic skins varying the orientation  : comparison with the literature.

38