1 Introduction

Dirac’s well-known equation [1, 2] characterises the relativistic dynamics of spin-\(\frac{1}{2}\) particles in the framework of quantum mechanics respecting also the principles of special relativity [3]. This equation has been very successful in its early days by providing a clear formalism for the spin of an electron as a point particle and also has led to the postulation of the existence of its anti-particle which was discovered shortly afterwards reassuring Dirac’s interpretation [4]. The Dirac equation also paved the way for quantum electrodynamics [5] and quantum field theory of electromagnetic interactions [6]. Nowadays, the Dirac equation is also an important tool for the description of the dynamics of charge carriers in carbon nanostructures like graphene [7] or more general in so-called Dirac electronic systems [8]. Hence, exact solutions of the Dirac equation are of great interest but are rare compared to its non-relativistic counterpart, the Schrödinger equation. Here various techniques for finding exact solutions like the factorisation methods have been developed during the last century and more recently received much attention in the context supersymmetric methods [9,10,11,12]. Also approximation methods like the WKB approach or perturbation methods are nowadays well-established for the non-relativistic quantum mechanics.

Supersymmetric (SUSY) quantum mechanics in the context of relativistic dynamics characterised by the Dirac equation was originally studied by Jackiw [13] and Ui [14]. See also the work by Cooper et al. [15] and the book by Thaller [3]. Due to their relevance to the current paper, let us also mention more recent studies focusing on the two-dimensional graphene in a magnetic field [16,17,18] and those related to \((1+1)\) dimensional Dirac electron interacting with a scalar or pseudo-scalar potential [19,20,21,22,23,24].

As the Dirac Hamiltonian with supersymmetry is known to be closely related to a non-relativistic Pauli–Schrödinger-type Hamiltonian [10, 13, 14], it is natural to investigate the possibility to derive exact solutions as well as approximations for a Dirac system by reducing this problem to a non-relativistic system. Whereas in the recent work [25], we focused on \((3+1)\)-dimensional systems, we will limit ourselves in the current paper to supersymmetric Dirac Hamiltonians in \((1+1)\) dimensions.

We will show in the next section that supersymmetry puts a very strict condition on the most general Dirac Hamiltonian in \((1+1)\) dimensions. In essence, we find that SUSY basically allows only for either a scalar or a pseudo-scalar potential. In both cases, the associated non-relativistic system is the well-studied Witten model of supersymmetric quantum mechanics [26]. After recalling in Sect. 3, some basic properties of supersymmetric Dirac Hamiltonians and the Witten model, we show in Sect. 4 how the eigenvalue problem of the Dirac system can be reduced to that of the Witten model. As an explicit example, the Dirac oscillator is considered and its spectral properties are derived from those of the standard harmonic oscillator on the real line. In Sect. 5, we will study the resolvent of the Dirac system and show how this can be expressed in terms of the resolvent of the non-relativistic Witten model. Again the Dirac oscillator is chosen as an explicit example and the corresponding Green’s function is derived in closed form. In Sect. 6, we will utilise the quasi-classical approximation of the Witten model to arrive at quasi-classical approximation of supersymmetric Dirac Hamiltonians. This approximation is known to respect the spectral symmetry implied by supersymmetry, and this property is also respected by the derived approximation formulas of the Dirac system. Finally in the conclusions, we give an outlook how that current approach can, for example, be applied to radial Dirac Hamiltonians and a Dirac particle in a box.

2 General Dirac Hamiltonian with supersymmetry

The most general \((1+1)\)-dimensional Dirac Hamiltonian, acting on Hilbert space \(\mathcal{H} = L^2(\mathbb {R})\otimes \mathbb {C}^2\), can be put into the form

$$\begin{aligned} H_\mathrm{D} := c p \sigma _1 + W(x) \sigma _2 + \left[ mc^2 + S(x)\right] \sigma _3 + eV(x)\mathbf{1}, \end{aligned}$$
(1)

where x and \(p=(\hbar /\mathrm{i})\partial _x\) are the position and momentum operators on \(L^2(\mathbb {R})\), respectively, \(\{\sigma _i|i=1,2,3\}\) are the Pauli matrices and \(\mathbf{1}\) is the \(2\times 2\) unit matrix acting on \(\mathbb {C}^2\). In the above, \(m>0\) and e stand for the mass and the charge of the Dirac particle moving along the real line \(\mathbb {R}\), respectively, and \(c>0\) represents the speed of light. This particle interacts with various potentials, namely a scalar potential S, a pseudo-scalar potential W (in \(3+1\) dimensional models this is also called a tensor potential) and an electro-static potential V, which is the 0-component of an electromagnetic vector potential (VU). The 1-component U of this vector potential has been omitted as it can be gauged away due the local u(1) gauge-invariance. In essence, in \((1+1)\) dimensions the presence of a vector potential would lead to a pure phase factor of the form \(\exp \{\pm \mathrm{i}(e/c)\int \mathrm{d}x\,U(x)\}\) in the upper and lower component of the 2-spinor, respectively.

Representing the Pauli matrices in the standard form

$$\begin{aligned} \sigma _1 = \left( \begin{array}{ll} 0 &{}\quad 1 \\ 1 &{}\quad 0 \end{array}\right) ,\qquad \sigma _2 = \left( \begin{array}{ll} 0 &{}\quad -\mathrm{i} \\ \mathrm{i} &{}\quad 0 \end{array}\right) ,\qquad \sigma _3 = \left( \begin{array}{ll} 1 &{}\quad 0 \\ 0 &{}\quad -1 \end{array}\right) , \end{aligned}$$
(2)

the above Dirac Hamiltonian can be put into a matrix form

