We investigate the dynamics of the electroweak phase transition within an extended Standard Model framework that includes one real scalar (Φ) and one complex scalar (S), both of which are SM gauge singlets. The global U(1) symmetry is softly broken to a \( {\mathcal{Z}}_3 \) symmetry by the S3 term in the scalar potential. After this U(1) symmetry breaking, the imaginary component of the complex scalar (S) acts as a pseudo-Nambu-Goldstone boson (pNGB) dark matter candidate, naturally stabilized by the \( {\mathcal{Z}}_2 \) symmetry of the scenario. Specially, the spontaneous breaking of the global U(1) symmetry to a discrete \( {\mathcal{Z}}_3 \) subgroup can introduce effective cubic terms in the scalar potential, which facilitates a strong first-order phase transition. We analyze both single-step and multi-step first-order phase transitions, identifying the parameter space that satisfies the dark matter relic density constraints, complies with all relevant experimental constraints, and exhibits a strong first-order electroweak phase transition. The interplay of these criteria significantly restricts the model parameter space, often leading to an underabundant relic density. Moreover, we delve into the gravitational wave signatures associated with this framework, offering valuable insights that complement traditional dark matter direct and indirect detection methods.
