Post-buckling of web-core sandwich plates based on classical continuum mechanics: success and needs for non-classical formulations

The paper investigates the post-buckling response of web-core sandwich plates through classical continuum mechanics assumptions. The compressive loading is assumed to be in the direction of the web plates. Equivalent Single Layer (ESL) plate formulation is used with the kinematics of the First order Shear Deformation Theory (FSDT). During the initial, membrane-dominated loading stages, it is observed that the effect of finite size of the periodic microstructure is barely influences the plate responses. At the higher loads, when bending is activated, the finite size of the microstructure activates secondary shear-induced bending moments at the unit cells of the plate. A method to capture the envelope of the maximum values of these bending moments is presented. The findings are validated with the shell element models of the actual 3D-geometry. Finally, the physical limits of the classical continuum mechanics are discussed in the present context.


Introduction
There is constant need to study new materials and structural configurations for thin-walled structures. Weight reduction and strength increase enable lighter and often more sustainable structural solutions that can be used in civil, naval and aeronautical applications. Structural efficiency in terms of strength-toweight and stiffness-to-weight is obtained by positioning of the material according to the load-carrying mechanism of the structure. In bending dominated applications, this often results in plates with periodic microstructure, for example, single-sided stiffened plates and sandwich panels with visibly discrete core; see Fig. 1. The periodic structure allows integration of functions to the panel(e.g. air-conditioning and cable tracks). Thus, unidirectional stiffening system with coarse spacing becomes attractive structural alternative; Refs. [1][2][3][4][5][6][7][8][9][10][11][12][13]. These unidirectional plates, especially in the form of steel sandwich plates, have been found to have excellent mechanical properties for buckling, bending and impact loads [14][15][16][17][18][19][20][21]. These structures are, however, challenging in terms of structural analysis. The direct inclusion of the microstructure to computational models of the entire structure leads to expensive pre-processing, analysis and post-processing times (i.e., computationally intensive). The problem compounds if the geometry changes and analysis is needed to be repeated several times(e.g., during structural optimization or reliability analysis [18,20]). Further complications are caused by the fact that often the failure of these panels occurs first locally in the microstructural level in the face and web-plates or at the welds by yielding or buckling, see Fig. 1. This means that certain level of accuracy is needed in the analysis of structural details.
Homogenization is the alternative method for a direct modeling of both micro-and macrostructure. In terms of plates, it has been discussed by several authors(see, e.g., Refs. [1,3,4,[11][12][13][14][20][21][22][23][24][25][26][27][28][29][30][31][32][33][34][35][36]). In plates, the main idea is to reduce the mathematical  [18]; full-scale ship collisions Ehlers et al. [37] and quasi-static plate punch experiments Körgesaar et al. [38]) description to Equivalent Single Layer (ESL) representation in which the deformations and stress resultants (forces and moments) are computed for known load and boundary conditions. In periodic sandwich panels, the closed unit cells, including the interacting face and web plates, form a structure that warps in bending due to the out-of-plane shear. With this phenomenon included, Libove and Hubka [1] presented an ESL-theory for corrugated-core sandwich plates which follows First order Shear Deformation Theory (FSDT). They derived the equivalent shear stiffness and provided series solutions to prevailing differential equations for bending. After this work, several papers have been written for other microstructures [26-30, 33, 37, 38] where the stiffness parameters change, but the differential equations remain the same. The approach has been also validated to geometrically moderate non-linear global deformations [26][27][28]. However, the problem arises in the assessment of stresses. As the local, microstructural length-scale can be close to that of macro-structural, the local oscillations of stresses within microstructure can have significant contribution to the overall stress state caused by macroscopic plate behavior. It has been experimentally shown that the out-of-plane shear damages the periodic structures from plate edges or close to point loads, thus the secondary stresses induced by shear need to be carefully assessed [20][21][22], see also Fig. 1. These findings are based on bending of plates. This situation can also occur during in-plane compression after buckling when the initially membrane-dominated responses are complemented with out-of-plane responses [39].
The objective of this study is to present the phenomena related to post-buckling response of web-core sandwich plates and thereby extending the work presented in Ref. [39]. Primary focus is on the bending responses due to the von Kármán strains. The ESL formulation based on the FSDT is utilized in order to show the connection between in-and out-ofplane deformations. In post-buckling regime, bending occurs, which activates secondary normal stresses of the microstructure. As the unit cell to plate length ratio, l micro /l macro , is not infinitesimal, the microfluctuations of stress at unit cell level contribute significantly to the total normal stress response of the face and web plates.
In order to show these aspects, we first formulate the geometrically non-linear differential equations of the problem at hand. We assume that the microstructure is linearly-elastic. Thus, at micro-scale, analytical, stiffness and strength formulations are possible and the prevailing physics can be explained in parametric form. Then we show how to eliminate the microstructural elements from the plate model. This elimination allows solution of the macro-scale problem with the classical ESL-FSDT formulations that can be found from numerous textbooks and the finite element implementations from several commercial codes. Then, at the post-processing stage we focus on assessment of stress concentrations at faces, webs and welds and reconsider the finite size of the microstructure. This shows the limitations of this paper and also the commercial FSDT-ESL FE-codes on capturing the prevailing post-buckling behavior. It also motivates extensions to the non-classical continuum mechanics formulations.