$$\begin{aligned} H_\mathrm{D} = \left( \begin{array}{ll} M_+&{}\quad A \\ A^\dag &{}\quad -M_- \end{array}\right) , \end{aligned}$$
(3)

where

$$\begin{aligned} M_\pm := mc^2 + S(x) \pm eV(x)\,\qquad \mathrm{and}\qquad A:= cp-\mathrm{i} W(x). \end{aligned}$$
(4)

In order to represent a supersymmetric Dirac Hamiltonian, these operators are required to obey the condition [3, 10]

$$\begin{aligned} AM_-=M_+A. \end{aligned}$$
(5)

This condition in essence guarantees that the squared Hamiltonian \(H^2_\mathrm{D}\) becomes block diagonal and also assures the existence of an exact Foldy–Wouthuysen transformation [27, 28]. Above condition, however, turns out to be very restrictive and leads to constraints on the potentials with which the Dirac particle interacts. To be more explicit, the constraints are \(V(x) = 0\) and \(S(x)=const.\), that is, the electro-static potential must vanish and the scalar potential has to be constant. Without loss of generality, we also set \(S=0\) as a constant S can always be absorbed by the mass term \(mc^2\).

Before studying the properties of such a supersymmetric Dirac Hamiltonian, let us briefly consider here alternative representations of the general Hamiltonian (1). In fact, the representation (2) of the Pauli matrices is the one where \(\sigma _3\) is in diagonal form. Obviously, one could choose a basis in the \(\mathbb {C}^2\) subspace where one of the other two Pauli matrices becomes diagonal. This would be equivalent to a cyclic permutation of the three Pauli matrices in (1) and keeping the representation (2) fixed. Note that the algebra obeyed by the Pauli matrices is invariant under such cyclic permutations. Hence, an alternative to (1) would be the general Hamiltonian

$$\begin{aligned} \tilde{H}_\mathrm{D} := c p \sigma _2 + W(x) \sigma _3 + \left[ mc^2 + S(x)\right] \sigma _1 + eV(x)\mathbf{1}. \end{aligned}$$
(6)

Here we use a tilde to distinguish this representation for the standard one. In the matrix representation (2), this Hamiltonian reads

$$\begin{aligned} \tilde{H}_\mathrm{D} = \left( \begin{array}{ll} \tilde{M}_+&{}\quad \tilde{A} \\ \tilde{A}^\dag &{}\quad -\tilde{M}_- \end{array}\right) . \end{aligned}$$
(7)

Hence, we may identify the operators \(\tilde{M}^\pm \) and \(\tilde{A}\) as follows.

$$\begin{aligned} \tilde{M}_\pm := eV(x)\pm W(x),\qquad \tilde{A}:= -\mathrm{i} cp + mc^2 + S(x). \end{aligned}$$
(8)

and the SUSY condition (5) leads us to the restriction \(W(x)=0\) and \(V(x)=const\). As the constant V in essence is a constant shift in the energy scale, we can set it to zero without loss of generality. Note that this corresponds to the supersymmetric representation of the Dirac Hamiltonian with a scalar field as discussed by Thaller [3] in his Section 5.5.1. See also Section 9.4 in Ref. [10] and the detailed discussion by Nogami and Toyama [19]. Here we note that the case of a scalar potential can be reduced to the previous case of a pseudo-scalar potential by identifying \(W(x) = mc^2 + S(x)\) and \(M_\pm =0\). See, for example, Ref. [22]. We conclude that there is no need for a separate discussion of this case in the context of the current work.

In principle there is a third representation (we use a hat to indicate this case) of the general Hamiltonian (1) given by

$$\begin{aligned} \hat{H}_\mathrm{D} := c p \sigma _3 + W(x) \sigma _1 + \left[ mc^2 + S(x)\right] \sigma _2 + eV(x)\mathbf{1}, \end{aligned}$$
(9)

which leads to the identification

$$\begin{aligned} \hat{M}_\pm := eV(x)\pm cp\,\qquad \mathrm{and}\qquad \hat{A}:= W(x)-\mathrm{i} (mc^2 + S(x)). \end{aligned}$$
(10)

Here the SUSY condition (5) implies \(mc^2 + S(x)=0\) or \(S(x)=-mc^2=const\). In addition, we have \(W(x)=const.\) Both imply that the operator \(\hat{A}=W=const.\), leading to a trivial SUSY structure as will become obvious in the next section. Hence, we will not pursue this case either. That is, from now on we will, without loss of generality, limit our discussion to the case of a pure pseudo-scalar W, which in turn is directly applicable to the second non-trivial case of a pure scalar potential.

3 Properties of supersymmetric Dirac Hamiltonians

As explicitly discussed in the previous section, without loss of generality, the most general supersymmetric Dirac Hamiltonian in \((1+1)\) dimensions can be put into the form

$$\begin{aligned} H_\mathrm{D}=\left( \begin{array}{ll} mc^2 &{}\quad cp-\mathrm{i} W(x)\\ cp +\mathrm{i} W(x) &{}\quad -mc^2 \end{array}\right) . \end{aligned}$$
(11)

That is, it is of the supersymmetric form (3) with

$$\begin{aligned} M_\pm = mc^2\,\qquad \mathrm{and }\qquad A=cp-\mathrm{i} W(x). \end{aligned}$$
(12)

The \(N=2\) SUSY structure becomes explicit by introducing the supersymmetric Hamiltonian

$$\begin{aligned} H_\mathrm{SUSY}:=\frac{1}{2mc^2}\left( H_\mathrm{D}^2-m^2c^4\right) =\left( \begin{array}{ll} H_+ &{}\quad 0 \\ 0 &{}\quad H_- \end{array}\right) \end{aligned}$$
(13)

with partner Hamiltonians

