In this work, we present an efficient numerical implementation of the finite element method for modal analysis that leverages various symmetry operations, including spatial symmetry in point groups and space-time symmetry in pseudo-hermiticity systems. We provide a formal and rigorous treatment, specifically deriving the boundary constraint conditions corresponding to symmetry constraints. Without loss of generality, we illustrate our approach via computing the modes of optical waveguides with complex cross-sections, accompanied by performance benchmark against the standard finite element method. The obtained results demonstrate excellent agreement between our method and standard FEM with significantly improved computational efficiency. Specifically, the calculation speed increased by a factor of 23 in the hollow-core fiber. Furthermore, our method directly classifies and computes the modes based on symmetry, facilitating the modal analysis of complex waveguides.