Introduction

In the recent time, physicists have developed much interest in searching for the exponential-type potentials. The reason is that most of the exponential-type potentials play an important role in physics, e.g. Yukawa potential is used in plasma, solid-state and atomic physics [1]. As a result, many authors have solved both relativistic and non-relativistic wave equations with these potentials. For instance, Zhang et al. [2] obtained approximate solutions of the Schrödinger equation with the Generalized Morse potential model including a centrifugal term. Jia et al. [3] solved six-parameter exponential-type potentials. Onate [4] obtained relativistic and non-relativistic solutions of the inversely quadratic Yukawa potential. Falaye et al. [5] obtained bound state solutions of the Schrödinger equation with Manning–Rosen potential. Hassanabadi et al. [6, 7] obtained Actual and general Manning–Rosen potential under spin and pseudospin symmetries of the Dirac equation; approximate solutions of the Schrödinger equation under Manning–Rosen potential in arbitrary dimensions via SUSY QM. Hamzavi et al. [8] obtained approximate spin and pseudospin solutions of the Dirac equation for inversely quadratic Yukawa potential and tensor interaction. Maghsoodi et al. [1] solved Dirac particles in the presence of Yukawa potential plus a tensor interaction in SUSY QM frame work. Ikhdair [9] obtained on the bound state solutions of the Manning–Rosen potential including an improved approximation to the orbital centrifugal term.

The solutions of the wave equations (either Schrödinger, Klein–Gordon, Dirac or D.K.P.) with any of these exponential-type potentials are obtained using different methods which include: asymptotic iteration method (AIM) [10–16], Nikiforov–Uvarov (N.U) method [17–20], exact/proper quantization rule [21], supersymmetric method [22–27], 1/N shifted expansion method [28], etc.

Motivated by the success in the exponential-type potentials, we attempt to investigate the solutions of the radial Schrödinger equation with a class of Yukawa potentials given as

$$V\left( r \right) = \frac{{ - br + rce^{ - \delta r} - ae^{ - 2\delta r} }}{{r^{2} }},$$
(1)

where a, b and c are potential strength and δ is the screening parameter. The potential is obtained by the addition of Hellmann potential and inversely quadratic Yukawa potential which is equivalent to the interaction of inversely quadratic Yukawa potential, Coulomb potential and Yukawa potential. This potential has its application where its components are useful. It is noted that the exact solution of the radial Schrödinger equation with potential (1) is not possible due to the presence of the inverse square term. Therefore, to obtain an approximate solutions, we employ a suitable approximation scheme. It is found that such approximation proposed by Greene and Aldrich [29]

$$\frac{1}{{r^{2} }} = \frac{{\delta^{2} }}{{\left( {1 - e^{ - \delta r} } \right)^{2} }}.$$
(2)

is a good approximation to the centrifugal/inverse square term which is valid for \(\delta \ll 1\) for a short potential range.

Our work is arranged as follows: in the following section, we obtain the bound state solutions. We present discussion and concluding remarks at the end of the article.

Bound state solutions

To study any quantum physical system, we solve the original Schrödinger equation given as [9, 30, 31]

$$\left( {\frac{{P^{2} }}{2\mu } + V\left( r \right)} \right)\psi_{n,\ell ,m} \left( r \right) = E_{n,\ell } \psi_{n,\ell ,m} \left( r \right),$$
(3)

where \(E_{n,\ell }\) is the energy, μ is the particle mass and V(r) is the potential. Setting the wave functions \(\psi_{n,\ell ,m} \left( r \right) = \frac{{R_{n,\ell } \left( r \right)Y_{\ell ,m} \left( {\theta ,\varphi } \right)}}{r}\), we obtained the following radial Schrödinger equation

$$\left[ {\frac{{{\text{d}}^{2} }}{{{\text{d}}r^{2} }} + \frac{2\mu }{{\hbar^{2} }}\left( {E_{n,\ell } - V\left( r \right)} \right) - \frac{{\ell \left( {\ell + 1} \right)}}{{r^{2} }}} \right]R_{n,\ell } \left( r \right) = 0.$$
(4)

With Eqs. (1) and (2), (4) becomes

$$\frac{{{\text{d}}^{2} R_{n,\ell } \left( r \right)}}{{{\text{d}}r^{2} }} = \left[ {V_{\text{eff}} - E_{\text{eff}} } \right]R_{n,\ell } \left( r \right),$$
(5)

