A singular one-dimensional bound state problem and its degeneracies

A Singular One-Dimensional Bound State Problem and its Degeneracies Fatih Erman1, Manuel Gadella2, Se¸cil Tunalı3, Haydar Uncu4 1Department of Mathematics, ˙ Izmir Institute of Technology, Urla, 35430, Izmir, Turkey 2Departamento de F´ısica Te´orica, At´omica y ´ Optica and IMUVA. Universidad de Valladolid, Campus Miguel Delibes, Paseo Bel´en 7, 47011, Valladolid, Spain 3Department of Mathematics, ˙ Istanbul Bilgi University, Dolapdere Campus 34440 Beyo˘glu, ˙ Istanbul 4Department of Physics, Adnan Menderes University, 09100, Aydın, Turkey E-mail: [email protected], man[email protected], [email protected], [email protected] June 8, 2017 Abstract We give a brief exposition of the formulation of the bound state problem for the one-dimensional system of Nattractive Dirac delta potentials, as an N×Nmatrix eigenvalue problem (ΦA=ωA). The main aim of this paper is to illustrate that the non-degeneracy theorem in one dimension breaks down for the equidistantly distributed Dirac delta potential, where the matrix Φ becomes a special form of the circulant matrix. We then give an elementary proof that the ground state is always non-degenerate and the associated wave function may be chosen to be positive by using the Perron-Frobenius theorem. We also prove that removing a single center from the system of Ndelta centers shift all the bound state energy levels upward as a simple consequence of the Cauchy interlacing theorem. Keywords. Point interactions, Dirac delta potentials, bound states. 1 Introduction Dirac delta potentials or point interactions, or sometimes called contact potentials are one of the exactly solvable classes of idealized potentials, and used as a pedagogical tool to illustrate various physically important phenomena, where the de Broglie wavelength of the particle is much larger than the range of the interaction. They have various applications in almost all areas of physics, see e. g., [1] and [2], and references therein. For instance, mutually non-interacting electrons moving in a fixed crystal can be modeled by periodic Dirac delta potentials, known as the Kronig-Penney model [3]. Another application is given by the model consisting of two attractive Dirac delta potentials in one dimension. This is used as a very elementary model of the chemical bound for a diatomic ion (H+ 2, for example) and has been discussed in [4, 5]. The interest of Dirac delta potentials and other one dimensional point potentials provides us with solvable (or quasi-solvable) models in quantum mechanics that give an insight for a better understanding of the basic features of the quantum theory. This makes them suitable for the purpose of teaching the discipline. In the recent pedagogical review [2], several interesting features of one-dimensional Dirac 1 delta potentials have been illustrated and multiple δ-function potential has been studied in Fourier space. Moreover, the bound state problem has been formulated in terms of a matrix eigenvalue problem. In this paper, we first give a brief review of the bound state spectrum of the NDirac delta potentials in one dimension by converting the time independent Schr¨odinger equation Hψ =Eψ for the bound states to the eigenvalue problem for an N×NHermitian matrix. This method is rather useful especially when we deal with large number of centers since the procedure that uses the matching conditions for the wave function at the location of the delta centers become cumbersome for large values of N. Once we formulate the problem as a finite dimensional eigenvalue problem, we show that there are at most Nbound states for Ncenters, using Feynman-Hellmann theorem (see page 288 in [6]). One of the main purposes of this paper is to show that this simple one-dimensional toy problem for more than three centers allows us to give an analytical example of the breakdown of the well-known non-degeneracy theorem for onedimensional bound state problems[7]. This shows us that we should not take the non-degeneracy theorem for granted particularly for singular interactions, where it may not be valid. This was first realized for the so-called one-dimensional Hydrogen atom [8], where the non-degeneracy theorem breaks down and has been studied for other one-dimensional singular potentials since then [9, 10, 11, 12, 13, 14]. In contrast to the degeneracies that appear in bound states, we give an elementary proof that the ground state is non-degenerate and the ground state wave function can be always chosen real-valued and strictly positive. In addition, we also show that all the bound state energies for Nattractive Dirac delta potentials increase if we remove one center from the system. All these results mentioned above become more transparent using some basic theorems from linear algebra, namely Perron-Frobenius theorem, the Cauchy interlacing theorem [15]. The simple proof of these theorems are given in Appendices so as not to interrupt the flow of the presentation. Our presentation is kept simple so that it is also accessible to a wide audience. 2 Bound States for NDirac Delta Potentials We consider a particle moving in one dimension and interacting with the attractive NDirac delta potentials located at aiwith strengths λi>0, where i= 1,2,··· , N. The time independent Schr¨odinger equation is then given by −~2 2m d2ψ dx2− N X i=1 λiδ(x−ai)ψ(x) = Eψ(x).(1) The above equation is actually a formal expression and its exact meaning can only be given by self-adjoint extension theory [16, 17, 18]. Here we follow a more traditional and heuristic approach used in the most quantum mechanics textbooks since the results that we obtain is completely consistent with the rigorous approach. As is well-known, the above equation can also be written as Hψ =Eψ (2) in the operator form, where H=P2 2m+Vand the potential energy operator Vfor the above particular case in the bra-ket formalism is V=− N X i=1 λi|aiihai|.(3) 2 Here |aiiis the position eigenket. In the coordinate representation, the action of Von the state vector |ψiis (V ψ)(x) = hx|V|ψi= N X i=1 λiδ(x−ai)ψ(ai) = N X i=1 λiδ(x−ai)ψ(x),(4) where we have used the fact δ(x−ai)ψ(ai) = δ(x−ai)ψ(x). This justifies the above formal potential operator (3) corresponds to the Schr¨odinger equation (1) with multiple Dirac delta potentials. Let us absorb the strengths λi’s into the bras and kets, i.e., √λi|aii=|fiiand similarly for bras. In terms of the rescaled bras and kets, the potential operator becomes V=PN i=1 |fiihfi|. Substituting this into (2) in the coordinate representation, we obtain hx|P2 2m|ψi− N X i=1hx|fiihfi|ψi=Ehx|ψi,(5) The rescaling is introduced to formulate the bound state problem in terms of an eigenvalue problem of a symmetric matrix, as we will see. Inserting the completeness relation Rdp 2π~|pihp|= 1 in front of |ψiand |fii, we obtain the following integral equation, which is actually the Fourier transformation: Z∞ −∞ dp 2π~ei ~px ˜ ψ(p)p2 2m−E= N X i=1 pλiZ∞ −∞ dp 2π~ei ~p(x−ai)φ(ai) (6) where hx|pi=ei ~px,hp|ψi=˜ ψ(p), and φ(ai) = hfi|ψi=√λiψ(ai). Since two functions with the same Fourier transforms are equal, equation (6) implies that: ˜ ψ(p) = N X i=1 pλi e−i ~pai p2 2m−Eφ(ai).(7) It is interesting to remark that this solution depends on the unknown coordinate wave function at aiand the energy E. If we use the relation between the coordinate and momentum space wave function through the Fourier transformation ψ(x) = Z∞ −∞ dp 2π~ei ~px ˜ ψ(p),(8) and insert (7) into above for x=ai, we obtain the following consistency relation ψ(ai) = N X j=1 pλjZ∞ −∞ dp 2π~ ei ~p(ai−aj) p2 2m−Eφ(aj).(9) Multiplying both sides of (9) by √λiand separating the j=ith term, we have "1−λiZ∞ −∞ dp 2π~ 1 p2 2m−E#φ(ai)−Z∞ −∞ dp 2π~ N X j=1 j6=ipλiλj"ei ~p(ai−aj) p2 2m−E#φ(aj)=0.(10) 3 This equation can be written as a homogeneous system of linear equations in a matrix form: N X j=1 Φij(E)φ(aj) = 0 ,(11) where Φij(E) =                1−λiZ∞ −∞ dp 2π~ 1 p2 2m−Eif i=j , −pλiλjZ∞ −∞ dp 2π~ ei ~p(ai−aj) p2 2m−Eif i6=j . (12) As usual, the matrix elements are denoted by Φij(E) and the matrix itself by Φ, so that Φ = {Φij(E)}. Let us first assume that E < 0, i.e., E=−|E|, so that there is no real pole in the denominators of the integrands. Let us now consider the integral in the off-diagonal part. The function under the integral sign has simple poles located at the points (p=±ip2m|E|) in the complex p-plane. In order to calculate this integral by the residue method, we have to take into account separately the situations ai< ajand ai> aj. We note that only the pole with sign plus (minus) lies inside the contour of integration for ai> aj(aj> ai). Due to the exponential function, the integral over the semicircle vanishes as its radius goes to infinite [19]. Then, the value of the integral is obtained multiplying by 2πi the residue at that point: Z∞ −∞ dp 2π~ ei ~p(ai−aj) p2 2m−E=           m ~p2m|E|exp −p2m|E|(ai−aj)/~,if ai> aj, m ~p2m|E|exp −p2m|E|(aj−ai)/~,if ai< aj. (13) The diagonal part of the matrix Φ can be evaluated similarly, so equation (12) becomes: Φij(E) =              1−mλi ~p2m|E|if i=j , −mpλiλj ~p2m|E|exp −p2m|E||ai−aj|/~if i6=j . (14) Equation (11) has only non-trivial solutions if det Φ(E) = 0. Therefore, the bound state problem is solved once we find the solution to the transcendental equation det Φ(E) = 0. After that, we can find the bound state wave functions in the coordinate representation through (8). Suppose that the bound state energy, say EB, is the root of det Φ(E) = 0, and we find φB(aj) = pλjψB(aj) from Eq.(11) associated with EB. Then, the bound state wave function at aiis ψB(ai) = 1 √λi φB(ai).(15) 4 Taking into account the above considerations, we use (15) into the bound state wave function in momentum space (7) so as to obtain ˜ ψB(p) = N X i=1 pλi e−i ~pai p2 2m−EB φB(ai).(16) Then, the bound state wave function in the coordinate space can be found by just taking the inverse Fourier transform of the above momentum space wave function ψB(x) = N X i=1 λiφB(ai)pm/2e−√2m|EB| ~|x−ai| ~p|EB|,(17) where φB(ai) is defined by PN j=1 Φij(−|EB|)φB(aj) = 0. Suppose now that E > 0. In this case, we have to find the wave function and contour integrals of the form Z∞ −∞ dp 2π~ ei ~p(ai−aj) p2 2m−E,(18) whose poles are now located at p=±√2mE on the real axis, and there are four different choices of contours, each of which gives different result [20]. It is easy to see that the wave function becomes now the linear combination of the complex exponentials e±i√2mE/~|x−ai|.(19) Such a function cannot be square integrable unless it is identically zero . Therefore, there is no bound state for E > 0. A similar analysis can be done for E= 0, where the wave function becomes divergent over the whole real axis. Therefore, we conclude that Emust be negative for bound states. From the physical point of view, the bound state energies are expected to be less than the values of the potential at asymptotes. For this reason, the bound state energies for finitely many Dirac delta potentials are negative. For a single center located at x= 0 with coupling constant λ, the matrix Φ is just a 1 ×1 matrix, i.e., a single function: Φ(E)=1−mλ ~√2m|E|. Now, the condition det Φ(E) = 0 means that Φ(E) = 0, so that the bound state energy is EB=−mλ2 2~2for a single center [6]. After having found the bound state energy, we can find the bound state wave function. For N= 1, φB(ai) is some constant, say C. Then, the bound state wave function in momentum space becomes ˜ ψB(p) = √λ1 p2 2m−EB C . (20) The constant Ccan be determined from the normalization constant: λ|C|2(2m)2Z∞ −∞ dp 2π~ 1 (p2+m2λ2 ~2)2= 1 .(21) The integral in Eq. (21) can also be evaluated using the residue theorem. However, in this case the residues are at p=±imλ ~and of order two. Taking the integral, we find C=√mλ ~. Now we can find the 5 wave function associated with this bound state in the coordinate space by taking its Fourier transform. We perform the integration exactly as we did in (13), and obtain [6] ψB(x) = √mλ ~e−mλ ~2|x|.(22) Let us first consider the special case of two centers, namely twin attractive (λ1=λ2=λ) Dirac δ potentials located at a1= 0 and a2=a. Then, the expression det Φ(E) = 0 yields to the following transcendental equation: e−a√2m|E| ~=± ~p2m|E| mλ −1!.(23) For convenience, we define κ≡√2m|E| ~. Suppose that κ+(κ−) corresponds to the solution of Eq. (23) with the positive (negative) sign in front of the parenthesis, i.e., e−aκ+=~2κ+ mλ −1,(24) or e−aκ−= 1 −~2κ− mλ .(25) The bound state energies correspond to non-zero solutions for κ±of the above equations (24) and (25). The first transcendental equation (24) always has one real root, which implies the presence of at least one bound state. This is clear from the following considerations: the left hand side of (24) is a monotonically decreasing function, which goes to zero asymptotically, while the right hand side is a monotonically increasing function without any asymptote. However, the second transcendental equation (25) may or may not have a real positive solution. One real root of Eq. (25) is expected for κ−= 0. However, this cannot correspond to a bound state. In order to obtain a non trivial root, we must impose the condition that the slope of the right hand side of (25) must be smaller than the slope of the left hand side in absolute value.  d dκ 1−~2κ mλ κ=0 < d dκ e−κaκ=0 .(26) This means that the distance between the centers must be greater than some critical value for two bound states: a > ~2 mλ .(27) Hence, we conclude that there are at most two bound states for attractive twin Dirac delta potentials. The first one appears unconditionally so that it corresponds to the ground state. On the other hand, the second bound state appears only if ais sufficiently large ( ~2 maλ <1). This corresponds to the excited state of the system. Actually, the explicit solutions to Eq.(24) and Eq.(25) can be easily found and then the bound state energies are E+=−λ 2+1 aWaλ 2e−aλ 22 , 6 E−=−λ 2+1 aW−aλ 2e−aλ 22 .(28) where Wis the Lambert Wfunction [21], defined as the solution of the transcendental equation y ey=z, i.e., y=W[z]. The above explicit solutions given in terms of Lambert Wfunction have been known in the literature, see for instance [22] and the recent work [23], where the non-linear generalization of the problem has been discussed. 3 Bound States as a Finite Dimensional Eigenvalue Problem In order to study location and properties of bound states more systematically, we consider the equation (11) as the particular case of an eigenvalue problem for the matrix Φ: Φ(E)A(E) = ω(E)A(E),(29) where ωis any of the eigenvalues of the matrix Φ. Then, the zeros of the eigenvalues of Φ are just the bound state energies. In other words, the roots of the equation ω(E) = 0 (30) give the bound state energies. Hence, the eigenvalues of the linear differential equation Hψ(x) = Eψ(x) are obtained through a non-linear transcendental algebraic problem, ω(E) = 0. Let us consider the N= 2 case. For twin centers, located at a1= 0, a2=a, the eigenvalues can be explicitly calculated: ω1= 1 + mλ ~p2m|E|−1−e−1 ~√2m|E|a ω2= 1 + mλ ~p2m|E|−1 + e−1 ~√2m|E|a(31) As shown in above Fig. 1, there are always two eigenvalues of the matrix Φ. However, for λ= 2 with ~= 2m= 1, there are two bound states only if the distance between the centers is greater than the critical value a= 1. Otherwise there is only one bound state, which is consistent with the result given in the previous part. When the centers are sufficiently close to each other, one of the bound states seems to disappear, since the zeros of the first eigenvalue seems to move to the negative real axis. We also observe from Fig. 1 that the bound state energies come closer and closer as the distance between them increases. This is not surprising since the eigenvalues (31) converge to 1 −mλ ~√2m|E|as a→ ∞ so that zeroes of these degenerate eigenvalues lead to degenerate bound states in the limiting case. Now, we shall show why our method is much easier to investigate the bound state spectrum as we increase the number of Dirac delta potentials. The number of bound states is an important characteristic of any system. There are several ways to determine it for some regular potentials [24]. It is noteworthy that we can determine the maximum number of bound states of this system from the behavior of the eigenvalues of the matrix Φ(E) through the Feynman-Hellmann theorem [25, 26]. 7 E¤ Ω a=0.8 2 4 6 8 10 -0.2 0.2 0.4 0.6 E¤ Ω a=1 2 4 6 8 10 -0.2 0.2 0.4 0.6 E¤ Ω a=2 2 4 6 8 10 -0.4 -0.2 0.2 0.4 0.6 E¤ Ω a=4 2 4 6 8 10 -0.2 0.2 0.4 0.6 Figure 1: The flow of the eigenvalues of the matrix Φ as a function of |E|for different values of a. Here λ= 2 and ~= 2m= 1. To find the behavior of the eigenvalues as a function of E, let us first take the derivative of Φij(E) with respect to E. We may interchange this derivative and the integral sign in (14), since all the matrix elements Φij are analytic functions on the half plane <(E)<0. Hence, we obtain dΦij dE =−pλiλjZ∞ −∞ dp 2π~ ei ~p(ai−aj) (p2 2m−E)2.(32) Now, let us make use of the Feynman-Hellmann theorem, which states that dω(E) dE =hAk|dΦij dE |Aki,(33) where Akis a given normalized eigenvector for ω(E). In other words, the Feynman-Hellmann theorem states that the derivative of the eigenvalue of a parameter dependent Hermitian matrix is equal to the expectation value of the derivative of the matrix with respect to its normalized eigenvector. The FeynmanHellmann theorem can be generalized for the degenerate states [27] but it does not change our conclusion that we will draw. Thus, we have dω(E) dE =−pλiλj N X i,j=1 (Ak i)∗Z∞ −∞ dp 2π~ ei ~p(ai−aj) (p2 2m−E)2Ak j =−Z∞ −∞ dp 2π~ 1 (p2 2m−E)2 N X i=1 e−i ~paipλiAk i 2 <0.(34) For E=−|E|,dω(E) d|E|>0. Since there are at most Ndistinct eigenvalues of the N×Nmatrix Φ and these eigenvalues are monotonically increasing functions of |E|, there must be at most Nbound states. This 8 conclusion would have been rather difficult to arrive just by following the standard method, in which the properties of the bound states are just determined by matching conditions at the locations of the delta centers. 4 Degeneracies in the Bound States for Periodically Distributed Centers Let us consider NDirac delta potentials located equidistantly, i.e., a0= 0, a1=a, a2= 2a, . . . , aN= (N−1)aand λ1=. . . =λN=λ. Then, the matrix Φ given in Eq. (14) takes the following form         c0c1··· cN−2cN−1 cN−1c0c1··· cN−2 cN−2cN−1....... . . . . .. . .......c1 c1c2··· cN−1c0        N×N (35) where c0= 1 −mλ ~√2m|E|and cj=cN−j=−mλ ~p2m|E|exp −p2m|E|ja/~(36) for all j= 1, . . . , N −1. The form of the matrix above (35) is usually known as the circulant matrix. By using the Fourier matrix, it can be diagonalized and its eigenvalues can be found easily [15]. However, showing this is the beyond the scope of the main aim of this paper. Nevertheless, it is a simple exercise to show that         c0c1··· cN−2cN−1 cN−1c0c1··· cN−2 cN−2cN−1....... . . . . .. . .......c1 c1c2··· cN−1c0                1 ζl ζ2l . . . ζ(N−1)l        =λl        1 ζl ζ2l . . . ζ(N−1)l        ,(37) where ζis the Nth root of the unity, i.e., ζ=e2πi/N and j= 0,1, . . . , N −1 and λl’s are the eigenvalues of the matrix Φ, given by ωj= N−1 X k=0 ckζjk .(38) Note that the above formula is reduced to (31) for N= 2. We realize that for N≥3, the eigenvalues are degenerate ωj=ωN−j(39) 9 For the third part, let us assume that λis degenerate. Hence, we can find two real orthonormal eigenvectors (uj) and (vj) associated with λ. Suppose that ui<0 for some i. From the addition of Eq. (B-1) and (B-4), we have λ(ui+xi) = X j aij (uj+xj)⇒λ(ui+|ui|) = X j aij (uj+|uj|).(B-5) Then, uj=−|uj|for every j. If we assume that ui>0 for some iand subtracting Eq. (B-4) from (B-1), we obtain uj=|uj|. That means uj=±|uj|and by applying the same procedure, we also have vj=±|vj|. Therefore, X j vjuj=±X j|vjuj|.(B-6) Since |uj|,|vj| 6= 0 for all j,|vjuj| 6= 0 which means that uand vcannot be orthogonal. Because of the contradiction with the first assumption, we conclude that λis non-degenerate. As for the last part, let (wj) be a normalized eigenvector associated with µsuch that µ<λ, X j aij wj=µ wi.(B-7) From the variational property and the non degeneracy of λ, we have λ > X ij aij |wi||wj|≥X ij aij wiwj=|µ|.(B-8) Appendix-C: A Proof of Cauchy Interlacing Theorem Here we give a simple proof of Cauchy interlacing theorem using intermediate value theorem. This proof is originally given in [35] and we give it here in order to be self-contained. Without loss of generality, the submatrix Boccupies rows 2,3, . . . , N and columns 2,3, . . . , N. Then, the matrix Ahas the following form: A=ay† yB(C-1) where †denotes the Hermitian conjugation. Since Bis also Hermitian, we can diagonalize it by a unitary transformation U: U†B U =D , (C-2) where D= diag(µ2, µ3, . . . , µN). For simplicity, let us define a new vector z= (z2, z3, . . . , zN)T:= U†y, where Tdenotes the transposition. Here we only give the proof for the special case, where µN< µN−1< ··· , < µ3< µ2and zi6= 0 for all i= 2,3, . . . , N. The complete proof can be found in [35]. Let V=10T 0U,(C-3) 16 where 0is the zero vector. Since Uis unitary, Vis also unitary. It is easy to see that V†A V =a z† z D .(C-4) Let us define the following function f: f(x) := det(xI −A),(C-5) where Idenotes the identity matrix. Since the determinant is invariant under unitary transformations, we have f(x) = det(xI −V†A V ), or explicitly f(x) =           x−a−z2∗−z3∗··· −zN−1∗−zN∗ −z2x−µ20 0 ··· 0 −z30x−µ30.... . . . . .. . .. . .. . ....0 −zN−10··· 0x−µN−10 −zN0 0 0 0 x−µN           .(C-6) If we expand this determinant along the first row, we get f(x)=(x−a)(x−µ2)···(x−µN)− N X i=2 fi(x),(C-7) where fi(x) = |zi|2(x−µ2)··· \ (x−µi)···(x−µN) for i= 2,3, . . . , N. Here the factor with a hat is deleted. Note that fi(µj) = 0 for j6=iand fi(µi)     >0if iis even , <0if iis odd . (C-8) Since f(µi) = −fi(µi), the sign of f(µi) is opposite to that of fi(µi). It is easy to see that f(x) is a polynomial of degree Nwith positive leading coefficient. Using this and the fact f(x) is the characteristic equation for the matrix A, and intermediate value theorem [36], we conclude that there exist Nroots λ1, λ2,··· , λNof the equation f(x) = 0 such that λN< µN< λN−1< µN−1<··· < λ2< µ2< λ1.(C-9) This is the result we wanted to show. References [1] Yu. N. Demkov and V. N. Ostrovskii, Zero-range Potentials and Their Applications in Atomic Physics (Plenum Press, 1988). 17 [2] M. Belloni and R. W. Robinett, Phys. Rep. 540, 25 (2014). [3] R. de L. Kronig and W. G. Penney, Proc. R. Soc. A 130, 499 (1931). [4] C. Cohen-Tannoudji, B. Diu and F. Laloe F, Quantum Mechanics, Vol. 1 (Wiley-Interscience, 2006). [5] I. R. Lapidus, Am. J. Phys. 38, 905-908 (1970). [6] D. J. Griffiths, Introduction to Quantum mechanics (Pearson Education, 2005). [7] L. D. Landau and E. M. Lifshitz, Quantum Mechanics (Pergamon Press, 1977). [8] R. Loudon, Am. J. Phys. 27, 649 (1959). [9] W. Kwong, J. L. Rosner, J. F. Schonfeld, C. Quigg and H. B. Thacker, Am. J. Phys. 48, 926 (1980). [10] J. M. Cohen and B. Kuharetz, J. Math. Phys. 34, 12 (1993). [11] K. Bhattacharyya and R. K. Pathak, Int. J. Quantum. Chem. 59, 219 (1996). [12] S. Kar and R. R. Parwani, Eur. Phys. Lett. 80, (2007). [13] S. De Vincenzo, Brazilian Journal of Physics, 38, no. 3A, (2008). [14] A. Dutt, T. Nath, S. Kar, and R. Parwani, Eur. Phys. J. Plus 127, 28 (2012). [15] C. D. Meyer, Matrix Analysis and Applied Linear Algebra (SIAM, 2000). [16] S. Albeverio, F. Gesztesy, R. Hoegh-Krohn, and H. Holden, Solvable Models in Quantum Mechanics, 2nd ed; (AMS, 2004). [17] G. Bonneau, J. Faraut, and G. Valent, Am.J.Phys. 69, 3 (2001). [18] V. S. Araujo, F. A. B. Coutinho, and J. F. Perez, Am.J.Phys. 72, 2 (2004). [19] R. B. Ash, W. P. Novinger, Complex Analysis (Dover, 2007). [20] P. Dennery, A. Krzywicki, Mathematics for Physicists (Dover, 1996). [21] R. M. Corless, G. H. Gonnet, D. E. G. Hare, D. J. Jeffrey, and D. E. Knuth, On the Lambert W Function, Advances in Computational Mathematics, 5No. 1, pp. 329-359 (1996). [22] S. Tunalı, Point interactions in quantum mechanics, Ms. Thesis, ˙ Izmir Institute of Technology, 2014. [23] H. Kovarik, A. Sacchetti, J. Phys. A: Math. and Theo. 43, (2010). [24] E. B. Manoukian, Quantum Theory, A wide spectrum (Springer, 2006). [25] R. P. Feynman, Phys. Rev. 56, 340 (1939). [26] H. Hellmann, Einf¨uhrung in die Quantenchemie (Franz Deuticke, Leipzig, 1937), p. 285. [27] S. R. Vatsya, Phys. Rev. B 69, 037102 (2004). 18 [28] J. Dorsey, C. R. Johnson, and Z. Wei, Spec. Matrices 2, 200 (2014). [29] J. W. Brown, R. V. Churchill, Complex Variables and Applications, Eighth edition (McGraw-Hill, 2004). [30] G. Palma, U. Raff, Can. J. Phys. 84, 787-800 (2006). [31] M. Reed and B. Simon, Methods of Modern Mathematical Physics, vol IV, (Academic Press, 1978). [32] F. Erman, O. T. Turgut, J. Phys. A: Math. Theor. 43, 335204 (2010). [33] F. Ninio, J. Phys. A: Math. Gen. 9, No. 8 (1976). [34] A. Sacchetti, J. Phys. A: Math. Theor. 49, 175301 (2016). [35] S-G. Hwang, The American Mathematical Monthly 111, No. 2, pp. 157 (2004). [36] J. Steward, Calculus: Early Transcendentals, 7th edition, (Brooks/Cole, Cengage Learning, 2012). 19