Through the introduced work a viable numerical scheme has been developed, using Numerov's method, to study the energy spectra of mesons and hadronic interactions. Upon the implementation of the non-relativistic quark model and the Quantum Chromodynamics theory elements which describe the phenomenological interactions between the charm-anticharm quarks via two proposed potentials, our model accurately predicts the mass spectra of charmed quarkonium as an example of mesonic systems. It was found that our results agree well with the experimental charmonium masses that were published via the Particle Data Group (PDG). The proposed model has been tested for numerical instabilities and found that it works well for both proposed potentials. Although this model is suited to study the dynamics of charmonium interactions, we found that the results for one of the potentials are in better agreement with the published data. In addition, we predicted the mass spectrum of new charmed multiplets using the same previously mentioned potentials.