where \(V_{\text{eff}} = \left[ {\frac{2\mu \delta }{{\hbar^{2} }}\left( {a\delta - b + c} \right) + \ell \left( {\ell + 1} \right)\delta^{2} } \right]\frac{{e^{ - \delta r} }}{{1 - e^{ - \delta r} }} - \left[ {\frac{{2\mu a\delta^{2} }}{{\hbar^{2} }} - \ell \left( {\ell + 1} \right)\delta^{2} } \right]\frac{{e^{ - \delta r} }}{{\left( {1 - e^{ - \delta r} } \right)^{2} }}\) and \(- E_{\text{eff}} = - \frac{2\mu \delta }{{\hbar^{2} }}\left( {E_{n,\ell } + b} \right) + \ell \left( {\ell + 1} \right)\delta^{2}\). Since we are dealing with a Schrödinger-like equation that we solved by means of SUSY QM [6, 32–35], the first step in the SUSY approach is finding the solution of the Riccati equation [6]. Using the shape invariance formalism, it can be easily seen that

$$W^{2} (r) - \frac{{{\text{d}}W(r)}}{{{\text{d}}r}} = V_{\text{eff}} - E_{\text{eff}} ,$$
(6)

whose solution is given as

$$W(r) = A - \frac{{Be^{ - \delta r} }}{{1 - e^{ - \delta r} }},$$
(7)

The ground state wave function \(R_{0,\ell } (r)\) as

$$R_{0,\ell } (r) = \exp \left( { - \int {W(r){\text{d}}r} } \right),$$
(8)

where W(r) is called the superpotential function in supersymmetric quantum mechanics [36–38]. Substitute Eq. (7) into Eq. (6), we deduce the following relations

$$A^{2} = - E_{\text{eff}} ,$$
(9)
$$B = \frac{2\mu a\delta }{{\hbar^{2} }} - \left( {\ell + 1} \right)\delta ,$$
(10)
$$A = \frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - B^{2} }}{2B}.$$
(11)

With Eq. (7), we can construct the two partner potentials as follows:

$$U_{ + } = W^{2} (r) + \frac{{{\text{d}}W(r)}}{{{\text{d}}r}} = A^{2} - \frac{{B(B + 2A)e^{ - \delta r} }}{{1 - e^{ - \delta r} }} + \frac{{B(B + \delta )e^{ - \delta r} }}{{\left( {1 - e^{ - \delta r} } \right)^{2} }},$$
(12)
$$U_{ + } = W^{2} (r) + \frac{{{\text{d}}W(r)}}{{{\text{d}}r}} = A^{2} - \frac{{B(B + 2A)e^{ - \delta r} }}{{1 - e^{ - \delta r} }} + \frac{{B(B - \delta )e^{ - \delta r} }}{{\left( {1 - e^{ - \delta r} } \right)^{2} }}.$$
(13)

Using the shape invariance technique [6, 39–42], it can readily be shown that the two partner potentials are shape invariant. Therefore, their relationship is written as

$$U_{ + } (r,a_{0} ) = U_{ - } (r,a_{1} ) + T(a_{1} ),$$
(14)

where a 0 is an old set of parameters in which the new set of parameters a 1 is obtained from and T(a 1) is a remainder that is independent of the variable r. Here, B = a 0 via mapping of the form a 1 → a 0 − δ,  a 2 → a 0 − 2δ,  a 3 → a 0 − 3δ. Thus, a generalization is drawn as a n  → a 0 − nδ. In terms of the parameters of the problem, we obtain the following relations:

$$T(a_{1} ) = \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{0}^{2} }}{{2a_{0} }}} \right] - \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{1}^{2} }}{{2a_{1} }}} \right],$$
(15)
$$T(a_{2} ) = \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{1}^{2} }}{{2a_{1} }}} \right] - \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{2}^{2} }}{{2a_{2} }}} \right],$$
(16)
$$T(a_{3} ) = \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{2}^{2} }}{{2a_{2} }}} \right] - \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{3}^{2} }}{{2a_{3} }}} \right],$$
(17)
$$T(a_{n} ) = \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{n - 1}^{2} }}{{2a_{n - 1} }}} \right] - \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{n}^{2} }}{{2a_{n} }}} \right],$$
(18)

The energy spectral can then be determine as follows:

$$E_{0,\ell }^{ - } = 0,$$
(19)
$$E_{n,\ell }^{ - } = E_{\text{eff}} + E_{0,\ell }^{ - } = \sum\limits_{k = 1}^{n} {T(} a_{k} ) = \left[ {\frac{{2\mu \delta \left( {b - c - a\delta } \right) - \ell \left( {\ell + 1} \right)\delta^{2} - a_{n}^{2} }}{{2a_{n} }}} \right].$$
(20)

