Effect of Membrane Components of Transverse Forces on Magnitudes of Total Transverse Forces in the Nonlinear Stability of Plate Structures

For an isotropic square plate subject to unidirectional compression in the postbuckling state, components of transverse forces in bending, membrane transverse components and total components of transverse forces were determined within the first-order shear deformation theory (FSDT), the simple first-order shear deformation theory (S-FSDT), the classical plate theory (CPT) and the finite element method (FEM). Special attention was drawn to membrane components of transverse forces, which are expressed with the same formulas for the first three theories and do not depend on membrane deformations. These components are nonlinearly dependent on the plate deflection. The magnitudes of components of transverse forces for the four theories under consideration were compared.


Introduction
In the mid-20th century, Reissner [1,2] presented a plate theory accounting for the transverse shear deformation effect. This is a stress-based approach. Mindlin [3] offered a theory based on a displacement approach, where transverse shear stress was assumed to be the same through the plate thickness, and the shear correction factor k 2 (the so-called Mindlin correction factor) appeared. For transversely inextensible plates and k 2 = 5/6, values of stresses are equal in the Reissner and Mindlin plate theories [4]. The theories are characterised by an equivalent approximation degree known as the Reissner-Mindlin plate theory. Their comparison is discussed in [5][6][7][8], etc. Theoretical considerations can be found in [9][10][11][12][13][14][15] for higher-order shear deformation theories as well.
In [4,5], the equations of Reissner and Mindlin plates, including the parameter, which allows for an interpretation of these theories for transversally inextensible plates, were derived. In [6], for Reissner, Mindlin and Reddy plate models, a solution to the rectangular transverse plate sinusoidally loaded and freely supported along all the edges is given. Vibrations in the plate-beam system, in which the Reissner-Mindlin plate model related to the Timoshenko beam model was applied, were analysed in [7]. In [11,16], it is shown that when the plate thickness is around zero, the solution for the Reissner-Mindlin plate becomes close to the solution within the Kirchhoff-Love plate theory (the so-called classical plate theory (CPT)). In [17], formulations of the mixed finite element were given based on the mechanism of the shear locking phenomenon and the general variation method for the Reissner theory, including the Lagrange multiplier method. An extensive literature survey devoted to the two primary plate theories, i.e., the Kirchhoff plate theory and the Reissner-Mindlin plate theory, is found in [18]. The main purpose of [18] was to present a history of refinement in the Reissner-Mindlin theory, showing the up-to-date state of knowledge in this field of research.
Endo and Kimura [19] first suggested the simple first-order shear deformation theory (S-FSDT). Employing Hamilton's principle, an alternative formula in which a deflection in bending is the primary variable instead of the angle of rotation in bending and, at the same time, some limitations on neglecting the Reissner boundary effects are imposed, was given [2,32]. In the S-FSDT, three equations in the original formulation are reduced to two and the boundary conditions are also subject to respective modifications; however, the way the system is modelled remains unaltered (for a more detailed analysis, see Appendix A). Moreover, in [26][27][28], when two independent variables ϕ and w s are considered, two differential equations with boundary conditions are attained. In the static analysis, two differential equations are uncoupled. The boundary conditions should be uncoupled as well, which is not possible in general. By introducing a bending relationship between the quantities, differences between the Reissner and Mindlin plate theories were investigated in [8]. In [24], the first two-variable shear deformation theory (FSDT) considering in-plane rotation, which allows one to correctly predict the response of plates for arbitrary boundary conditions in the analysis of buckling and vibrations of isotropic plates, was presented.
Interesting critical remarks to the above-mentioned plate theories can be found in [32][33][34][35][36]. To consider the Reissner boundary effect, a rotary potential, which is a fast-varying solution to the boundary layer, should be applied apart from the function ϕ. The boundary effect covers only some boundary conditions (e.g., pure tension or contact problems). The author of [32] suggested to refer to the presented version of the theory as a modern form of the CPT.
The finite element method employs the FSDT [9,18,33,[37][38][39][40][41]. The effect of shear locking in finite shell elements and a loss in accuracy was explained in terms of the occurrence of solutions to the boundary layer. A shear locking problem occurs in the FEM, as shape functions cannot approximate a fast-variable solution to the boundary layer. [33,39]. Shear locking does not cause membrane deformations. In the majority of cases, the effect of shear deformation on displacements should be considered only. Solutions to the boundary layer are neglected.
In composites widely applied nowadays, the behaviour of individual layers can be affected considerably by transverse shear deformation [42]. These materials show low shear characteristics beyond the plane, which should be accounted for in numerous analyses (e.g., [15,20,21,43]).
In the works devoted to the plate dynamics, the effect of transverse rigidity on shear deformation and the influence of rotary inertia on frequencies and modes of plate vibrations were investigated (e.g., [3,5,15,17,27,30,31]). In [44], for the plate jointly supported along the whole circumference and subjected to free vibrations, it was shown that deflections in bending and shear deformation vibrated in phase (i.e., the total deflection is equal to their sum) in the first branch, whereas in the second branch, shear deformation and bending deflections vibrate in antiphase, where the deflection in shear deformation is predominant, i.e., the total deflection is in phase along with the deflection in shear deformation.
In the literature, apart from vibrations, the stability of individual rectangular plates is analysed (e.g., [15,16]). In [22], a buckling analysis of isotropic and orthotropic plates employing a two-variable refined plate theory is presented, whereas in [43], the S-FSDT for composite plates with four unknowns is discussed. The authors do not know any works devoted to the nonlinear stability of plates accounting for the transverse shear deformation effect.
All the works under discussion refer to the cases when there are no membrane forces in the plate structure. These forces appear in thin-walled structures (where h/a < 0.05) for loads exceeding the critical loads, that is to say, for postbuckling equilibrium paths. In the literature known to the authors of the present paper, there is a lack of works devoted to the nonlinear stability of thin-walled plates, in which transverse shear deformation is considered.
In the classical theory of thin plates (CPT), total equivalent transverse Kirchhoff forces were introduced only in [36,45] employing the variational method as far as the literature known to the authors is concerned. In the CPT, a notion of equivalent Kirchhoff transverse forces is introduced to satisfy a proper number of boundary conditions. In the variational approach to the CPT, there is no need to introduce the notion of Kirchhoff forces. These forces "emerge themselves" from the theory in such an approach. As shown in [46], it is necessary to introduce the notion of total equivalent Kirchhoff transverse forces, which results from Stokes' theorem concerning a change of the surface integral for the equilibrium equations into a plate circumference-oriented integral, that is to say, for the boundary conditions. In total Kirchhoff forces, two components of transverse forces appear; one of them is a derivative of internal moments, and the second is a projection of membrane forces on the transverse direction. The membrane forces also appear in other nonlinear problems, such as the deflection of thin-walled transversely loaded plates.
In the present work, the authors have decided to deal with the influence of these additional membrane components on the magnitudes of total Kirchhoff forces within the CPT for isotropic square plates subject to compression in the postbuckling state. These limitations were taken to facilitate an interpretation of the obtained results. For verification purposes, solutions to the Reissner theory (FSDT) and the Mindlin theory within the S-FSDT approach, i.e., after the introduction of two independent functions of displacements along the z-axis (i.e., the total lateral displacement w and the bending deflection φ), are presented. In the FSDT and the S-FSDT, the Reissner boundary condition was neglected. The governing equations within the three theories under consideration were derived with variational methods, allowing one to indicate two different components of transverse forces resulting from internal moments and membrane forces. For transversally inextensible plates, the membrane shear forces are independent of membrane deformation. For these three theories, the results for membrane forces and total forces were presented.
In composite materials, transverse shear deformation substantially affects the delamination of composites. In the failure criteria of composites, the impact of transverse components of membrane forces (i.e., in compression) is neglected. In the authors' opinion, these components are predominant in the postbuckling state and should be considered in composite failure criteria. The main aim of the paper is to draw attention to the theoretical background for membrane components of transverse forces in the expressions for transverse forces in the theory of thin plates, which are not accounted for in FEM shell elements.

Formulation of the Problem
The nonlinear stability of a square isotropic plate freely supported along the whole circumference and subject to compression along the x-axis ( Figure 1) is analysed. The plate material is assumed to obey Hooke's law. In this study, for postbuckling equilibrium paths, transverse shear forces are analysed in detail for the transversally inextensible plate. The analysis is conducted within three theories of thin plates, namely the classical plate theory CPT (i.e., the Kirchhoff plate theory), the simple first-order shear deformation theory (S-FSDT) in a two-variable refined plate version and the Reissner plate theory (FSDT).
The governing equations of the three theories under consideration are presented in Appendix A. The equations were derived within a variational approach, which allows the equilibrium equations and the boundary conditions to be expressed explicitly. The solutions to the nonlinear problem of stability of the square plate for the three theories are presented in Appendix A. Instead of a system of two equations of equilibrium in the central plate plane (i.e., after an introduction of the function of Airy forces F, the system is satisfied identically), an equation of inseparability of deformations was derived (Appendix A).
According to the considerations presented in the Appendix A, transverse shear forces have two components (compare: the FSDT (A13), the S-FSDT (A24), the CPT (A33), respectively).
The first components depend on the derivatives of internal moments on the plate. Thus, the components can be referred to as transverse shear forces in bending. These forces are expressed with the following relationships: (A45) for the CPT and (A55) for the S-FSDT, correspondingly. The forces have a very similar structure. A difference lies only in the reduction factor 1/(1 + η) in (A55). Moreover, for the FSDT in (A64), a difference with respect to the S-FSDT occurs in the numerical coefficient two instead of (3-ν) for the S-FSDT. A change in the numerical coefficient results from different boundary conditions for the FSDT and the S-FSDT. A more detailed analysis can be found in the Appendix A.
The second components depend on projections of membrane transverse forces on the direction perpendicular to the central plate plane. These components can be referred to as transverse shear forces in compression. For the three theories under consideration, these forces are expressed with identical formulas (compare (A46) for the CPT, (A56) for the S-FSDT and (A65) for the FSDT in the Appendix A). It is caused by the fact that the effect of shear deformation is not accounted for, as the forces are determined on the basis of the displacement w and the function of Airy forces F.

Analysis of the Calculation Results
A detailed analysis was conducted for a steel square plate ( Figure 1) of the following dimensions: a = 100 mm, h = 1 mm and the material constants: The ideal plate is supported freely along all edges and subjected to uniform compression with the stress p along the x-axis. The boundary conditions for the three theories under consideration (i.e., the FSDT, the S-FSDT and the CPT) are given in detail in the Appendix. The analytical results attained were verified with the commercial ANSYS software [47] employing the FEM (details to be found in Appendix A).
In the detailed analysis, the postbuckling state (or the so-called postbuckling equilibrium path) was dealt with, as only then plate deflections appear for the perfect plate. It is accompanied by the appearance of two transverse components of shear forces, that is to say, in bending and compression (the so-called membrane components).
The following index symbols are introduced in the study: C for the CPT, S for the S-FSDT, F for the FSDT, and A for ANSYS (FEM), respectively.
Firstly, for the three theories, the corresponding bifurcation loads (or the so-called critical loads), listed in detail in Table 1, were determined. According to the Appendix, values of the bifurcation loads for the FSDT (A62) and the S-FSDT (A54) are identical and slightly lower by the factor 1/(1 + η) than the CPT. For the data assumed in the analysis, according to (A63), we have η = 0.000564, which corresponds to 1/(1 + η) = 0.9994. As can be seen, according to (A63), corrections for the S-FSDT and the FSDT are very inconsiderable when compared to the CPT for the assumed ratio of (h/a = 0.01). The results obtained within the three theories are in conformity with the FEM outcomes. The determined value of the critical stress (A43) for the CPT was introduced into the relationship for the total equivalent Kirchhoff forceQ C x (A47a) and then the componentQ C x when p = p C cr , is equal tô DWcos πx a sin As can be easily noticed, in (1) there is a minus sign at both terms of the right-hand side. However, mutual relations depend on the relationships of products of trigonometric functions. The first term Q C x p = p C cr in (1) attains extreme values for x = 0; a and y = a/ 2 and the second respectively minimum for x = 0; a and y = a/ 2. Attention should be drawn to the fact that when x = 0 and y = a/ 2, the first term has a minus sign, and the second term has a plus sign. The opposite situation takes place when x = a and y = a/ 2, i.e., a plus sign is in the first term and a minus sign is in the second. The extreme valueŝ Q C x are attained inside the square plate. A further analysis dealt with postbuckling states. According to (A39c), the force component N xy equals zero, and one of two membrane components of transverse forces (according to (A12), (A23) and (A32)) vanishes as well.
Next, maximal absolute values of components of transverse forces in bending, membrane components (or in compression) or total forces for the three theories, determined according to the formulas given in Appendix A for five overload values of critical load, i.e., 1.2 ≤ p θ /p θ cr ≤ 2.0 (where the index θ = C, S, F), are listed in Table 2. In this table, values of the dimensionless deflection W/h and 1/(1 + η) are also presented.
For the CPT, the equivalent Kirchhoff forces Q C x , Q C y are equal according to (A45). On the other hand, values of membrane components of the transverse forces Q C x , Q C y differ depending on the overload p θ /p θ cr . For the overload equal to 1.2, the ratio of maximal absolute components Q C x /Q C y equals almost 5, whereas for the overload equal to 2, the ratio of membrane components is 1.4. The membrane forces Q C x , Q C y are independent of membrane deformation. It results from the fact that the membrane force Q C x has a term linearly dependent on deflection and in the third power, which for Q C y is in the third power only. A detailed analysis can be found in Appendix A. Components of the total equivalent Kirchhoff forceQ C y are always higher for the range of loads under analysis thanQ C x . It follows from the term that is linear with respect to W, dependent on the overloadQ C x . In Figure 2, the maximal absolute values of components of transverse Kirchhoff forces for the CPT versus p C /p C cr , listed in Table 2, are presented. For the S-FSDT, the maximal absolute values of components of the transverse forces Q S x max and Q S y max are the same in practice as for the CPT, which results from a very low value of the correction η. For the FSDT, force components in bending are 1.35 times lower for the CPT and the S-FSDT (cf. Formulas (A64) and (A45)). Similarly as for the S-FSDT, components of the total transverse forcê Q F y are always larger thanQ F x . When total transverse forces are accounted for in the CPT, the S-FSDT and the FSDT, they yield higher values than the equivalent Kirchhoff force by approx. 1.5 times for the component with respect to the x-axis (i.e., with a lower index x) and more than 2 times for the component with respect to the y-axis.
For the FEM, the values of components of the transverse forces Q A x max are higher than Q A y max . At the overload equal to 1.2, the ratio Q A x max / Q A x max is 1.05, but for the overload of 2.0, it is equal to 1.36, respectively. The values Q A x max and Q A y max are closest to the equivalent Kirchhoff force Q C x , Q C y . Thus, the transverse forces Q A x max , Q A y max determined within the FEM have a different character than the total transverse forces for the CPT, the S-FSDT and the FSDT, determined on the basis of components in bending and compression. It can originate from the fact that membrane components were neglected in the FEM analysis.
In Figure 3, the maximal absolute values of transverse resultant forces for the CPT, the S-FSDT and the FSDT (of which the values are listed in Table 2) versus overload are collected. The results for ANSYS are shown as well. The present study was primarily aimed at drawing attention to a necessity to consider the effect of membrane components on total transverse forces for 1.2 ≤ p θ /p θ cr ≤ 2.0. To illustrate the effect of components of transverse Kirchhoff forces for the CPT, their distributions are presented in subsequent figures (Figures 4-9). The components are shown in the contour drawings (denoted as a) and plane drawings (denoted as b) for the whole square plate and the overload p C /p C cr = 2 in the following sequence: Q C x (Figure 4), Q C x ( Figure 5),Q C x (Figure 6), Q C y (Figure 7), Q C y ( Figure 8) and Q C y (Figure 9). The distributions of Q C x and Q C y presented in Figures 4 and 7 are the same according to (A45). When the distributions of membrane forces are compared, the distribution for Q C x ( Figure 5) is more complex than Q C y (Figure 8). However, in the authors' opinion, the total transverse forcesQ C x ( Figure 6) andQ C y (Figure 9), which should be employed in failure criteria, are the most important. The componentsQ C y are larger thanQ C x , opposite to what happens in the FEM (Table 2). It should be underlined once again that for the square plate under analysis, one of two membrane force components, which is dependent nonlinearly on the deflection W, equals zero (i.e., for N xy = 0).

Conclusions
The effect of membrane components of transverse forces on total transverse forces within the three theories: the CPT, the S-FSDT and the FSDT, was discussed. When membrane components are taken into consideration, an increase can be observed in transverse forces equal to 1.5 times, at least for the square plate, freely supported along the whole circumference under consideration. It results from the fact that membrane components of transverse forces depend nonlinearly on the plate deflection. The results were compared to the FEM. The membrane transverse forces are independent of membrane deformation.
In composite materials, a failure of the structure resulting from delamination exerts a considerable effect on its integrity and load-carrying capacity. Transverse shear effects significantly influence the behaviour of composites. In the composite failure criteria, the impact of transverse force components in compression is neglected. From the authors' viewpoint, these components prevail in the postbuckling state, which was proven in this study and should be considered in the failure criteria of composites, for instance, the Hashin failure criterion for 3D, LaRC04(3D), Matrix Failure under the additional condition that σ 33 = 0. The following geometrical relationships for the plate component are assumed [36,44,46] And where u, v are components of the plate displacement vector along the axis x, y, correspondingly, w is the total displacement vector along the z-axis and ψ x , ψ y are the rotation angles of a transverse normal due to bending of the axes x, y, respectively, whereas the plane x − y overlaps the central plane before its buckling.
In this approach to the plate theory when transverse shear is accounted for, it is assumed that the total rotation angles of the normal to the central plane in two planes are, respectively [19,44] w ,x = ψ x + β x w ,y = ψ y + β y (A3) where β x , β y are the transverse shear angles. Internal cross-sectional forces can be expressed in the form [19,44] The coefficient k 2 occurring in Formulas (A6) is known as the Mindlin correction factor. The total potential energy Π of a thin rectangular plate of the following dimensions: × b × h can be written as [36] where U is the elastic strain internal energy, W is the work of external forces and p 0 (y) is the plate prebuckling external load in the central plane along the x-axis. [Q x,x + Q y,y + (N x w ,x + N xy w ,y ) ,x + N xy w ,x + N y w ,y ) ,y δwdxdy = 0 (A8c) N xy w ,x + N y w ,y + Q y δwdx y = const = 0 (A10c) 0 M xy δψ x dx y = const = 0 (A10d) 0 M y δψ y dx y = const = 0 (A10e) Equation (A8) is a system of equilibrium equations and relationships, and (A9) and (A10) are boundary conditions for x = const and y = const, respectively. Equations (A8)-(A10) hold for the FSDT. The above equations were determined from variational methods for which, according to Stokes' theorem, the surface integral can be transformed into a circumference-oriented integral. The form the equations are presented follows from it.
According to (A6) and (A8d,e), transverse forces are expressed with the relations where the upper index F was introduced for the FSDT. These are components of transverse forces dependent on derivatives of internal moments.
In Conditions (A9c) and (A10c) for the components of transverse forces, forces of the same characters are added, i.e., additional components of transverse forces that depend on projections of membrane forces, that is to say [36,45] These forces do not affect membrane deformations. When (A11) and (A12) are considered, according to (A9c) and (A10c), the total transverse forceŝ Next, the independent functions ψ x , ψ y (A2) were substituted with a potential function φ(x, y) [19,[32][33][34]44] This approach consists of the fact that the function φ is treated as deflections in bending. Because the function φ is introduced, the Reissner boundary conditions [32], not considered in the present study, are ignored.
If β x = w s,x , β y = w s,y hold for (A3), the shear deflection is obtained from the relationship w s = w − φ. From (A3) and (A14), it follows that When (A15) is accounted for, the internal forces (A5) and (A6) are written as From the variation of the total energy Π (A7) with respect to the displacement variable components u, v, w, φ, we obtain [Q x,x + Q y,y + (N x w ,x + N xy w ,y ) ,x + N xy w ,x + N y w ,y ) ,y δwdxdy = 0 (A18c) for the plate corner, i.e., for x = const and y = const 2M xy δφ|x = const y = const = 0 (A21) Equations (A18)-(A21) correspond to the S-FSDT, i.e., to a two-variable refined plate theory. According to (A6), (A19d) and (A20d), the transverse forces dependent on the variable φ are expressed as where the upper index S was introduced for the S-FSDT. The membrane components of the transverse forces dependent on the variable w, according to (A19c) and (A20c), take the form Analogously to (A13), the total transverse forcesQ S x andQ S y for the S-FSDT are as followŝ Comparing the formulas for transverse forces in pairs of (A11) and (A22), as well as (A13) and (A24), it can be easily noticed that in the case of the S-FSDT, we have a coefficient two at a derivative of the torque M xy , which for the FSDT is equal to one.
In the classical theory of plates (CPT), transverse forces are neglected (A6) and, moreover, in (A2) and (A3), it should be [36] w ,x = ψ x w ,y = ψ y (A25) Taking into consideration (A25) in (A5), we have [46] When the above-mentioned relations are accounted for and it is assumed that Q x = Q y = 0 in (A7), the following equations are obtained [46] • Equations of equilibrium for y = const 0 N y δvdx y = const = 0 (A29a) 0 N xy δudx y = const = 0 (A29b) 0 M y,y + 2M xy,x + N xy w ,x + N y w ,y δwdx y = const = 0 (A29c) 0 M y δw ,y dx y = const = 0 (A29d) for the plate corner, i.e., for x = const and y = const In the history of the CPT, for the first two components in (A28c) and (A29c), a term of equivalent Kirchhoff transverse forces, was introduced and defined as where the upper index C was used for the CPT. By analogy to the FSDT and the S-FSDT, according to (A28c) and (A29c), the following components were assumed Q C x = N x w ,x + N xy w ,y (A32a) The above-mentioned components of transverse forces result from the projection of membrane forces. Thus, they can be referred to as equivalent Kirchhoff membrane forces.
Taking into account the two above-mentioned systems, the total equivalent Kirchhoff transverse forcesQ C x andQ C y for the CPT were written as [36,45,46] Q C y = Q C y + Q C y = M y,y + 2M xy,x + N y w ,y + N xy w ,x (A33b) Formulas (A24) and (A33) have the same structure. It should be remembered that for the S-FDST, these equations are for the two variables w, φ, whereas for the CPT, only for the variable w.

Appendix A.1.4. Shear Forces
When we compare relationships (A13) for the FSDT, (A24) for the S-FDST and (A32) for the CPT, we can see that the expressions for total transverse forces are identical in practice. In the FSDT, at the derivative of the torque M xy there is a coefficient equal to one and not two, as it takes place for the S-FSDT and the CPT. It is caused by additional boundary conditions (A9c) and (A10c) for the FSDT.
For postbuckling states, there are membrane forces in the central plane that yield simultaneously projections for additional components of transverse forces (A12), (A23) and (A32) for the FSDT, the S-FSDT and the CPT, correspondingly. In the case of postbuckling equilibrium paths, these additional components are larger than the components of transverse forces (A11), (A22) and (A31).
Appendix A.2. Solution to Governing Equations within the FSDT, the S-FSDT and the CPT A square isotropic plate freely supported along all edges, compressed along the x-axis (Figure 1), was analysed. The plate of the dimensions a and the thickness h was assumed to have the following material constants: Young's modulus E and Poisson's ratio ν. The considerations were limited to an elastic range.

Appendix A.2.1. Continuity Equation of Deformations (the So-Called Inseparability Equation of Deformations)
On the assumption that relations (A1) and (A4) hold for the three theories under consideration (i.e., the FSDT, the S-FSDT and the CPT), two first equations of equilibrium for each theory are the same (cf. (A8a,b), (A18a,b) and (A27a,b)). To solve them, a function of Airy forces F was introduced. It is defined as follows [36,45] N x = σ x h = F ,yy (A34a) Taking into consideration (A34) in (A8a,b), (A18a,b) and (A26a,b), it was found that both the equations were identically zero. That system of equations was substituted by one equation referred to as a continuity equation of deformations or an inseparability equation of deformations.
For this purpose, a system of Equations (A1) was rewritten as ε x,yy + ε y,xx − γ xy,xy = w 2 ,xy − w ,xx w ,xy (A35) When we take account of relations (A1), (A4), (A34) in the above-mentioned system, we obtain a continuity equation of deformations in the form [36,45] ∇∇F ≡ F ,xxxx + 2F ,xxyy + F ,yyyy = E w 2 ,xy − w ,xx w ,xy (A36) The equation is linear with respect to F and nonlinear with respect to w. The deflection of a square plate freely supported is approximated in the following way [36] w = Wsin πx a sin πy a (A37) which satisfies the following boundary conditions When (A37) is substituted into (A36), a function of the membrane forces F is determined, and then components of the membrane forces are defined as [36] where p is the stress along the x-axis (Figure 1). The functions of forces (A39) fulfil the following boundary conditions [36,46] u( One of the solutions is Relationship (A39). The second solution is the equation of Equilibrium (A27c) with respect to w, which was solved with the Galerkin-Bubnov method [36].
The assumed function of deflection in (A37), according to (A26), fulfils the boundary conditions of free support To solve the problem, Relationships (A37) and (A39) are introduced into Equation (A27c). Finally, we obtain an equation that describes the postbuckling equilibrium path for the CPT Relationship (A42) can also be expressed as [36] 1 Transverse forces in the CPT are referred to as equivalent Kirchhoff transverse forces (A31) and, according to (A37), expressed with the following relations Wcos πx a sin πy a (A45a) Wsin πx a cos πy a (A45b) From (A32), for equivalent Kirchhoff forces in compression, taking account of (A3) and (A39), we have For transversally inextensible membrane plates, the shear forces are independent of membrane deformation.
The relations for total equivalent Kirchhoff transverse forces result from the two above-mentioned systems of equations and (A33) As one can easily notice, the first terms on right-hand sides in (A47) are positive, whereas the remaining terms are negative. The two first terms in (A47a) and the first term in (A47b) are a linear function of W, whereas the remaining terms are nonlinear with respect to W. To solve the problem, the deflection W should be determined from Equation (A44) for the given load p, and next the values of transverse forces (A45)-(A47).
The first componentQ C x in (A47) depends formally on the value of the compressive load p, whereas the secondQ C y does not. However, it should be remembered that for the given load p > p C cr , we have the deflection W, and thus, indirectly, the equivalent Kirchhoff transverse forcesQ C x andQ C y depend on p and W.

Appendix A.2.3. Solution to the Nonlinear Problem of Stability for the S-FSDT
A solution to the nonlinear stability problem for the S-FSDT consists of three equations: (A36) and (A18c,d). The solutions to (A36) are relationships (A39). The solution to the remaining two equations, when we take into account (A37), is predicted in the form which, according to (A16) fulfils, the boundary conditions of free support Having substituted (A37) and (A48) into (A18d), a linear relation between amplitudes of the functions W and Φ was determined as Nonlinear Equation (A18c) was solved within the Galerkin-Bubnov method, as for the CPT. For this purpose, Relations (A16), (A17), (A39) and (A50) were introduced into the above-mentioned equation, and after some transformations, an equation describing the postbuckling equilibrium path was obtained Relation (A52) can be expressed in an analogous way to (A44) as From the formal point of view, Relationships (A42) and (A52) have an identical structure due load and the deflection W, and a difference lies in critical forces only. From a comparison of (A43) and (A53), we have The transverse components of forces in the S-FSDT, according to (A16), (A22) and (A50), are in the form Comparing relations (A47) for the CPT and (A57) for the S-FSDT, it can be easily seen that the formulas differ only by the factor 1/(1 + η) in the first term of right-hand sides of expressions in (A57). Moreover, the first two components of the right-hand side of (A57a) and the first term in (A57b) depend linearly on W, whereas the remaining ones are in the third power for W.
Appendix A.2.4. Solution to the Nonlinear Problem of Stability for the FSDT For the FSDT, a solution to the nonlinear stability problem consists of a function of Airy forces F (A39) and a system of equilibrium equations (A8c-e). The solutions to these equations are the function w (A37) and the functions ψ x , ψ y , which were assumed in the form The function w (A37) fulfils the boundary conditions (A38), whereas the functions ψ x , ψ y satisfy the conditions of free support ψ x (y = 0) = ψ x (y = a) = 0 (A59a) After the substitution of (A58) into two equations (A8d,e), the following linear dependencies between the functions W, Ψ x , Ψ y are attained where relation (A51) holds. A further step was a solution to the nonlinear Equation (A51), which describes the postbuckling equilibrium path. This solution was obtained employing, similarly as for the CPT and the S-FSDT, the Galerkin-Bubnov method. The postbuckling equation of equilibrium has the form As can be easily seen in (A62), the values of critical loads for the S-FSDT and the FSDT are identical. The postbuckling equilibrium paths for the S-FSDT (A52), (A52a) and the FSDT (A61), (A61a) are the same, and, additionally, they differ from the CPT only by the factor 1/(1 + η).
In a comparison of the Kirchhoff theory (i.e., CPT) with the Mindlin (the S-FSDT version) and Reissner (FSDT) theories, the reduction factor η plays a very important role (A51).
For the CPT, the factor η, according to [3,46], can be written for the case when G → ∞ , which leads to the relation η = 0.
In the last stage, the determined components of transverse forces in the FSDT, according to (A60), (A37) and (A8), were written as When Relations (A47) for the CPT and (A57) for the S-FSDT are compared, one can easily see that they differ in the factor 1/(1 + η) only in the first term of right-hand sides of expressions (A57).
Like for the CPT and the S-FSDT, the first two terms of the right-hand side (A66a) and the first term in (A66b) depend linearly on W, whereas the remaining terms are in the third power for W.
Appendix A.2.5. Solution to the Nonlinear Problem of Stability in the FEM (ANSYS) In the given system of coordinates ( Figure A1), a square plate jointly supported was investigated. The support was achieved taking the degrees of freedom away along the following directions: x for x = 0 (displacement u = 0), y for y = 0 (displacement v = 0), z for z = 0 (displacement w = 0, perpendicular to the central plane of the plate). Numerical calculations were performed with the commercial ANSYS ® software [47] based on the FEM. In Figure A1, a compression force along the x-direction, applied to the node, is indicated. The displacements u and v for x = a and y = a (a-plate length/width), respectively, were assumed to be constant, using the Couple Dof's function, which allows for controlling a group of nodes (the so-called Slaves) with one main node (the so-called Master). In the numerical model, a Shell181 finite element with six degrees of freedom was used. The FEM model had 2500 elements, 15,606 degrees of freedom, and the element size was 2 mm. In the nonlinear problem, an option of Update geom, i.e., a possibility to impose an initial shape imperfection, which is a solution to the plate postbuckling state, was involved. In the computations, the imperfection amplitude equal to w0 = 0.0001 h was assumed. The problem was solved within the Newton-Raphson method for the system of equilibrium equations, which is a system of algebraic equations, for 1000 load steps within the range 1.001 ≤ P/P cr ≤ 2.000.
In the Shell181 element description [46], the transverse forces are denoted as Q 13 and Q 23 . For uniform notations in the presented study, it is assumed that where the index A refers to ANSYS ® [47] based on the FEM. In FEM procedures, the manufacturer issues certain supplements and corrections with respect to the traditional version of the assumed theory, for instance, for the FSDT.