Notations
The plate is assumed to consist of structural elements with small thickness representing the face and web plate(s). This justifies the use of Kirchhoff hypothesis locally at the face and web-plates. The thicknesses of the top and bottom face plates are denoted by t t and t b, respectively, and these are positioned in the xy-plane. The web plates are in the xz-plane and have thickness t w , spacing s, and height h c . The plate has two coordinate systems, namely: global xyz and local x l y l z l (see Fig. 2). The origin of the global coordinate system is located at the geometrical mid-plane of the plate and the origin of the local coordinate is located at the geometrical mid-planes of face or the web plates under consideration.

Classical, homogenized FSDT for periodic plates
The deformation of the periodic plate is composed of global bending deflection of the mid-surface and local deflection due to warping of face and web plates. Thus, the total deflections along the three coordinate directions of the faceplate can be represented as u x; y; z; x l ; y l ; z l ð Þ¼u 0 x; y ð Þþz/ g;x x; y ð Þ þ z l / x;l x; y; x l ; y l ð Þ ð 1Þ v x; y; z; x l ; y l ; z l ð Þ¼v 0 x; y ð Þþz/ g;y x; y ð Þ þ z l / y;l x; y; x l ; y l ð Þ ð 2Þ w x; y; z; x l ; y l ; z l ð Þ¼w g x; y ð Þþw l x; y; x l ; y l ð Þ : where subscript 0 denotes the in-plane membrane displacements at the geometrical mid-plane, as shown in Fig. 2. The subscripts g and l denote bending actions at the geometrical mid-plane of the entire sandwich plate and face plate respectively. The global deflection of webs is defined completely by the face plate deflection at their intersection, w web = w g . In addition, due to perfect connection with the faces, the webs bend in the y-direction due to the local deflection of the faces, w l . This deflection of the webs is defined completely by the face plate rotations at the web-toface-interface. Thus, additional kinematical variables are not needed for the web deformations. The rotation of the mid-plane of the entire plate is denoted by / g , the in-plane displacements by u 0 and v 0 and deflection by w g . The global and local rotations are which means the sandwich plate behaves according to the FSDT. The relative strains, expressed in column vector, are where the underlined terms denote the von Kármán strains, which are accounted at the macroscale. The plates are assumed to be made of isotropic material and follow Hooke's law. Thus, the elasticity matrices of the different layers (subscripts: t = top face, b = bottom face, w = web) of the sandwich plate are where the web properties are smeared equally over the unit cell width by the rule-of-mixtures. The stresses per layer are: The stress resultants are obtained by throughthickness integration and are Àh=2 r xx r yy r xy Due to the periodic microstructure, smearing the shear modulus of the core by an integrating over the thickness of the sandwich does not result in correct stiffness values. In the x-direction this is due to the shear flow of the thin-walled section. In the y-direction this is due to the higher order warping deformations of the unit cells. Instead, we incorporate the shear stiffness that accounts directly these two effects and write: The equilibrium of a plate element is governed by the following equations: where the underlined terms are associated with local shear-induced warping response. It should be also noticed that due to the von Kármán strains, the in-plane membrane and out-of-plane shear forces are coupled, by Eq. (20). The differential equations in terms of displacements are then: The membrane, membrane-bending coupling and bending stiffnesses are The shear stiffness is obtained by unit cell analysis of the periodic structure. In x-direction (see Ref. [21]) this is: where the shear correction factor is obtained by integration of shear flow through top, web and bottom faces; for details see Ref. [21]. In the y-direction the shear stiffness is computed by unit frame analysis (see Ref. [21] and Fig. 3b) and is given by where k h is the rotation stiffness of the laser-stake weld and k-parameters model the relative stiffness of faces and webs. For symmetric sandwich panels the parameter k Q = 1/2 indicating that equal amount of shear is carried out by top and bottom face plates.

