<p>This paper presents new shell elements that require only the shell midsurface to be discretized, utilize complete three-dimensional constitutive laws, and contain only seven translational degrees of freedom along each nodal normal. The middle surface of such a shell element is approximated by a variable-order quadrilateral spectral element, which is scaled along the thickness to represent the full shell geometry. Along the nodal normal, the in-plane displacement fields are assumed to be linear by interpolating the displacements at the top and bottom points, while the transverse displacement field is quadratic by interpolating the displacements at the top, middle, and bottom points. Lagrange polynomials are employed for the interpolation and the coefficients of the linear and quadratic terms are respectively scaled by the shell thickness and its square to reduce the condition number of the global stiffness matrix particularly for thin shells. The inverse Jacobian matrix is approximated as a quadratic matrix polynomial using Neumann expansion to guarantee the solution accuracy for thick curved shells. The assumed natural strain method is applied to alleviate transverse shear locking, membrane locking, and curvature thickness locking. A general procedure for recovering the transverse stresses is proposed by analytically integrating the static equilibrium equations along the thickness direction, and the involved functional derivatives with respect to the in-plane coordinates are evaluated via the differential quadrature method. The capabilities of the developed shell elements are sufficiently demonstrated through numerical examples.</p>