Accurate, efficient and versatile computational methods for functionally graded piezoelectric (FGP) shells are still restricted. This contribution proposes new FGP shell finite elements with superior comprehensive performance. In the formulation, a shell element is treated as a three-dimensional continuum and its middle surface is represented with a quadrilateral spectral element. Along the thickness, the shell geometry is described by scaling the middle surface while the displacements and electric potential are approximated via a consistent quadratic Lagrange interpolation. The electric potential is scaled by a factor λ to avoid numerical instability and the assumed natural strain method is applied to eliminate numerical locking. The piezoelectric material properties are assumed to vary continuously through the shell thickness obeying a simple power-law function. Binomial theorem is utilized to perform analytical integral along the thickness in evaluating the stiffness and mass matrices and nodal load vector. The validity of the developed approach is demonstrated by solving piezoelectric or FGP plate problems with reference solutions. A simply supported FGP square plate, a fully clamped FGP cylindrical shell and a partly clamped FGP hyperbolic paraboloid shell are systematically analyzed. The influence of the power-law index and span-to-thickness ratio on the static and free vibration behaviors of the FGP structures is investigated. The optimal value of λ for general FGP shells is also determined.