This gives energy equation as

$$E_{n,\ell } = \delta \left( {\frac{{\delta \hbar^{2} \ell (\ell + 1)}}{2\mu } - b} \right) - \frac{{\delta^{2} \hbar^{2} }}{2\mu }\left[ {\frac{{\frac{2\mu }{{\delta \hbar^{2} }}\left( {b - c - a\delta } \right) - \left( {n + \ell + 1} \right)^{2} - \ell \left( {\ell + 1} \right)}}{{2\left( {\ell + n + 1} \right)}}} \right]^{2} .$$
(21)

Wave function

To obtain the un-normalized wave function, we define a variable of the form \(y = e^{ - \delta r}\) and substitute it into Eq. (5) to have

$$\frac{{{\text{d}}^{2} R_{n,\ell } (r)}}{{{\text{d}}y^{2} }} + \frac{1}{y}\frac{{{\text{d}}R_{n,\ell } (r)}}{{{\text{d}}y}} + \frac{{Ny^{2} Qy + P}}{{\left( {y\left( {1 - y} \right)} \right)^{2} }}R_{n,\ell } (r) = 0,$$
(22)

where

$$\begin{aligned} P & = \frac{{2\mu (b\delta + E_{{n,\ell }} )}}{{\delta ^{2} \hbar ^{2} }} - \ell (\ell + 1),\; \\ Q & = - \frac{{2\mu (b\delta + c\delta + E_{{n,\ell }} )}}{{\delta ^{2} \hbar ^{2} }}, \\ N & = \frac{{2\mu (a\delta + c\delta ^{2} + 2E_{{n,\ell }} )}}{{\delta ^{2} \hbar ^{2} }} \\ \end{aligned}$$
(23)

Analyzing the asymptotic behavior of Eq. (22) at origin and at infinity, it can be tested that when r → 0(y → 1) and when r → ∞(y → 0) Eq. (22) has a solution

$$R_{n,\ell } (y) = y^{\eta } \left( {1 - y} \right)^{\varepsilon } ,$$
(24)

where

$$\eta = \sqrt { - \frac{{2\mu (b\delta + E_{n,\ell } )}}{{\delta^{2} \hbar^{2} }} + \ell (\ell + 1)} ,$$
(25)
$$\varepsilon = \frac{1}{2} + \frac{1}{2}\sqrt {1 - \frac{2\mu a}{{\delta \hbar^{2} }} + \ell (\ell + 1)} .$$
(26)

Now, taking a trial wave function of the form \(R_{n,\ell } (y) = y^{\eta } \left( {1 - y} \right)^{\varepsilon } f(y)\) and substitute it into Eq. (22), we have