$$\begin{aligned} H_+ := \frac{1}{2mc^2}AA^\dag ,\qquad H_- := \frac{1}{2mc^2}A^\dag A \end{aligned}$$
(14)

and the SUSY charges

$$\begin{aligned} Q:= \frac{1}{2mc^2}\left( \begin{array}{ll} 0 &{}\quad A \\ 0 &{}\quad 0 \end{array}\right) ,\qquad Q^\dag = \frac{1}{2mc^2}\left( \begin{array}{ll} 0 &{}\quad 0 \\ A^\dag &{}\quad 0 \end{array}\right) . \end{aligned}$$
(15)

These operators obey the \(N=2\) SUSY algebra

$$\begin{aligned} H_\mathrm{SUSY}=\{Q,Q^\dag \},\quad Q^2 = 0 =(Q^\dag )^2,\quad \{Q,\sigma _3\}=0, \quad [ H_\mathrm{SUSY}, \sigma _3 ] = 0,\quad \sigma _3^2 = 1, \end{aligned}$$
(16)

where \(\sigma _3\) plays the role of the grading (or Witten) operator. In fact, the partner Hamiltonians have the explicit form

$$\begin{aligned} H_\pm =\frac{p^2}{2m} +\varPhi ^2(x)\pm \frac{\hbar }{\sqrt{2m}}\varPhi '(x), \end{aligned}$$
(17)

which is the well-known one-dimensional non-relativistic Witten model of supersymmetric quantum mechanics [10, 26]. Note that in the above, we have rescaled the pseudo-scalar potential \(W(x)=: \sqrt{2mc^2}\, \varPhi (x)\), where \(\varPhi \) may now be identified with the SUSY potential of the Witten model.

This Witten model has been studied extensively in the last 25 years and finds many applications in various fields of physics. Here let us summarise the most essential properties of this model. Obviously, both partner Hamiltonians are non-negative and furthermore are essential isospectral, that means, there strictly positive eigenvalues are identical. Let us assume that both have a pure discrete spectrum represented by the eigenvalues \(\varepsilon _n\) with associated eigenstates \(\phi _n^\pm \). That is, we have

$$\begin{aligned} H_\pm \phi ^\pm _n = \varepsilon _n\phi ^\pm _n,\qquad \varepsilon _n >0,\qquad n=1,2,3,\ldots . \end{aligned}$$
(18)

In addition, in case of an unbroken SUSY, there exists a zero-energy ground state, which for convenience is then assumed to belong to \(H_-\). This can always be achieved with a change of sign in the SUSY potential. Hence, in the case of unbroken SUSY, we have in addition the ground state \(\phi ^-_0\) associated with \(\varepsilon _0=0\) defined via \(A\phi ^-_0=0\) leading to the explicit form

$$\begin{aligned} \phi ^-_0(x) = N\exp \left\{ -\frac{\sqrt{2m}}{\hbar }\int \mathrm{d} x\,\varPhi (x)\right\} = N\exp \left\{ -\frac{1}{\hbar c}\int \mathrm{d} x\,W(x)\right\} \end{aligned}$$
(19)

with N denoting a normalisation constant. The eigenstates (18) belonging to the strictly positive spectrum are related to each other via the SUSY transformations

$$\begin{aligned} A\phi ^-_n = \sqrt{2mc^2\varepsilon _n}\phi ^+_n,\qquad A^\dag \phi ^+_n = \sqrt{2mc^2\varepsilon _n}\phi ^-_n,\qquad n=1,2,3,\ldots . \end{aligned}$$
(20)

Before studying the spectral properties of \(H_\mathrm{D}\), let us briefly comment on the supersymmetric representation which covers the case of a scalar potential S. As pointed out in the previous section, the case of a pure scalar potential can be covered by above discussion when we identify

$$\begin{aligned} W(x) := mc^2 + S(x),\quad \tilde{A} := -\mathrm{i} \left( cp+\mathrm{i}W(x)\right) ,\quad H_\mathrm{SUSY} = (\tilde{H}^2_\mathrm{D}-e^2V^2)/2mc^2, \end{aligned}$$
(21)

where \(V=const.\)

It is now also obvious that the third case mentioned in Sect. 3 does lead to a trivial SUSY structure as in this case \(\hat{A}=W=const.\), which leads to trivial partner Hamiltonians \(H_\pm = W^2/2mc^2=const.\)

4 Spectral properties of the Hamiltonian

In this section, we will review the spectral properties of the most general supersymmetric Dirac Hamiltonian in \((1+1)\) dimensions, which is fully characterised by a pseudo-scalar potential as shown in the previous section, cf. Eq. (11). Such supersymmetric Dirac Hamiltonians are known to be block-diagonalisable via a unitary transformation separating positive and negative energy eigenspaces. Indeed, it is possible to show that there exists a unitary matrix U, see for example Ref. [3, 10], which transforms the Dirac Hamiltonian (11) into the form

$$\begin{aligned} H_\mathrm{FW}:= U H_\mathrm{D} U^\dag =\left( \begin{array}{ll} \sqrt{2mc^2H_++m^2c^4} &{}\quad 0 \\ 0 &{}\quad -\sqrt{2mc^2H_-+m^2c^4} \end{array}\right) . \end{aligned}$$
(22)

Hence, the positive and negative energy spectrum of \(H_\mathrm{FW}\) and hence also of \(H_\mathrm{D}\) is fully determined by that of \(H_+\) and \(H_-\), respectively. As both partner Hamiltonians are essential isospectral, the spectrum of \(H_\mathrm{D}\) is in fact symmetric about the origin,

$$\begin{aligned} E^\pm _n = \pm \sqrt{2mc^2\varepsilon _n + m^2c^4},\qquad n>0. \end{aligned}$$
(23)