Unit cell analysis for out-of-plane shear
The local response, opposite to stiffener direction, is assumed to have the same values as the global response at the unit cell edges in terms of displacement, w l (x,y) = w g (x,y), w l (x,y ? s) = w g (x,y ? s), as shown in Fig. 3a and b. However, along the unit cell the local displacement field in the y-direction can differ from global due to shear-induced warping (see Fig. 3b). Due to the classical continuum assumptions the average shear strain of horizontal and vertical sliding, c yz = (c yz,h ? c yz,v )/2, is utilized. The local warping deflection due to shear force Q yy is given as (see Ref. [21]): with the k-factors defined in Eq. (30). The curvature is The local bending moments and shear forces are The volume averages are Thus, with these average values the underlined terms in equilibrium Eqs. (18) and (19) and differential Eqs. (24) and (25) based on displacements become zero. The shear force Q yy is produced by the unit cell warping, see Fig. 3b. This gives us the homogenized differential equations that are derived based on the classical continuum mechanics assumption, that l micro / l macro = 0. These equations can be found from numerous text books of composite materials(see, e.g., Ref. [40]). In these books the corresponding analytical solutions and finite element formulations are also presented. The corresponding finite element implementations can be found from most of the commercial FE codes.

Calculation of envelope of periodic stresses from homogenized response
The periodic shear-induced bending moments and normal stresses (see Fig. 3c) are important as the with these included the total stress always exceeds the averaged, homogenized solution values. The characteristic length, l micro , of these normal stresses is equal to web plate spacing, s, l micro = s. These shear-induced bending moments and normal stresses have zero mean as shown in Fig. 3d and Eq. (36). Next we assume that, as the homogenized structure is in equilibrium, the periodic structure must be in equilibrium too. Then, the periodic strains and stresses are (see Ref. [21]): where Eqs. (34) and (38) together create periodic strain field from homogenized, smooth, strain field. From strength viewpoint shape of this micro-fluctuation is not important, but the maximum value is. Thus, for linear distribution within unit cell, (see Eq. (34) and Fig. 3c and d) It should be noted that when l micro = s=0, these bending moments and resulting shear induced stresses become zero and result is the same as produced by the homogenized solution without any additional postprocessing.

General
The example presented here is taken from Ref. [39]; it is extended here in order to explain the main postbuckling phenomena. A square plate with length and width of L = B=3.60 m is considered. Thickness of the face and web plates are t t = 2.5 mm and t w = 4.0 mm, respectively. Core height is h c = 40 mm and the web plate spacing is s = 120 mm giving l micro /l macro = s/ L = 1/30. The interface between web and face plate is assumed to be rigid in order to simplify the analysis. It has been shown in Refs. [21,36] that non-local plate formulations are needed in the cases were the laserstake weld is assumed to be flexible or if l micro /l macroratio would be significantly larger. In cases where the l micro /l macro -ratio would be significantly smaller, the shear-induced stress fluctuations would be less important. Material is assumed to be linear-elastic with Youngs modulus 206 GPa and Poisson ratio 0.3.
The FSDT-ESL problem is solved numerically by using Finite Element Method. The non-linear analysis is carried out in two steps. The first eigenmode is first computed. It is used as the shape of the initial out-ofplane imperfection and is given the magnitude of 0.01% of the plate length, L. Then, the geometrically non-linear analysis is carried out to trace the postbuckling path. Abaqus software, version 6.9, is used. A subspace iteration solver is used for the eigenvalue analysis and the modified Riks procedure for the postbuckling path. In order to secure converged results in FSDT-ESL, a mesh of 50 9 50 S4R shell elements are used. Simply supported boundary conditions are considered, with the loaded edges kept straight and the unloaded edges free to pull in.

