Soon after Dirac published his relativistic wave equation  in 1928, it was tried out with a variety of potentials, from which a gradual awareness of some of its strange properties arose. From the beginning, of course, it was evident that solutions with arbitrarily negative energies could exist, but as time went on, it was found that they could not always be cleanly separated from the positive energy solutions.
Even so, the harmonic oscillator never seems to have been a specific object of study. It was one of the cases included in a systematic review of power law potentials, mostly negative powers, which might induce singularities at the origin of the second order differential equation which one wanted to quantize.
Nikolsky  examined it in 1930, Postepska  five years later, after which it seems to have laid dormant until Sewell's perturbation approach  in 1949. The latter immediately attracted Titchmarsh's attention  because the spectrum was actually a continuum with sharp resonances amenable to his general theory of eigenfunction expansions based on the m-function of Weyl's thesis .
In recent years another differential equation has come to be called the ``Dirac harmonic oscillator,'' with applications to quark physics. What has happened is that Dirac's version of a relativistic wave equation does not exhaust all the possibilities for simple multicomponent wave equations.
In one dimension, there are four linearly independent matrices which are usefully referred to a quaternion basis. To keep real matrices, a somewhat Lorentzian version of the quaternions is required wherein the squares of two of them are while only the third is ; all of Hamilton's quaternions square to .
Whatever the preferred symbolism, both the Schrödinger equation and the Dirac equation have coefficient matrices which are linear combinations of and . There is no reason to omit from the construction of the coefficient matrix; indeed such a term would arise from separating Dirac's equation had an angular momentum term been present.
From the point of view of symmetry and accidental degeneracy, the new variant not only has a discrete spectrum, but an interesting symmetry group. It also reflects properties ascribed to the quark model of nuclear matter, so it is hardly surprising that it has become an object of study in its own right. Nevertheless, the present discussion is primarily concerned with resonance and continuum wave functions, so it is confined to the ``classical'' Dirac Harmonic Oscillator.
Most of the figures included in this report were obtained with the CalComp plotter attached to the PDP-10 computer at the Instituto Naciónal de Energía Nuclear in Salazar, using a program collection called SERO together with the graphics package <PLOT>. In the interim these programs have been adapted to a series of microcomputers, and to the types of visual display which they employ; some of them have also been used.