In addition, in the case of an unbroken SUSY the eigenvalue \(E^-_0=-mc^2\) also belongs to the spectrum of \(H_\mathrm{D}\). The corresponding eigenstates are given by

$$\begin{aligned} \psi _n^\pm = U^\dag \tilde{\psi }^\pm _n,\qquad \mathrm{with}\qquad \tilde{\psi }^+_n=\left( \begin{array}{c}\phi _n^+ \\ 0 \end{array} \right) ,\qquad \tilde{\psi }^-_n=\left( \begin{array}{c}0 \\ \phi _n^- \end{array} \right) . \end{aligned}$$
(24)

In other words, the spectral properties of \(H_\mathrm{D}\), that is,

$$\begin{aligned} H_\mathrm{D} \psi ^\pm _n = E^\pm _n\psi ^\pm _n, \end{aligned}$$
(25)

are fully characterised by those of \(H_\mathrm{SUSY}\). In case of unbroken SUSY in addition we have

$$\begin{aligned} H_\mathrm{D} \psi ^-_0 = -mc^2\psi ^-_0,\qquad \psi ^-_0 = U^\dag \left( \begin{array}{c}0 \\ \phi _0^- \end{array} \right) = \left( \begin{array}{c}0 \\ \phi _0^- \end{array} \right) . \end{aligned}$$
(26)

Note that \(U=1\) on \(\ker H_-\) as will be shown below.

Let us now study the unitary transformation matrix U in more detail. According to the general SUSY approach [3, 10], this matrix is given by

$$\begin{aligned} U := a_+ + \sigma _3 \,\mathrm{sgn}\, Q_1 \, a_-\, \qquad \mathrm{with} \qquad a_\pm := \sqrt{\frac{|H_\mathrm{D}|\pm mc^2}{2 |H_\mathrm{D}|}}. \end{aligned}$$
(27)

Here the self-adjoint supercharge \(Q_1\) is defined as follows

$$\begin{aligned} Q_1 := \sqrt{2mc^2}\left( Q+Q^\dag \right) =\left( \begin{array}{ll} 0 &{}\quad A \\ A^\dag &{}\quad 0 \end{array}\right) \end{aligned}$$
(28)

and the definition of the sign function is

