Abstract In this work, we have solved the radial Schrö-dinger equation for the Woods–Saxon potential together with coulomb (r>Rc), centrifugal terms and spin orbit in-teraction by using a new type of Nikiforov-Uvarov (NU) method. This approach is based on the solution of the Second-Order Linear Differential Equations (SOLDE). The mandatory specific choices of the required parameters in this technique restricts the application of this method to the Schrödinger equation with complicated potential pro-files, which means that the NU method cannot efficiently be employed to solve more realistic physical systems. Due to the mentioned difficulties in evaluating the equivalent second order algebraic equation in the NU method, the analytical NU method has to be extended to more effi-cient version that is combined with numerical methods (that leads to a semi-analytical method). We have solved it by combination of the NU method with the numerical fitting schema. The numerical fitting schema helps us to find the mentioned second order algebraic equation. Oth-erwise, complicated changes of variables or overwhelming algebraic treatments to deriving the energy eigenvalues and the wavefunctions are required. The current approach is simpler, more flexible and efficient. This technique can also be developed to be suitable for the equations other than the Schrodinger one. The Woods–Saxon potential is also a short-range interaction in the potential model for nuclear physics and has predictions for the nuclear shell model and distribution of nuclear densities. We have ob-tained a semi-analytical energy eigenvalues and eigen-functions for various values of n, l, and j quantum num-bers. Agreement of 5/2+ and 1/2+ wavefunctions with the published works is also obtained which also shows the ac-curacy of our method.