3D Finite Shell Element Analyses for Validation
In order to validate the FSDT-ESL approach, a 3D model of the actual periodic plate is used; see Fig. 4 for details. The 3D-plate is modelled using shell elements (S4) that follow the Kirchhoff hypothesis (assumption in FSDT-ESL model). Concentrated nodal forces act at web plates in the nodes in the geometrical mid-plane. Six shell elements per web plate height and between webs are used as this has been shown to produce converged results in buckling and bending problems [20,41]. Simply supported boundary conditions are considered, with the loaded edges kept straight by constraint equations and the unloaded edges free to move in-plane. The transverse deflection is zero only at the nodes at the geometric mid-plane. This allows the rotation of the plate around the mid-plane edge.

Results
The load-end-shortening, load-deflection and loadout of plane shear behavior is presented in Fig. 5 Fig. 6 for load points A and C. Figure 7 present the corresponding normal stresses with shear-induced stresses included and excluded. From Fig. 5 it is clear that the non-linear response is predicted very accurately with the FSDT-ESL in comparison to the 3D-model of the actual geometry (3D-FEA). Both load-end-shortening and load-deflection curves overlap until local buckling occurs at the unit cells at point C. It is also clear that this point is well beyond the panel level buckling, around point B. After this point B, the out-of-plane shear, Q yy , increases rapidly due to the von Kármán non-linearity, see Eq. (20). Figure 6 shows that the spatial distribution of out-of-plane shear is accurately predicted with the FSDT-ESL. Even at the point of local buckling, the FSDT-ESL and 3D-FEA results overlap each other in Fig. 4 Description of the constrains implementation of 3D-FEA. Reproduced from Ref. [39] average sense. Figure 7 shows that the membrane stresses of the face plates, predicted by both FSDT-ESL and 3D-FE methods, overlap at the load point A. In post-buckling, due to global bending, they start to differ significantly from each other. It is also seen that the membrane stress in the faces is significantly lower than the total stress at the surface of the face plates which is magnified by the secondary bending. It is seen that the stresses from 3D-FEA are within the maximum stress envelope curves defined by Eq. (41). This highlights the importance of taking the finite size, l micro /l macro , of the microstructure into account when computing the stress response. The stress jumps at the location of webs close to plate edges indicate significant bending in the webs and welds due to shearinduced warping of the unit cells.