$$\begin{aligned} \mathrm{sgn}\,x := \left\{ \begin{array}{c} +1 \\ 0 \\ -1 \end{array}\right. \qquad \mathrm{for} \qquad \begin{array}{c} x > 0 \\ x = 0 \\ x <0 \end{array}. \end{aligned}$$
(29)

To make things a bit more explicit, let us first look into the spectral properties of \(Q_1\). It is straight forward to show with the help of the SUSY transformations (20) that the states

$$\begin{aligned} \chi _n^\pm := \frac{1}{\sqrt{2}}\left( \begin{array}{c}\phi ^+_n \\ \pm \phi ^-_n \end{array}\right) \end{aligned}$$
(30)

are eigenstates of \(Q_1\), that is,

$$\begin{aligned} Q_1\chi _n^\pm = \pm \sqrt{2mc^2\varepsilon _n}\chi _n^\pm . \end{aligned}$$
(31)

In addition, for unbroken SUSY, \(Q_1\) has a zero eigenvalue as

$$\begin{aligned} Q_1\chi _0^- = 0\qquad \mathrm{for}\qquad \chi _0^- := \left( \begin{array}{c}0\\ \phi ^-_0 \end{array}\right) , \end{aligned}$$
(32)

which implies \(U= \mathbf{1}\) on \(\mathrm{ker}\, H_-=\mathrm{span}\,|\phi _0^-\rangle \langle \phi _0^-|\). Having this in mind, it is obvious that the operator \(\mathrm{sgn}\,Q_1\) can explicitly be written as

$$\begin{aligned} \mathrm{sgn}\,Q_1 = \left( \begin{array}{ll} 0 &{}\quad \left( AA^\dag \right) ^{-1/2}A \\ \left( A^\dag A\right) ^{-1/2}A^\dag &{}\quad 0 \end{array}\right) . \end{aligned}$$
(33)

Note that \(\mathrm{sgn}\,Q_1\chi _n^\pm =\pm \chi _n^\pm \).

Together with the diagonal matrix

$$\begin{aligned} |H_\mathrm{D}|=\left( \begin{array}{ll} \sqrt{AA^\dag +m^2c^4} &{}\quad 0 \\ 0 &{} \quad \sqrt{A^\dag A+m^2c^4} \end{array}\right) \end{aligned}$$
(34)

and the dimensionless operator \(a:=A/mc^2\), we have

$$\begin{aligned} a_\pm =\left( \begin{array}{ll} \displaystyle \sqrt{\frac{\sqrt{aa^\dag +1} \pm 1}{2\sqrt{aa^\dag +1}}} &{}\quad 0 \\ 0 &{}\quad \displaystyle \sqrt{\frac{\sqrt{a^\dag a + 1} \pm 1}{2\sqrt{a^\dag a + 1}}} \end{array}\right) \end{aligned}$$
(35)

leading to the explicit but rather complicated expression

$$\begin{aligned} U=\frac{1}{\sqrt{2}}\left( \begin{array}{ll} \displaystyle \sqrt{\frac{\sqrt{aa^\dag +1} + 1}{\sqrt{aa^\dag +1}}} &{}\quad \displaystyle (aa^\dag )^{-1/2}a \sqrt{\frac{\sqrt{a^\dag a + 1} - 1}{\sqrt{a^\dag a + 1}}}\\ \displaystyle -(a^\dag a)^{-1/2}a^\dag \sqrt{\frac{\sqrt{aa^\dag +1} - 1}{\sqrt{aa^\dag +1}}} &{}\quad \displaystyle \sqrt{\frac{\sqrt{a^\dag a + 1} + 1}{\sqrt{a^\dag a + 1}}} \end{array} \right) . \end{aligned}$$
(36)

Noting that

$$\begin{aligned} a\phi _n^- = \sqrt{\frac{2\varepsilon _n}{mc^2}} \phi _n^+,\qquad a^\dag \phi _n^+ = \sqrt{\frac{2\varepsilon _n}{mc^2}} \phi _n^-, \end{aligned}$$
(37)

which implies

$$\begin{aligned} a^\dag (aa^\dag )^{-1/2}\phi _n^+=\phi _n^-,\qquad a(a^\dag a)^{-1/2}\phi _n^-=\phi _n^+, \end{aligned}$$
(38)

we obtain explicit expressions for the eigenfunctions (24) of the original Dirac Hamiltonian (11) in terms of the eigenfunctions of the corresponding non-relativistic Witten model

$$\begin{aligned} \begin{array}{l} \psi _n^+ = \displaystyle \frac{1}{\sqrt{2}}\left( \begin{array}{c} \sqrt{1+\left( 1+2\varepsilon _n/mc^2\right) ^{-1/2}}\phi _n^+ \\ \sqrt{1-\left( 1+2\varepsilon _n/mc^2\right) ^{-1/2}}\phi _n^- \end{array} \right) = \frac{1}{\sqrt{2}}\left( \begin{array}{c} \sqrt{1+\frac{mc^2}{E_n^+}}\phi _n^+ \\ \sqrt{1+\frac{mc^2}{E_n^-}}\phi _n^- \end{array} \right) ,\\ \psi _n^- = \displaystyle \frac{1}{\sqrt{2}}\left( \begin{array}{c} - \sqrt{1-\left( 1+2\varepsilon _n/mc^2\right) ^{-1/2}}\phi _n^+ \\ \sqrt{1+\left( 1+2\varepsilon _n/mc^2\right) ^{-1/2}}\phi _n^- \end{array} \right) = \displaystyle \frac{1}{\sqrt{2}}\left( \begin{array}{c} - \sqrt{1-\frac{mc^2}{E_n^+}}\phi _n^+ \\ \sqrt{1-\frac{mc^2}{E_n^-}}\phi _n^- \end{array} \right) ,\\ \psi _0^- = \left( \begin{array}{c} 0 \\ \phi _0^- \end{array} \right) . \end{array} \end{aligned}$$
(39)

Note that these states constitute a complete orthonormal set on \(\mathcal{H}\) due to the orthonormality of \(\langle \phi _m^+|\phi _n^+\rangle = \langle \phi _m^-|\phi _n^-\rangle = \delta _{mn}\) and the relation \(E^-_n=-E^+_n\). Hence, we have reduced the complete eigenvalue problem of a general supersymmetric Dirac Hamiltonian (11) to that of the associated non-relativistic Witten model (17). Above results cover several special cases discusses before. For example, they agree with previous work exclusively dedicated to the Dirac oscillator [21, 29]. The special case of a vanishing mass characterising a Dirac electron in graphene under the influence of a magnetic field has in essence been discussed by Kuru, Negro and Nieto [16] and agrees with above result when taking the limit \(m\rightarrow 0\).

4.1 The Dirac oscillator in (1 \(+\) 1) dimensions

As an illustrative example, let us consider the \((1+1)\)-dimensional version of the Dirac oscillator [30] which is characterised by the pseudo-scalar potential

$$\begin{aligned} W(x) = mc\omega x, \qquad \omega > 0. \end{aligned}$$
(40)

Obviously the two partner Hamiltonians are represented by a shifted harmonic oscillator

$$\begin{aligned} H_\pm = \frac{p^2}{2m}+\frac{m}{2}\omega ^2x^2\pm \frac{\hbar \omega }{2}. \end{aligned}$$
(41)

The corresponding spectral properties of these Hamiltonians are well known and given by

$$\begin{aligned} \varepsilon _n = \hbar \omega n,\qquad \phi _n^+ = |n-1\rangle , \qquad \phi _n^-=|n\rangle . \end{aligned}$$
(42)

Here \(n\in \mathbb {N}\) for \(H_+\), whereas \(n\in \mathbb {N}_0\) for \(H_-\), and \(|n\rangle \) denotes a standard harmonic oscillator eigenstate which obeys the relations

$$\begin{aligned} b|n\rangle = \sqrt{n}|n-1\rangle ,\qquad b^\dag |n\rangle = \sqrt{n+1}|n+1\rangle , \end{aligned}$$
(43)

with \(b:=\mathrm{i} p/m\omega +x\) and \(b^\dag = -\mathrm{i} p /m\omega +x\) being the standard annihilation and creation operators of the harmonic oscillator. They are related to the dimensionless operators a and \(a^\dag \) introduced above as follows:

$$\begin{aligned} a=-\mathrm{i}\sqrt{\frac{2\hbar \omega }{mc^2}}\,b,\qquad a^\dag = \mathrm{i}\sqrt{\frac{2\hbar \omega }{mc^2}}\,b^\dag . \end{aligned}$$
(44)

Hence via (23) and (39) we immediately find that the eigenvalues and eigenstates of the Dirac oscillator Hamiltonian

$$\begin{aligned} H_\mathrm{D}=\left( \begin{array}{ll} mc^2 &{}\quad cp-\mathrm{i} mc\omega x \\ cp+ \mathrm{i} mc\omega x &{}\quad -mc^2 \end{array} \right) \end{aligned}$$
(45)

are given by

$$\begin{aligned} \begin{array}{l} E_n^+= mc^2\sqrt{1+2n\frac{\hbar \omega }{mc^2}},\qquad \psi _n^+={\displaystyle \frac{1}{\sqrt{2}}}\left( \begin{array}{c} \sqrt{1+\left( 1+2n\frac{\hbar \omega }{mc^2}\right) ^{-1/2}}\,|n-1\rangle \\ \sqrt{1-\left( 1+2n\frac{\hbar \omega }{mc^2}\right) ^{-1/2}}\,|n\rangle \end{array} \right) , \\ E_n^-= -mc^2\sqrt{1+2n\frac{\hbar \omega }{mc^2}},\qquad \psi _n^-={\displaystyle \frac{1}{\sqrt{2}}}\left( \begin{array}{c} -\sqrt{1-\left( 1+2n\frac{\hbar \omega }{mc^2}\right) ^{-1/2}}\,|n-1\rangle \\ \sqrt{1+\left( 1+2n\frac{\hbar \omega }{mc^2}\right) ^{-1/2}}\,|n\rangle \end{array} \right) , \end{array}\end{aligned}$$
(46)

where \(n=0\) is only allowed for the second line with \(\psi ^-_0\) having only a lower component, cf. last line in Eq. (39). As mentioned before, these results are not new and agree with previous work dedicated to the Dirac oscillator in \((1+1)\) dimensions [21, 29]. Let us note that we can also recover from this the zero-mass case discussed by Kuru, Negro and Nieto [16] when setting in the above result \(2m\omega \rightarrow \omega \) and taking the limit \(m\rightarrow 0\).

5 Spectral properties of the resolvent

In this section, we will study the resolvent of a general supersymmetric Dirac Hamiltonian (11) defined by

$$\begin{aligned} G_\mathrm{D}(z):=\frac{1}{H_\mathrm{D}-z},\qquad z\in \mathbb {C}\backslash \mathrm{spec}\, H_\mathrm{D} \end{aligned}$$
(47)

which can be expressed in terms of the so-called iterated resolvent g as follows

$$\begin{aligned} G_\mathrm{D}(z) = \left( H_\mathrm{D} +z\right) g(z^2), \qquad g(z^2):= \frac{1}{H^2_\mathrm{D}-z^2}. \end{aligned}$$
(48)

As \(H^2_\mathrm{D}\) is block-diagonal,

$$\begin{aligned} H^2_\mathrm{D} =2mc^2\left( \begin{array}{ll} H_+ &{}\quad 0 \\ 0 &{}\quad H_- \end{array}\right) +m^2c^4, \end{aligned}$$
(49)

so is

$$\begin{aligned} g(z^2)=\frac{1}{2mc^2}\left( \begin{array}{ll} G_+(\zeta (z)) &{}\quad 0 \\ 0 &{}\quad G_-(\zeta (z)) \end{array}\right) , \end{aligned}$$
(50)

where

$$\begin{aligned} G_\pm (\zeta ) := \frac{1}{H_\pm -\zeta } \end{aligned}$$
(51)

denotes the resolvent for the SUSY partner Hamiltonians (17) and we have introduced

$$\begin{aligned} \zeta (z):= \frac{z^2}{2mc^2}-\frac{mc^2}{2}. \end{aligned}$$
(52)

Hence we can express the resolvent (47) in terms of the non-relativistic resolvents (51) as follows

$$\begin{aligned} G_\mathrm{D}(z)=\frac{1}{2mc^2}\left( \begin{array}{ll} (z+mc^2)G_+(\zeta (z)) &{}\quad AG_-(\zeta (z))\\ A^\dag G_+(\zeta (z))&{}\quad (z - mc^2)G_-(\zeta (z))\end{array}\right) . \end{aligned}$$
(53)

Utilising the spectral representations

$$\begin{aligned} G_\pm (\zeta )=\sum _{\varepsilon _n\ge 0}\frac{|\phi _n^\pm \rangle \langle \phi _n^\pm |}{\varepsilon _n -\zeta }, \end{aligned}$$
(54)

where the \(\varepsilon _0=0\) term is only present in \(G_-\) in case of an unbroken SUSY, we arrive with the help of the SUSY transformations (20) at the spectral representation of (47)

$$\begin{aligned} \begin{array}{l} G_\mathrm{D}(z)= \left( \begin{array}{ll} 0 &{}\quad 0\\ 0 &{}\quad \displaystyle \frac{|\phi ^-_0\rangle \langle \phi ^-_0|}{z+mc^2}\end{array}\right) \\ \quad \displaystyle + \sum _{\varepsilon _n>0}\frac{1}{2mc^2\varepsilon _n +m^2c^4-z^2} \left( \begin{array}{ll} \left( z+mc^2\right) ~|\phi ^+_n\rangle \langle \phi ^+_n|~~~&{}\quad \left( 2mc^2\varepsilon _n\right) ^{1/2} |\phi ^+_n\rangle \langle \phi ^-_n|\\ \left( 2mc^2\varepsilon _n\right) ^{1/2} |\phi ^-_n\rangle \langle \phi ^+_n| ~~~ &{}\quad \left( z-mc^2\right) ~ |\phi ^-_n\rangle \langle \phi ^-_n|\end{array}\right) . \end{array} \end{aligned}$$
(55)

Again the first term, which has a pole at \(E_0=-mc^2\), is only present in case of an unbroken SUSY. The poles of the second term reflect the energy eigenvalues as given in (23) as expected.

Note that besides the spectral representation, the diagonal form (50) also allows to derive a path integral representation of the Green’s function following, in essence, the same procedure as presented in Ref. [25]. Below we will derive the Green’s function of the Dirac oscillator in closed form directly from that of the standard harmonic oscillator.

5.1 The Green’s function of the Dirac oscillator

Let us reconsider the Dirac oscillator in \((1+1)\) dimensions. As we have seen in the previous section, the partner Hamiltonians of the Dirac oscillator are related to the standard harmonic oscillator Hamiltonian

$$\begin{aligned} H_0 := \frac{p^2}{2m}+\frac{m}{2}\omega ^2x^2 \end{aligned}$$
(56)

via constants shifts given by \(H_\pm = H_0 \pm \hbar \omega /2\). Hence, the two partner resolvents (54) can be obtained from the usual harmonic oscillator resolvent

$$\begin{aligned} G_0(\zeta ):=\frac{1}{H_0-\zeta } \end{aligned}$$
(57)

via the relation \(G_\pm (\zeta )=G_0(\zeta \mp \hbar \omega /2)\). Following a recent work by Glasser and Nieto [31], the coordinate representation of (57) can be given in closed form

$$\begin{aligned} \begin{array}{rl} G_0(x'',x';\hbar \omega \epsilon )&:=\langle x''|G_0(\hbar \omega \epsilon )|x'\rangle \displaystyle =\sqrt{\frac{m}{\pi \omega \hbar ^3}}\Gamma \left( \textstyle \frac{1}{2}-\epsilon \right) D_{\epsilon -1/2}(\mu x_+)D_{\epsilon -1/2}(-\mu x_-). \end{array} \end{aligned}$$
(58)

In the above, \(D_\nu \) stands for the parabolic cylinder function, \(\Gamma \) is Euler’s gamma function, \(\mu :=\sqrt{2m\omega /\hbar }\) and \(x_\pm \) stands for the maximum and minimum of \(x''\) and \(x'\), respectively, that is, \(x_+:=\max (x'',x')\) and \(x_-:=\min (x'',x')\). For simplicity, we have introduced the dimensionless parameter \(\epsilon := \zeta /\hbar \omega =\frac{z^2}{2mc^2\hbar \omega }-\frac{mc^2}{2\hbar \omega }\). With this explicit expression, it is now also possible to find the corresponding closed-form expression for the Green’s function of the Dirac oscillator. In fact, with the help of the relations

$$\begin{aligned} D'_\nu (y)+(y/2)D_\nu (y)= \nu D_{\nu -1}(y),\quad -D'_\nu (y)+(y/2)D_\nu (y)= D_{\nu +1}(y) \end{aligned}$$
(59)

one finds following useful relations for the creation and annihilation operators (43),

$$\begin{aligned} b D_\nu (\mu x)=\nu D_{\nu -1}(\mu x),\quad b^\dag D_\nu (\mu x) = D_{\nu +1} (\mu x), \quad b^\dag b D_\nu (\mu x) = \nu D_{\nu } (\mu x). \end{aligned}$$
(60)

From these relations, it is straight-forward to obtain the closed-form expressions

$$\begin{aligned} \begin{array}{l} \displaystyle \frac{1}{2mc^2}\langle x''|AG_-(\zeta )|x'\rangle = \frac{\mathrm{i}}{\hbar c}\frac{{\mathrm{sgn}} (x''-x')}{\sqrt{2\pi }} \Gamma (1-\epsilon )D_{\epsilon -1}(\mu x_+)D_\epsilon (-\mu x_-)\\ \displaystyle \frac{1}{2mc^2}\langle x''|A^\dag G_+(\zeta )|x'\rangle = \frac{\mathrm{i}}{\hbar c}\frac{{\mathrm{sgn}} (x''-x')}{\sqrt{2\pi }} \Gamma (1-\epsilon )D_{\epsilon }(\mu x_+)D_{\epsilon -1}(-\mu x_-)\\ \displaystyle \frac{z+mc^2}{2mc^2}\langle x''|G_+(\zeta )|x'\rangle =\frac{1}{\hbar c}\frac{z+mc^2}{\sqrt{2\pi mc^2{\hbar {\omega }}}} \Gamma (1-\epsilon )D_{\epsilon }(\mu x_+)D_{\epsilon }(-\mu x_-)\\ \displaystyle \frac{z-mc^2}{2mc^2}\langle x''|G_-(\zeta )|x'\rangle =\frac{1}{\hbar c}\frac{z-mc^2}{\sqrt{2\pi mc^2{\hbar {\omega }}}} \Gamma (-\epsilon )D_{\epsilon -1}(\mu x_+)D_{\epsilon -1}(-\mu x_-) \end{array} \end{aligned}$$
(61)

