Accurate modeling of composite structures is important for their safe application under different loading conditions. To provide accurate predictions of three-dimensional (3D) stress fields in an efficient computational framework, in this study, a modeling approach that builds upon the recently developed hierarchical Serendipity Lagrange finite elements (FEs) is employed. The approach provides layerwise (LW) and equivalent single-layer (ESL) models for analyzing constant- and variable-stiffness laminated beam structures. To enhance the capability of the ESL model, two zig-zag (ZZ) functions, namely,Murakami’s ZZ function (MZZF) and the refined ZZ theory function (RZT), are implemented. For constant-stiffness laminated and sandwich beams, the RZT ZZ function predicts the structural response more accurately than the MZZF. However, for variable-stiffness laminated structures, the applicability of RZT remains unknown and its accuracy is therefore tested within the present formulation. Results obtained are validated against 3D closed-form and 3D FE solutions available from the literature. For similar levels of accuracy, significant gains in computational efficiency are achieved over 3D FE and LW models by using the ESL approach with RZT ZZ functions.