Discussion
The investigation presented above indicates that when sandwich structures, with visible periodic core, are homogenized, special attention must be paid on the stress assessment, even in the case of buckling assessment where loading is of membrane type. In linear regime this loading is carried out mostly by pure membrane actions of face and web plates as Fig. 7 at point A indicates. As the post-buckling takes place, the membrane actions start to interact with the out-ofplane deformation, due to the von Kármán nonlinearity. As the core is visibly discrete, l micro /l macro-= 1/30, and the unit cells warp in shear, the magnitude of shear-induced secondary normal stresses at the faces and webs becomes significant. When these are added to the membrane stresses, the total stress can be significantly higher than the membrane stress only. This effect becomes larger when the role of out-ofplane shear c yz or Q yz increases and when the unit cell size increases. So, this phenomenon can be considered as type of size effect. Usually in homogenization theories we assume that these are infinitely far apart, i.e. l micro /l macro = 0. When departing from this assumption, very soon we end up dealing with another type of size effect, the assumptions of non-classical continuum mechanics.
This issue is seen through following example. Postbuckling analysis requires two stage analysis with initial analysis to define the initial deformation shape which is then followed by geometrically non-linear analysis. In present formulation, the bifurcation buckling load for simply supported rectangular plate is given as [40]: where m and n are used to denote the number of halfwaves in directions x and y respectively. In present case, the loading is assumed to be uniaxial, thus the load ratio factor k = 0, orthotropy ratio for shear, D Q1 / D Q2 -ratio, is very high (e.g. 100-1000). Thus, the buckling load minimum is obtained when m = n = 1. Problems with this equation occur when the load is turned to be along y-axis or biaxial compression is considered. In this case, the high orthotropy in shear causes a situation where minimum does not converge for m = n = 1, but decreases as function of n, while m = 1. In this case the buckling load in finite element solution becomes mesh size dependent; the smaller is Top: x-normal stresses due to membrane action for load-levels A; Center at load-level C and Bottom: at load-level C with envelopes of the shear induced secondary normal stress the mesh size, the smaller is the buckling load. This is another type of size effect which results from classical continuum mechanics assumptions. The decrease is unphysical and can be corrected by incorporation nonclassical continuum mechanics assumptions into our FSDT-ESL model. As shown by Romanoff and Varsta [21], Jelovica and Romanoff [42] by thick-face plates sandwich theory and by Karttunen et al. [36] by micropolar theory, the shear deformations of sandwich panels can have only finite wave-lengths, that is, finite n-values. When using thick-faces effect or micropolar solution, physically correct behavior is obtained. This is a result of the fact that we can split the out-of-plane shear strain to symmetric and antisymmetric parts [36]. Therefore, the present investigation should be extended in future to compression in transverse direction and investigations based on non-classical continuum formulations. This calls for micropolar plate elements, and recently such study has been reported by Nampally et al. [43]. In order to assess the strength of real welded structures following issue must be handled. The stress values seen in point C of the case study are very close to the material yield point. Typically, the steel used for these panels has a yield point at 355 MPa for the faces and 235 MPa or 355 MPa for the web plates. Thus, as Fig. 7, shows it is crucial to recover the microstructural stress if first fiber yield is to be assessed. The assessment of laser-stake welds is a bit more challenging. Jutila [44] carried out experiments with Digital Image correlation system of pull-out strength of laser-stake welds. Strength values up to 1000 MPa for tension for steel faces and webs of 355 MPa and 235 MPa, respectively are reported. The difference is  [43,44] also seen in the hardness of the welds, see Fig. 8. Due to rapid changes in the hardness, the weld deforms only moderately in softer faces and webs and microrotation is seen between these two, in the heat affected zone. This is indicating that non-classical continuum mechanics are needed in the weld modelling. However, there is another effect that requires careful investigation. This is the contact between the faces and webs at the laser-stake welds when being bended. Due to the contact, the stresses are redistributed and significant variations are seen in the moment carrying capacity of the welds due to small variations in weld geometrical properties. Measured strength values are M cr = 1400Nm/m (see Ref. [45] for details), which corresponds Q yz = 21.5kN/m in shear. Thus, the strength of welds would not be reached yet in present context for individual stress components. Proper modelling of these phenomenon in 3D FEA is therefore a significant challenge which we cannot solve in the present context. Also, an experimental study is needed to gain more understanding of prevailing failure modes in real structures.
In the present case it is indicated that the microstructure can buckle during the deformation, which is known to reduce the stiffness globally. As the local and global deformations interact, there is a need for coupled models where the geometrical nonlinearity in one length-scale can be mapped correctly to the another one. Such example in context of classical continuum mechanics has been presented in Reinaldo Goncalves et al. [46] with one way coupling and in for example Geers et al. [21] and Rabzcuk et al. [30] for two-way coupling. However, this work should be extended to the non-classical continuum mechanics due to the reasons mentioned above.

Conclusions
The paper presented, a phenomena related to postbuckling of web-core sandwich plates. Equivalent Single Layer (ESL) formulation with First order Shear Deformation Theory (FSDT) was used to identify the problem parameters in closed form. During the axial load increase, these plates have multiple load-carrying mechanisms that change due to von Kármán nonlinearity. During the membrane-dominated loading stages, the effect of finite size of the periodic microstructure is barely present, while at higher loads when bending is activated, the finite size of the unit cells activates secondary shear-induced bending moments. Due to this effect, the normal stress levels become significantly higher than the homogenized plate theory would predict. A method to capture the envelope of the maximum values of these stresses is presented and validated with 3D-Finite Elementmodels of the actual geometry. It is also discussed that classical continuum mechanics has its limits when failure of welded web-core sandwich structures is concerned. The welds experience micro-rotation close to failure point. Due to this also the antisymmetric outof-plane shear strain activates at the plate level which calls for non-classical continuum mechanics. These extensions are left for future work. SSAB, Deltamarin, Koneteknologiakeskus Turku and Meyer Turku Shipyard. The paper was finalized in the Academy of Finland project (#310828) called ''Ultra Lightweight and Fracture Resistant Thin-Walled Structures through Optimization of Strain Paths''. This support is gratefully acknowledged. Thanks to IT Centre for Science for providing the computational capabilities.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.