which constitute the four components of the resolvent (53) in the coordinate representation, that is \(\langle x''|G_\mathrm{D}(z)|x'\rangle \), for the Dirac oscillator (45). Obviously the pole at \(\epsilon = 0\) leads to the eigenvalue \(z=E^-_0=-mc^2\) whereas the other poles at \(\epsilon = n \in \mathbb {N}\) result in the eigenvalues as given in (46).

6 Quasi-classical approximation

In the previous section, we have seen that the supersymmetric Dirac problem in \((1+1)\) dimensions can be reduced to the corresponding non-relativistic Witten model. That is, the spectral properties of \(H_\mathrm{D}\) are fully deduced from those of the partner Hamiltonians \(H_\pm \). In particular, all exactly solvable Witten models, e.g. those which are shape-invariant, immediately result in the exact solutions of the corresponding Dirac system. Even closed-form expression of the resolvent \(G_\mathrm{D}\) can be obtained when the corresponding resolvents \(G_\pm \) are given in closed form.

However, the above general discussion of Sect. 3 is not limited to only exactly solvable Witten models. It may also be applied to approximation methods. Here we consider the so-called quasi-classical approximations to the energy eigenvalues of the Witten model. In fact, this quasi-classical approximation is in general better than the usual WKB approximation as it respects the SUSY-induced symmetry between the eigenvalues of the two partner Hamiltonians and is exact for the ground-state energy in case of unbroken SUSY [10]. To be more explicit, in case of unbroken SUSY with the ground state belonging to \(H_-\) the eigenvalues \(\varepsilon _n\) for the SUSY partners \(H_\pm \) are obtained via the so-called CBC formulaFootnote 1

$$\begin{aligned} \int _{x_\mathrm{L}}^{x_\mathrm{R}}\mathrm{d} x \sqrt{2m\left( \varepsilon -\varPhi ^2(x)\right) } = \hbar \pi n,\qquad n\in \mathbb {N}_0, \end{aligned}$$
(62)

where the left and right turning points \(x_\mathrm{L} \le x_\mathrm{R}\) are defined by \(\varPhi ^2(x_\mathrm{L/R})=\varepsilon \). Note that as before \(n=0\) is only allowed for \(H_-\) and obviously implies \(\varepsilon _0 = 0\), which is exact. For \(n\in \mathbb {N}\) above CBC formula provides a quasi-classical approximation to the jointed eigenvalues of \(H_\pm \). In case of shape-invariant systems, this approximation is even exact for any n.

In the case of a broken SUSY, the corresponding quasi-classical approximation is given by the so-called EIJ formulaFootnote 2 and reads

$$\begin{aligned} \int _{x_\mathrm{L}}^{x_\mathrm{R}}\mathrm{d} x \sqrt{2m\left( \varepsilon -\varPhi ^2(x)\right) } = \hbar \pi \left( n-\frac{1}{2}\right) ,\qquad n\in \mathbb {N}. \end{aligned}$$
(63)

Again this formula reflects the strict iso-spectral property of the partner Hamiltonians as well as the strict positivity \(\varepsilon _n>0\) for all n. Furthermore, this approximation results in exact eigenvalues in case of shape-invariant systems.

Let us now utilise the relation (23) between the eigenvalues of the Dirac Hamiltonian and the Witten Hamiltonians. For unbroken SUSY, the relativistic version of the quasi-classical approximation then reads

$$\begin{aligned} \int _{x_\mathrm{L}}^{x_\mathrm{R}}\mathrm{d} x \sqrt{E^2 - m^2c^4 - W^2(x)} = c\hbar \pi n,\qquad n\in \mathbb {N}_0, \end{aligned}$$
(64)

where the turning points \(x_\mathrm{L} \le x_\mathrm{R}\) are now given by the relation \(W^2(x_\mathrm{L/R})=E^2-m^2c^4\). Again for \(n=0\) this results in the exact ground-state energy \(E_0^-=-mc^2\). The approximate eigenvalues are given by \(E_n^\pm =\pm \sqrt{E^2}\) for \(n\in \mathbb {N}\). Similar, for broken SUSY, we have

$$\begin{aligned} \int _{x_\mathrm{L}}^{x_\mathrm{R}}\mathrm{d} x \sqrt{E^2 - m^2c^4 - W^2(x)} = c\hbar \pi \left( n-\frac{1}{2}\right) ,\qquad E_n^\pm =\pm \sqrt{E^2},\qquad n\in \mathbb {N}. \end{aligned}$$
(65)

Again, as in the non-relativistic case, whenever the SUSY potential W is shape-invariant both formulas (64) and (65) reproduce the exact spectrum of the corresponding Dirac Hamiltonian. We leave it as an exercise to the reader to verify that for \(W(x)=mc\omega x\) the formula (64) reproduces the spectrum of the Dirac oscillator as given in (46).

7 Conclusions

In the current paper, we have studied the most general supersymmetric Dirac Hamiltonian in \((1+1)\) space-time dimensions. It has been shown that the spectral properties of the relativistic Hamiltonian \(H_\mathrm{D}\) as given in Eq. (11) can be reduced to the non-relativistic spectral problem (18) of the Witten model. Hence the main results are represented by Eqs. (23) and (39). In addition, we have shown that the relativistic resolvent kernel (47) can be expressed in terms of those of the Witten model (51) via Eq. (53). As an explicit example, we have chosen the Dirac oscillator for which the closed-form expression (61) of the Green’s function was obtained for the first time. Obviously, any of the shape-invariant non-relativistic SUSY models, which exhibit exact solutions, immediately result in exact solutions of the corresponding relativistic model.

With the discussion of Sect. 6, we may also consider not exactly solvable systems by applying the supersymmetric quasi-classical approximations. Here for example, one could study an harmonic oscillator systems being characterised by a SUSY potential \(\varPhi (x)= |x|^d\) and \(\varPhi (x)=\mathrm{sgn}\, x |x|^d\) exhibiting broken and unbroken supersymmetry, respectively. It shall be noted that in the limit \(d\rightarrow \infty \) theses systems simulate a particle in the box with various boundary conditions and the quasi-classical approximation is known [10] to become exact in that limit, too. Hence, this will provide another route towards the study of a Dirac particle in a box, a topic being still of interest [35,36,37,38,39].

Despite the fact that the current paper focuses on the \((1+1)\) dimensional Dirac systems, some of the present results are valid in the more general case of arbitrary supersymmetric Dirac Hamiltonians. For example, the results of Sect. 4 and in particular the explicit expression (36) for the unitary transformation matrix U are valid for an arbitrary SUSY Hamiltonian (3) as long as \(M_+=M_-=mc^2\), which is the case for almost all supersymmetric Dirac Hamiltonians [10].

Finally, let us mention that the present results also apply to radial symmetric Dirac Hamiltonians being supersymmetric. In fact, due the spherical symmetry the radial Dirac Hamiltonian is of the same form as given in (11), cf. eq. (9.102) in [10], with a pseudo-scalar potential being of the form \(W(r)=\kappa /r+\varPhi (r)\) now acting on the positive half-line \(r\in \mathbb {R}^+\) and \(\kappa \) denotes the eigenvalues of the spin-orbit operator.