On a sequence of approximate solutions: bending of a simply supported square plate

The solution from Kirchhoff’s theory for a simple textbook problem of bending of a simply supported square plate under doubly sinusoidal vertical load is re-examined concerning the exact solution of the three-dimensional problem. Well-known deficiencies in widely used classical theories such as Reissner’s theory and first-order shear deformation theory (FSDT) along with Kirchhoff’s theory are highlighted. These theories are based on the stationary property of relevant total potential. They do not provide proper initial solutions for the generation of a sequence of approximate solutions converging to the exact solution. The aim of the present work is to eliminate deficiencies in these classical theories. A sixth-order system consisting of three second-order equations is initially formulated with the aid of a normal strain ignored in FSDT. These second-order equations, uncoupled from transverse deflection in Kirchhoff’s theory, govern in-plane displacements and transverse stresses. A supplementary problem consisting of a fourth-order system of equations is included to rectify inadequacies of polynomial expansions in the thickness-wise distribution of approximate solutions. More realistic and practical edge conditions are considered. An auxiliary problem is formulated, and its solution is used in the generation of a proper sequence of two-dimensional problems. Higher-order corrections to primary solutions are proposed through an iterative procedure.


Introduction
One finds several theories reported in the literature for the analysis of bending of plates within the classical small deformation theory of elasticity. These theories are generally based on making suitable assumptions about thicknesswise distribution of displacements and/or stresses (or strains) to derive two-dimensional plate equations from three-dimensional equations. They are primarily intended to overcome inherent limitations in the theory of Kirchhoff (1850) which is simple and continuously used to obtain design information. For a brief survey of most of these theories, one may refer to the following reviews: Reissner (1985), Lo et al. (1977Lo et al. ( , 1978, Lewinski (1986Lewinski ( , 1987, Blocki (1992), Kienzler (2002), and Batista (2010).
Kirchhoff's theory consists of a single-variable model in which in-plane displacements are expressed in terms of gradients of vertical deflection w 0 (x, y) from the assumption of zero transverse shear strains. Solution for w 0 is such that zero face shear conditions and edge support conditions are satisfied. It is governed by a fourth-order equation associated with two edge conditions instead of three edge conditions required in a 3-D problem. A consequence of this lacuna is the well-known Poisson-Kirchhoff boundary conditions paradox unresolved satisfactorily over the past 16 decades. Reissner (1985) in his article felt that inclusion of transverse shear deformation effects was directly (or by implication) presumed to be the key in resolving Poisson-Kirchhoff's boundary conditions paradox. In his pioneering work, Reissner (1944) brought out the use of transverse shear stresses in deriving a sixth-order plate theory. Equations governing (average) displacements in his sixth-order theory (Reissner 1945) consist of Kirchhoff's fourth-order equation and a second-order equation coupled through edge conditions. They correspond to plate element equations governing stress resultants consistent in the reduction of 3-D equations through calculus of variations. After Reissner (1944Reissner ( , 1945, several sixth-and higher-order shear deformation theories were reported. Governing equations derived either from using calculus of variations or from statics correspond to errors in satisfying in-plane equilibrium equations and/or constitutive relations. The assumption of zero transverse shear strains is discarded, resulting in a three-variable model of the widely used first-order shear deformation theory (FSDT), based on Hencky (1947). In this theory, w 0 and in-plane displacements [u, v] = z[u 1 (x, y), v 1 (x, y)] are coupled in governing differential equations and boundary conditions. In Reissner's sixth-order theory (Reissner 1945) and FSDT, applied in-plane shear is combined with transverse shear resulting in a 'torsion-type' problem in flexure instead of in-plane shear combining with reactive vertical shear implied in Kelvin and Tait's (Love 1934) physical interpretation of contracted boundary condition in Kirchhoff's theory. In fact, torsion problem is associated with flexure problem (Vijayakumar 2011b), whereas flexure problem (unlike that directly or indirectly implied in the energy methods) is independent of the torsion problem. This torsion problem is, however, an approximation to the associated torsion problem in the bending of plates and results in unrealistic vertical deflection of the plate. In fact, a solution of the associated torsion problem is to nullify the effect due to applied or reactive edge stress τ xy in the flexure problem (In this connection, a pertinent observation is that Kirchhoff's theory and FSDT are valid for the analysis of hard and soft, simply supported plates, respectively). This can be inferred from numerical results reported by Lewinski (1990) in which higher-order theories give decreasing values of vertical deflection. Recently, Kirchhoff's theory is modified such that the resulting sixth-order theory of bending and St. Venant's theory of torsion are mutually exclusive to each other (Vijayakumar 2009).
In the present author's recently developed theory designated as 'Poisson's theory of plates in bending' , reactive transverse shear stresses are initially determined independent of material constants retaining their thicknesswise parabolic distributions. Dependence on material constants is through higher-order corrections to these stresses. These higher-order corrections are determined through an iterative procedure based on the equilibrium of a 3-D infinitesimal element instead of on the use of the stationary property of total potential in the energy methods. Dependence of analysis on vertical deflection w 0 (x, y) is eliminated (one should note that the coupling of analysis with w 0 in energy methods is due to the work done by the applied transverse stresses during deformation and also the root cause for Poisson-Kirchhoff's boundary conditions paradox). Also, the coupling between flexure problem and associated torsion problem is eliminated through an iterative method (Vijayakumar 2011a). The problem at each stage of iteration is defined by a sixth-order system of equations using a solution at the preceding stage of iteration. In the earlier investigations (Vijayakumar 2009(Vijayakumar , 2011a, in-plane displacements are expressed in terms of gradients of two functions, ψ and φ, in which ψ is related to w and φ is required to decouple bending and torsion problems. The function φ is an auxiliary plane harmonic function from zero rotation (v, x − u, y ) denoted by ω z about z-axis. The function φ was introduced earlier by Reissner (1947) as a 'stress function' governed by a 'wave equation' with imaginary velocity of propagation. Kirchhoff 's theory, FSDT, Reissner's theory, and other shear deformation theories and reported higher-order theories based on energy principles are definitely not useful (explained later) in the generation of a proper sequence of 2-D problems converging to a 3-D problem.
In Kirchhoff's theory and shear deformation theories including the widely used FSDT, vertical deflection w 0 (average displacement in Reissner's sixth-order theory) is same in the face parallel planes. However, it is known physically that neutral plane deflection has to be higher than face deflection. This is due to the fact that an elastic medium is on either side of the neutral plane, whereas it is on one side of each of the top and bottom faces. Proper estimation of neutral plane deflection is dependent on normal strain ε z ignored in Kirchhoff's theory and shear deformation theories. In an earlier investigation (Vijayakumar 2011a), it is shown that the expansion of displacements in polynomials of thickness coordinate z is not adequate to properly estimate face and neutral plane deflections. This fact is overlooked in the analysis of even isotropic homogeneous plates through the widely used FSDT and other shear deformation theories. A solution to a supplementary problem is required to obtain neutral plane deflection. It is, however, observed that the error in the estimation of face deflection is much higher than that of neutral plane deflection. However, it is desirable to provide uniform approximation to deformations through thickness of the plate. It appears to the author's knowledge that no suitable 2-D modeling which gives more or less the same percentage of approximation to thickness-wise distribution of a displacement variable is reported until now. Proper higher-order theories are to reduce only maximum error in the estimation of a physical variable to a desirable level.
Variational methods based on the principle of virtual work consist of virtual increments in displacement variables. Due to arbitrary virtual increments, one derives static equations governing real displacements independent of these virtual displacements which do not participate in these equations. If a 3-D variable, for example, u(x, y, z) is expressed in infinite series with sum on (i, j, k = 0, 1, 2, ••••••) u x; y; z ð Þ¼ ∑α i x i À Á ∑b j y j À Á ∑c k z k À Á then determination of real constants (a i , b j , c k ) from relevant algebraic equations generated in the numerical methods such as finite difference methods, finite element methods, etc. is presumably valid for solutions of such 3-D variables. However, expansion of displacements based on polynomials of thickness coordinate z is not adequate as mentioned earlier in obtaining proper solutions of displacements. Moreover, sequence of 2-D problems derived from variational principles do not converge to 3-D problems since Hence, it is necessary to consider equilibrium equations of a 3-D infinitesimal element (instead of plate element) to generate proper sequence of 2-D problems governing 2-D displacement variables.
The condition zero ω z decoupling flexure from the torsion problem is satisfied in Kirchhoff 's theory. By imposing this condition in FSDT, the sum of in-plane normal strains (ε x + ε y ) in isotropic plates is governed by a second-order equation due to applied transverse loads. Using this sum of in-plane strains together with condition ω z = 0 gives uncoupled second-order equations governing in-plane displacements. Reactive transverse shear stresses are gradients of (ε x + ε y ) uncoupled from w 0 .
In Poisson's theory mentioned earlier, in-plane distribution of (ε x + ε y ) is denoted by ψ(x, y) so that the gradients ψ ,x and ψ ,y give in-plane distributions of transverse shear stresses independent of material constants. The function ψ is related to normal strain ε z , and ψ = 0 implies ε z = 0 along the supported edge. The condition ε z = 0 in place of w 0 = 0 is mathematically valid and more so at point supports along the wall of the plate.
In the bending problem of a simply supported square plate under doubly sinusoidal vertical load, accuracy achieved in the estimation of face deflection (Vijayakumar 2011a) is much less than that of neutral plane deflection at the end of the first stage of iteration. This is likely to persist even in higher-order approximations. This deficiency is rectified earlier (briefly described and used later in the present analysis) through solution of an auxiliary bending problem (Vijayakumar 2013). This auxiliary problem is based on consideration of the influence of σ z = ± zq/2 satisfying face load conditions and associated constant transverse shear stresses on the preliminary solutions of primary bending problems. These transverse stresses are treated as virtual quantities since they do not participate in the static in-plane equilibrium equations as the abovementioned virtual incremental displacements in the application of principle of virtual work. Constant transverse shear stresses, however, contribute in the integrated equilibrium equations, and linear σ z contributes in the in-plane normal stress-strain relations. These transverse stresses are absent in Kirchhoff's theory.

Bending of a simply supported square plate
A square plate bounded within 0 ≤ X, Y ≤ a, Z = ±h planes with reference to the Cartesian coordinate system (X, Y, Z) is considered. The material of the plate is homogeneous and isotropic with elastic constants E (Young's modulus), ν (Poisson's ratio), and G (Shear modulus) that are related to each other by E = 2(1 + ν) G. For convenience, coordinates X, Y, Z and displacements (U, V, W) in non-dimensional form x = X/a, y = Y/a, z = Z/h, (u, v, w) = (U, V, W)/h and half-thickness ratio α = (h/a) are used. With the above notation and ⇄ indicating interchange, equilibrium equations in stress components are (with 3-D stress components as functions of coordinates x, y, and z) in which the suffix after ' , ' denotes the partial derivative operator.
The plate is subjected to asymmetric load σ z = ± (q 0 /2) sin (πx/a) sin (πy/a) and zero shear stresses along z = ±1 faces. In a simply supported plate, conditions along x (and y) constant edges are In displacement-based models, stress components are expressed in displacements via six stress-strain constitutive relations and six strain-displacement relations. These relations within the classical small deformation theory of elasticity are as follows: Constitutive relations: Strain-displacement relations: Methods Analysis consists of exact solutions of 3-D equations of bending problem with and without torsion, 2-D equations from Poisson's theory and extended Poisson's theory.

Exact solutions
With reference to the proper sequence of 2-D problems converging to a 3-D problem, an exact solution of the above-mentioned problem of a simply supported square plate presented briefly earlier (Vijayakumar 2011a) is described here to assess solutions of the present 2-D problems. Equilibrium equations in displacements in the 3-D problem are From satisfying the above equilibrium equations, one Zero shear stresses and vertical load condition along faces give From the above solutions, one obtains with ν = 0.3 and α = 1/6 (i.e., plate thickness ratio 2h/a = 1/3) The above estimates correspond to those from the associate torsion problem since ω z ≠ 0 in Equations 7 and 8. In the case ω z = 0 decoupling bending and torsion problems, equilibrium equations are in terms of (u, v, σ z ) given by Estimates from Kirchhoff's theory and FSDT The above face deflection parameter from Kirchhoff's theory is 2.27. It is the same in all face parallel planes. Correction due to coupling with torsion is 1.45 (Vijayakumar 2011b), whereas it is 1.42 from FSDT and other sixthorder theories (Lewinski 1990). With reference to numerical values reported in Lewinski's article, the above correction is less than 1.23 in Reissner's 12th-order and other higher-order theories. It clearly shows that shear deformation theories and other higher-order theories do not lead to the solutions of bending problems. It is to be noted that vertical deflection w 0 from FSDT is usually considered to be w 0N , but its estimated value (3.69) of the deflection parameter is much nearer the w 0F .
In-plane variables (u 1 , v 1 ) uncoupled from w 0 are basic variables as in earlier investigations (Vijayakumar 2009(Vijayakumar , 2011a. Due to the condition ω z = 0 required to decouple bending and torsion problems, one obtains reactive transverse stresses from thickness-wise integration of equilibrium equations in which f 2 = (1 − z 2 ) / 2 and f 3 = (z − z 3 /3) / 2. Constitutive relation gives ε z = −μ f 1 e 1 with f 1 = z and e 1 = ε z1 . From satisfying face load condition, one gets the equation governing e 1 as In Poisson's theory, e 1 in the above equation is replaced by ψ 2 (x, y), and it is solved with edge condition ψ 2 = 0 in the simply supported plate problem.
As mentioned earlier, the function ψ 2 is related to normal strain ε z and ψ 2 = 0 implies ε z = 0.
Equations governing u 1 and v 1 are given by so that Homogeneous equations show that u 1 and v 1 are conjugate harmonic functions coupled through conditions (2) on u 1 and v 1 along x (and y) constant edge In the above equations, T x and T xy are prescribed stress distributions.
From zero face shear conditions, one obtains w 0 (x, y) in terms of known [u 1 , v 1 ] in the form The vertical deflection thus obtained corresponds to face deflection (note that zero face shear conditions do not participate in the 3-D domain equations unlike in Kirchhoff's theory). Its zero value along the edge requires simply a support to prevent vertical deflection of intersection of the face with the wall of the plate. Such a support also ensures zero deflection of the neutral plane along its edge. It is because the face and neutral plane deflections (w 0F , w 0N ) are one and the same since w(x, y, z) with (w 2 + ε z1 )= 0 can be expressed as However, f 2 (z) of the second-order correction w 2 to the vertical deflection w is parabolic; thereby, the applied or reactive transverse shear stress along an edge is parabolic.
Equations 15 and 17 form a sixth-order system to determine ψ 2 , u 1 , and v 1 . Relevant three edge conditions correspond to those required by Poisson. Hence, the theory is designated as 'Poisson's theory of plates in bending' , parallel to Kirchhoff's theory. It is significant to note that the condition w 0 = 0 prescribed along segments of the edge has no effect on the displacements [u 1 , v 1 , w 0 ] obtained from the present Poisson's theory. It is to be noted that Poisson's theory gives an additional value of 0.27 to w(x, y, 1) due to ε z1 to the earlier w 0F estimate from Kirchhoff's theory.

Numerical results and discussion
After the first stage of iteration involving f 3 (z)[u 3, v 3 ] in the earlier work (Vijayakumar 2011a), higher-order correction to w 0 uncoupled from torsion is 1.26 so that total correction to the value from Kirchhoff's theory is about 1.54. Neutral plane deflection w 0N is corrected from the solution of a supplementary problem. Here the total correction over face deflection is about 0.66 giving a value of 4.46 which is very close to the exact value 4.49. Hence, it is safe to conclude that second-order corrections in the displacements and transverse stresses from the iterative method serve the purpose of assessing data from Kirchhoff 's theory and FSDT.
Solution for w 0 consists of satisfying (i) zero face shear conditions and (ii) edge support conditions in w 0 . In Kirchhoff's theory, both of these requirements are met through the single variable w 0 governed by a fourthorder equation. In Reissner (1945), the average vertical displacement satisfies condition (ii), and in FSDT, condition (i) is not satisfied but its effect is included through shear energy correction factor. In any case, both of them are approximations to associated torsion problem (Vijayakumar 2011b). The present analysis provides some clarity with regard to the preliminary solution of primary flexure problem. It shows that determining in-plane displacements, bending stresses, and reactive transverse stresses from integration of equilibrium equations is uncoupled from w 0 . However, edge condition w 0 = 0 or contracted stress resultant V x = V xo in Kirchhoff's theory is different from the condition along x (and y) = constant edge Error in the estimated value (=3.81) of w 0F is relatively high compared to the accuracy achieved in neutral plane deflection. In the earlier work (Vijayakumar 2011a), estimation of w 0F is further improved by modifying in-plane displacements (u* 3, v* 3 ) such that (τ xz2 , τ yz2 ) are independent of ε z1 . In the present example, correction to the face deflection changes to 1.43 so that face deflection value is 3.97 (=2.27 + 0.27 + 1.43) which is under 4.7% from the exact value. The correction 1.43 to the face deflection is due to ε z3 from the constitutive relation. Determination of ε z3 in terms of (σ z3 , u 3 , v 3 ) involves lengthy algebra and arithmetical work.
The above analysis consists of a basic sixth-order system and a supplementary fourth-order system. Since the second-order corrections are mainly due to inclusion of σ z1 in the in-plane constitutive relations, it is much simpler to find its effect from the following extended Poisson's theory without consideration of higher-order displacement components (u 3 , v 3 ).

Improved estimates
The analysis presented here eliminates the assumptions in Kirchhoff 's theory. Preliminary solution consists of a combination of solutions of four problems presented elsewhere (Vijayakumar 2013) in the analysis of bending of symmetric laminates with anisotropic plies in which the analysis of the present title problem is included as illustrative example. The first one is a solution of an auxiliary problem to account for thickness-wise linear distribution of σ z = z q(x, y) / 2 neglected in the constitutive relations in Kirchhoff's theory and shear deformation theories. The second one is concerned with determining [u 1c , v 1c ] due to σ z = z q(x, y) / 2 in the inplane constitutive relations. In the case of a simply supported plate, [u 1c , v 1c ] are governed by a fourth-order system with homogeneous in-plane edge conditions. The third one consists of solution of a fourth-order system of finding [u 1 , v 1 ] b due to transverse stresses obtained in the second problem. Finally, the fourth one consists of solution of a fourth-order system of a supplementary problem to obtain a more or less uniform approximation of thickness-wise distributions of deformations in the plate.

Auxiliary problem
One finds that σ z1 is neglected in shear deformation theories. Moreover, transverse shear stress conditions along the edges are specified with reference to reactive (statically equivalent) stresses in the form of parabolic distributions, though mathematically valid but not practical. Kirchhoff's theory is based on assuming identically zero transverse shear strains and neglecting σ z1 in the constitutive relations. The former assumption results in highly constrained deformations of the plate.
In view of the above observations, transverse shears along faces of the plate are initially assumed in the form [τ xz0 , τ yz0 ] = [T xz (x, y), T yz (x, y)]. It is more reasonable and practical to assume that they are the same in face parallel planes instead of w 0 (x, y). Correspondingly, [τ xz0 , τ yz0 ] along constant x (and y) edges are They are more realistic edge conditions and more practical if they are independent of in-plane coordinates. In view of Equation 1b, [τ xz0 , τ yz0 ] are expressed as so that ψ 0 (x, y) is governed by in which Δ is the plane Laplace operator (∂ 2 /∂x 2 + ∂ 2 /∂y 2 ) and σ z1 = q/2 with σ z = z σ z1 . The above equation is to be solved with edge condition either ψ 0 = 0 or its normal gradient equal to the applied transverse shear stress along the edge. It is to be noted that the assumed face stresses [T xz , T yz ] correspond to [τ xz0 , τ yz0 ] in Equation 21 consistent with applied or reactive vertical shears along the edges.
The condition ψ 0 = 0 is different from the usual condition w 0 = 0 in Kirchhoff's theory. The function ψ 0 is related to normal strain ε z which is proportional to Δw 0 . If the plate is free of applied transverse stresses, that is, the plate is subjected to bending and twisting moments only, ψ 0 (x, y), thereby, Δw 0 ≡ 0 from Equation 2 (In Kirchhoff's theory, the tangential gradient of Δw 0 is proportional to the corresponding gradient of applied τ xy along the edge of the plate implied from Kelvin and Tait's physical interpretation of the contracted transverse shear condition). This Laplace equation is not adequate to satisfy two in-plane edge conditions. One needs its conjugate harmonic functions implying v 1 , x = u 1 , y to express in-plane displacements in the form The function φ 0 was introduced earlier by Reissner (1945) as a stress function in satisfying Equation 2.
After finding w 0 and φ 0 from solving in-plane equilibrium equations (Equation 1a,b), correction w 0c to w 0 from zero face shear conditions is given by One should note here that w(x, y) = w 0 + w 0c is the face deflection. It does not satisfy the prescribed w = 0 edge condition. As mentioned in Poisson's theory, its zero value along the edge requires only a support to prevent the vertical deflection of intersection of the face with the wall of the plate. Correction w 0c is obviously due to the tangential gradient of τ xy in vertical shear V x (y) . In Kirchhoff's theory, w(x, y) is a bi-harmonic function. It is obtained in the simply supported (hard type) plate satisfying w = 0 along with in-plane condition along the x (and y) constant edge In Kirchhoff 's theory, w(x, y) obtained above corresponds to neutral plane deflection w 0N with edge condition w 0N = 0. Obviously, this edge condition is not required to obtain [u 1 , v 1 ] and face deflection. With reference to the present analysis, it is relevant to note the following observations: In the preliminary solution, σ z , α[τ xz , τ yz ], α 2 [σ x , σ y , τ xy ] are of O(1). As such, estimation of in-plane stresses, thereby, in-plane displacements, are not dependent on w 0 . The support condition w 0 = 0 along the plate edge does not play any role in the estimation of in-plane displacements. In fact, Kirchhoff's theory gives only a lower bound for α for the validity of small deformation theory. Validity of its application is a questionable proposition for small values of α even slightly higher than this lower bound since the theory gives highly constrained solutions for displacements. Shear deformations are included in FSDT in which face shear conditions are not satisfied and Reissner's theory in which edge condition on w 0 is not satisfied. Both of these theories give more or less the same maximum vertical deflection. It shows that satisfying either edge condition on w 0 or face shear conditions has equal effect on deformations. In each of these theories, both conditions are satisfied only in the limit. Limiting process, however, nullifies the effect due to the applied or reactive edge stress τ xy in the bending problem. Moreover, linear variation of σ z in z is neglected in constitutive relations. The present analysis shows that this linear σ z results in more flexible plate deformations in conjunction with zero face shear conditions.
In 3-D problems of homogeneous plates, w 0 (x, y, 0) and w 0 (x, y, 1) are the neutral plane deflection w 0N and face deflection w 0F , respectively. Both of them are 2-D functions and they are from the thickness-wise integration of ε z in a 2-D problem. Determination of w 0N used as an unknown variable is dependent on edge support condition and that of w 0F on zero face shear conditions. The present analysis is with reference to finding w 0F and the need to satisfy that edge condition on w 0N does not exist.
[u 1c , v 1c ] due to transverse stresses [τ xz0 , τ yz0 , zσ z1 ] With reference to the second problem mentioned earlier, the in-plane displacements to be determined are z(u 1 , v 1 ) c . Here transverse stresses are expressed in the form in which (τ xz2 , τ yz2 , σ z3 ) c yet to be determined are corrections due to (u 1 , v 1 ) c . Since the added terms are dependent on material constants, the above transverse shear stresses are no longer applied stresses along the edges.
In order to keep σ z3c as a free variable, the function f 3 (z) in Equation 29 is replaced by f 3 *(z) = (f 3 − z/3) so that f 3 *(±1) = 0. Then the face load condition is satisfied with In further analysis, it is convenient to express σ z with β 1 = 1/3 in the form In-plane distributions [u 1 , v 1 ] c are modified in the form with The introduction of w 0 (x, y) term in the above equation is to avoid its participation in the static in-plane equilibrium equations. Modified transverse shear stresses from stress-strain and strain-displacement relations are Substitution of the above shear stresses in Equation 1b with σ z = z[q/2 − β 1 σ z3c ] and e 1c = (ε x1 + ε y1 ) c gives For the use of [u 1 *, v 1 *] c in the integration of equilibrium equations (Equation 1a,b), displacements [u 1 , v 1 ] c due to the requirement of v 1c,x = u 1c,y to decouple bending and torsion problems are expressed in the form in the case of hard and/or soft simply supported edge conditions. Contributions of ψ 1c and w 0 in [u 1 , v 1 ] c * are one and the same in giving corrections to w(x, y, z) and transverse stresses (in fact, the contribution of w 0 is through strain-displacement relations in static equilibrium equations and through constitutive relations in throughthickness integration of equilibrium equations). Since w 0 does not participate in the in-plane static equations, its contribution is through [u 1 , v 1 ] c in the integrated in-plane equilibrium equations. Hence, w 0 in [u 1 , v 1 ] c * is replaced by ψ 1c (so as to be independent of w 0 used in strain-displacement relations) so that [u 1 , v 1 , ε x1 , ε y1 , γ xy1 ] c * are From the integration of Equation 1a using the above strains in Equation 12a along with v 1 , x = u 1 , y , reactive transverse stresses are One equation governing in-plane displacements (u 1 , v 1 ), noting that σ z3 from Equation 29 is negative of the one from Equation 34 due to (f 3 , zz + f 1 ) = 0, is given by By equating σ z3s (from the integration of Equation 1a,b with s variables) with β 1 σ z3 (from static equations with * variables), one equation governing [u 1s , v 1s By expressing [u 1s , v 1s ] = −α [ψ 1s , x , ψ 1s , y ] due to the second equation v 1s , x = u 1s , y , Equation 58 becomes a fourth-order equation in ψ 1s to be solved with two inplane conditions along constant x (and y) edges in the present example ψ 1s is related to e 1c by the equation c 1s is related to c 1c by the equation so that c 1s = 0.55(2q 0 /E) with ν = 0.3 and α = 1/6. Neutral plane deflection w 0N = 4.36(2q 0 /E), (=2.27 + 0.15 + 1.39 + 0.55) is underestimated by 2.90% from the exact value (4.49). Evaluation of neutral plane deflection involves u 1c and u 1s associated with f 2 (z) and cos (πz/2) distributions of transverse shear stresses, respectively. As such, the solution of the present auxiliary problem provides secondorder corrections to the solutions of the primary problem. It is to be noted that the wide variation of percentage of errors from 0.67 to 4.80 in w along the vertical line from the neutral plane to the face plane at the first stage of iteration using (u 3 , v 3 ) in Poisson's theory is reduced from 2.16 to 2.90 from the face plane to the neutral plane without using (u 3 , v 3 ) in extended Poisson's theory. One iteration using (u 3 , v 3 ) is expected to reduce this gap even further along with decreased percentage of errors. A significant implication of this observation is that solution of the present auxiliary problem is necessary so as to obtain a more or less uniform approximation to thickness-wise distributions of deformations.
The fourteenth-order system of equations in the analysis of a simply supported plate presented here forms the initial set in the iterative procedure. Moreover, the analysis clearly establishes that the second-order corrections in Kirchhoff 's theory are mainly due to inclusion of σ z = z σ z1 in the constitutive relations and not due to consideration of shear deformations presumed in FSDT and Reissner's theory. It is possible to decrease the percentage of errors in the estimation of w(x, y z) by considering the influence of in-plane displacements thus obtained in estimation of transverse stresses. Successive iterative steps consist of a pair of fourth-order system of equations forming a proper sequence of sets of 2-D problems converging to a 3-D problem.

Conclusions
Estimation of transverse stresses in the preliminary solution is independent of material constants through Poisson's theory in the analysis of primary bending problem defined from Kirchhoff 's theory. Hence, it is unaltered in the analysis of bending of homogeneous and laminated plates with orthotropic and anisotropic materials. Application of transverse shear stress independent of thickness coordinate z along the edges of the plate is more practical instead of stress resultants which have no unique z-distributions. Coupling of w 0 (x, y) in obtaining in-plane displacements which are the root cause for the Poisson-Kirchhoff boundary condition paradox is eliminated. In addition, inclusion of the tangential gradient of applied inplane shear stress in transverse shear in Kirchhoff's theory is removed. Vertical stress σ z linear in z neglected in Kirchhoff's theory and shear deformation theories gives primary second-order corrections instead of shear deformations. Thickness-wise distributions of displacements in terms of polynomials in z are not adequate to reduce 3-D problems into a sequence of 2-D problems.