In the high Deborah (De) or Weissenberg (Wi) number flow, i.e., highly elastic flow, the presence of the elastic instability for the benchmark problem of viscoelastic fluid, 4:1 abrupt contraction flow, has been reported. The flow instability leads to the melt fracture and shark skin feature of fiber spinning, i.e., the polymer spinning. Furthermore, this instability results in quality degradation during the transportation of the battery slurry, which is a representative viscoelastic material in chemical processes. To solve the exponentially increasing stress at high De, log conformation reformulation (LCR) method, which preserves the positive definiteness of conformation tensor without stabilization methods, has been proposed. The elastic instability, or flow instability, has been reported only when the finite volume method (FVM) was used at high De with simple viscoelastic constitutive equation, Oldroyd-B model. In contrast, a stable solution can be obtained up to a relatively low De number with the finite element method (FEM) without any instability. In this study, we compared FEM-based numerical results about the instability to FVM-based results for validation and expanded to 3D contraction flow instabilities of viscoelastic fluid by implementing the LCR method to FEM. Furthermore, we analyzed the instabilities for complex viscoelastic constitutive equations such as Giesekus and Leonov models. Here, TR-FDI with adaptive time stepping was used to prevent accumulated numerical errors from time integration. Above some critical De number, flow changes unsteady with oscillating vortex, and we tracked the time evolution of vortex length. We compared FEM-based results to our experimental data using the 4:1 abrupt contraction pipe system and model fluids from low to high De number flow to examine the transient vortex size and pressure drop before and after the contraction.