$$f^{''} (y) + f^{'} (y)\frac{2\eta + 1 - y(2\eta + 2\varepsilon + 1)}{y(1 - y)} + f(y)\frac{{\left( {\eta + \varepsilon } \right)^{2} + N}}{y(1 - y)} = 0.$$
(27)

Equation (27) is a differential equation satisfied by the hypergeometric function. Thus, its solution is obtained as:

$$f(y) = {}_{2}F_{1} \left( { - n,n + 2(\eta + \varepsilon );2\eta + 1,y} \right).$$
(28)

Replacing the function f(y) with the hypergeometric function and write a complete radial wave function as:

$$R_{n,\ell } (y) = y^{\eta } \left( {1 - y} \right)^{\varepsilon } {}_{2}F_{1} \left( { - n,n + 2(\eta + \varepsilon );2\eta + 1,y} \right).$$
(29)

Results and discussion

Some special cases of interest are studied here, when a = 0,  our potential (1) reduces to Hellmann potential

$$V\left( r \right) = \frac{{ - b + ce^{ - \delta r} }}{r},$$
(30)

which reduces the energy Eq. (21) to

$$E_{n,\ell } = \delta \left( {\frac{{\delta \hbar^{2} \ell (\ell + 1)}}{2\mu } - b} \right) - \frac{{\delta^{2} \hbar^{2} }}{2\mu }\left[ {\frac{{\frac{2\mu }{{\delta \hbar^{2} }}\left( {b - c} \right) - \left( {n + \ell + 1} \right)^{2} - \ell \left( {\ell + 1} \right)}}{{2\left( {\ell + n + 1} \right)}}} \right]^{2} .$$
(31)

When a = c = 0, the potential (1) reduces to Coulomb potential

$$V\left( r \right) = \frac{ - b}{r},$$
(32)

and the energy Eq. (21) turns to

$$E_{n,\ell } = \delta \left( {\frac{{\delta \hbar^{2} \ell (\ell + 1)}}{2\mu } - b} \right) - \frac{{\delta^{2} \hbar^{2} }}{2\mu }\left[ {\frac{{\frac{2\mu b}{{\delta \hbar^{2} }} - \left( {n + \ell + 1} \right)^{2} - \ell \left( {\ell + 1} \right)}}{{2\left( {\ell + n + 1} \right)}}} \right]^{2} .$$
(33)

When a = b = 0, the potential (1) reduces to Yukawa potential

$$V\left( r \right) = \frac{{ - ce^{ - \delta r} }}{r},$$
(34)

and the energy Eq. (21) turns to

$$E_{n,\ell } = \frac{{\delta^{2} \hbar^{2} \ell (\ell + 1)}}{2\mu } - \frac{{\delta^{2} \hbar^{2} }}{2\mu }\left[ {\frac{{\frac{ - 2\mu c}{{\delta \hbar^{2} }} - \left( {n + \ell + 1} \right)^{2} - \ell \left( {\ell + 1} \right)}}{{2\left( {\ell + n + 1} \right)}}} \right]^{2} .$$
(35)

When b = c = 0, potential (1) reduces to inversely quadratic Yukawa potential

$$V\left( r \right) = \frac{{ - ae^{ - 2\delta r} }}{{r^{2} }},$$
(36)

and the energy Eq. (21) turns to

$$E_{n,\ell } = \frac{{\delta^{2} \hbar^{2} \ell (\ell + 1)}}{2\mu } - \frac{{\delta^{2} \hbar^{2} }}{2\mu }\left[ {\frac{{\frac{ - 2\mu a}{{\hbar^{2} }} - \left( {n + \ell + 1} \right)^{2} - \ell \left( {\ell + 1} \right)}}{{2\left( {\ell + n + 1} \right)}}} \right]^{2} .$$
(37)

In Table 1, we numerically reported the energy eigenvalues for a class of Yukawa potential (\(a = 1 \times 10^{ - 5} ,b = 2 \times 10^{0} , \quad c = - 1 \times 10^{0} )\), Hellmann potential (\(a = 0, b = 2 \times 10^{0} , \;c = - 1 \times 10^{0} )\), Coulomb potential (\(a = 0 = c = 0,\quad b = 3)\) and Yukawa potential (\(a = b = 0, \;c = - 3 \times 10^{0} )\) for 2p, 3p, 3d, 4p, 4d and 4f. In Table 2, we have reported the energy eigenvalues of these potentials for \(n = 0,\quad \in \ell = 0;\) \(n = 1,\quad \ell = 0,1;\) \(n = 2,\; \ell = 0, 1, 2;\) \(n = 3, \; \in \ell = 0, 1, 2, 3\) and \(n = 4,\; \ell = 0, 1, 2, 3, 4.\) In Table 1, energy increases as the screening parameter increases. In Table 2, energy increases as n increases. In Table 3, we compared our result for the Hellmann potential with the result from two other methods.

Table 1 Bound states energy for the class of Yukawa potential \(E_{n,\ell }^{\text{CYP}}\), Hellmann potential \(E_{n,\ell }^{\text{HP}}\), Coulomb potential \(E_{n,\ell }^{\text{CP}}\) and Yukawa potential \(E_{n,\ell }^{\text{YP}}\)
Full size table
Table 2 Bound states energy for the class of Yukawa potential \(E_{n,\ell }^{\text{CYP}}\), Hellmann potential \(E_{n,\ell }^{\text{HP}}\), Coulomb potential \(E_{n,\ell }^{\text{CP}}\) and Yukawa potential \(E_{n,\ell }^{\text{YP}}\)
Full size table
Table 3 Ro-vibrational energy spectrum \(\left( { - E_{n,\ell } } \right)\) for the Hellmann potential with \(2m = \hbar = 1,\) \(b = 2\) and \(c = 1\)
Full size table

Conclusions

We have obtained approximate solutions of the Schrödinger equation by combining inversely quadratic Yukawa potential, Yukawa potential and Coulomb potential. We deduced that from the energy equation of these combined potentials, the energy equation of Hellmann potential as well as these individual potential can be obtained. In Tables 1 and 2, we have numerically reported the equivalence of the energy of these potentials.