Abstract Epidemiological models feature reliable and valuable insights into the prevention and transmission of life-threatening illnesses. In this study, a novel SIR mathematical model for COVID-19 is formulated and examined. The newly developed model has been thoroughly explored through theoretical analysis and computational methods, specifically the continuous Galerkin–Petrov (cGP) scheme. The next-generation matrix approach was used to calculate the reproduction number ( R 0 ) ({R}_{0}) . Both disease-free equilibrium (DFE) ( E 0 ) ({E}^{0}) and the endemic equilibrium ( E ⁎ ) ({E}^{\ast }) points are derived for the proposed model. The stability analysis of the equilibrium points reveals that ( E 0 ) ({E}^{0}) is locally asymptotically stable when R 0 < 1 {R}_{0}\lt 1 , while E ⁎ {E}^{\ast } is globally asymptotically stable when R 0 > 1 {R}_{0}\gt 1 . We have examined the model’s local stability (LS) and global stability (GS) for endemic equilibrium \text{&#x00A0;} and DFE based on the number ( R 0 ) ({R}_{0}) . To ascertain the dominance of the parameters, we examined the sensitivity of the number ( R 0 ) ({R}_{0}) to parameters and computed sensitivity indices. Additionally, using the fourth-order Runge–Kutta (RK4) and Runge–Kutta–Fehlberg (RK45) techniques implemented in MATLAB, we determined the numerical solutions. Furthermore, the model was solved using the continuous cGP time discretization technique. We implemented a variety of schemes like cGP(2), RK4, and RK45 for the COVID-19 model and presented the numerical and graphical solutions of the model. Furthermore, we compared the results obtained using the above-mentioned schemes and observed that all results overlap with each other. The significant properties of several physical parameters under consideration were discussed. In the end, the computational analysis shows a clear image of the rise and fall in the spread of this disease over time in a specific location.