PHYSICAL REVIEW A 87, 054501 (2013)

Instability of tripositronium

Sergiy Bubin and Oleg V. Prezhdo
Department of Chemistry, University of Rochester, Rochester, New York 14627, USA

Kálmán Varga
Department of Physics and Astronomy, Vanderbilt University, Nashville, Tennessee 37235, USA

(Received 29 March 2013; published 13 May 2013)

The stability of tripositronium, a system consisting of three electrons and three positrons, has been investigated
systematically by varying the repulsion strength between like-charged particles. The possibility of the existence
of a Ps3 bound state that is stable against dissociation appears utterly unlikely based on the results of variational
calculations employing all-particle explicitly correlated Gaussian basis functions.

DOI: 10.1103/PhysRevA.87.054501 PACS number(s): 31.15.ac, 31.15.xt, 36.10.Dr

I. INTRODUCTION

Since the early work of Wheeler [1] on quantum mechanical
entities composed of electrons and positrons (which he called
polyelectrons) there has been a significant interest in finding
particle-antiparticle systems stable against dissociation. This
interest has taken place on both theoretical and experimental
sides. While Mohorovičić [2], Ruark [3], and Wheeler [1]
envisaged the existence of the positronium atom (Ps), Deutsch
[4] provided experimental evidence for it in 1951. In his work
Wheeler also predicted the stability of a larger, three-particle
system e+e−e− (Ps− or positronium negative ion). The first
experimental observation of Ps− was made three decades later
by Mills [5], who also measured its annihilation decay rate [6].
Soon after the work of Wheeler, Hylleraas and Ore [7] showed
rigorously, using the variational method, that even a larger
four-particle Ps2, which is a matter-antimatter analog of the
hydrogen molecule where positrons play the role of “nuclei,”
is also stable. The experimental confirmation of its existence
was made just recently when Cassidy and Mills performed
experiments with intense positron bursts implanted into a thin
film of porous silica [8,9]. This discovery has created a new
chapter in the experimental positronium physics and chemistry
and invigorated new theoretical studies. Further advances
in positron technology may permit much higher density of
positrons and raise questions of whether larger polyelectrons
can exist in nature [10]. Eventually, one might reach the point
where a macroscopic number of Ps atoms can undergo a phase
transition to form a Bose-Einstein condensate. This quantum
matter-antimatter system would be of great interest for many
reasons, in particular for possible production of an annihilation
γ -ray laser [8,11]. Because of the very short wavelength such a
laser could be used to probe objects as small as atomic nuclei.

Theoretical discovery and investigation of new matter-
antimatter systems remains a challenging task [12,13]. This
particularly concerns systems composed of a few electrons
and positrons. There are several reasons for this. First, these
systems (if they exist) are expected to be weakly bound. This
is illustrated by the fact that the binding energies of both
Ps− and Ps2 are only 5% and 3% of their total energies.
Second, there is no separation into heavy “nuclei” and light
electrons; thus all particles must be treated equally. Finally,
traditional quantum chemical ansatzes based on expansions in
terms of single-particle atomic orbitals are known to converge

extremely slowly for systems containing positrons due to
strong attraction between oppositely charged particles and
formation of clusters [14–17]. In order to be able to describe
such systems accurately it is necessary to approximate the
wave function with an expansion that depends explicitly on all
interparticle coordinates.

In this work we have performed a systematic study of
the stability of Ps3 (e+e+e+e−e−e− or positronium trimer).
In the framework of the variational method complemented
with the use of explicitly correlated Gaussian (ECG) basis
functions we demonstrated that this system is unlikely to
have bound states. While there has been previous work [18]
where ECG-type calculations of Ps3 were attempted, that work
did not provide a systematic study of the stability. Due to
the weak binding and complex structure of Ps3 (if it exists)
the convergence of the trial wave function in the stochastic
variational method (SVM) might have lead to a wrong
configuration (a dissociated system) corresponding to a local
minimum of the total energy. Capturing the right configuration
using a rough initial guess for the nonlinear parameters of
ECGs is generally not guaranteed when the system is very close
to dissociation. In this work we remedied this issue by slowly
varying the Coulomb repulsion strength and maintaining the
system in a bound state. As our calculations have shown, in the
limit of the exact interaction strength Ps3 becomes dissociated.

II. METHOD

It has been demonstrated in a number of applications
[18–21] that the variational method in conjunction with the
use of ECG basis functions provides a powerful framework
for studying bound states of Coulomb few-body systems with
very high accuracy. In fact, many reliable predictions of the
stability of quantum systems containing positrons were first
made using ECG expansions [22–25].

In the present work we are concerned with computing
the lowest S-state of a six-particle system with unit charges.
For definiteness, let us assume that particles 1, 2, and 3 are
positrons while particles 4, 5, and 6 are electrons. The ECG
basis functions suitable for calculations of zero total orbital
angular momentum states have the following form:

φk = exp
( − αk

12R
2
12 − αk

13R
2
13 − · · · − αk

N−1,NR2
N−1,N

)
, (1)

054501-11050-2947/2013/87(5)/054501(4) ©2013 American Physical Society

http://dx.doi.org/10.1103/PhysRevA.87.054501


BRIEF REPORTS PHYSICAL REVIEW A 87, 054501 (2013)

where Rij = |Ri − Rj | are interparticle distances and αk
ij are

nonlinear parameters. Superscript k in αk
ij reflects the fact that

the nonlinear parameters are unique for each basis function.
N is the total number of particles (six in our case).

ECG basis functions explicitly include all possible in-
terparticle distances, and each distance comes with its own
parameter. These nonlinear parameters can be tuned in order
to achieve a compact and accurate representation of the wave
function of the system. It is this flexibility due to a large number
of parameters that gives the ECG basis the ability to describe
very accurately the correlated motion of particles in essentially
arbitrary few-body systems.

Functions (1) are both rotationally and translationally
invariant. Effectively, they only describe the intrinsic motion of
the particles. The motion of the system as a whole is eliminated
from consideration. More details describing technical aspects
of calculations with ECG basis functions can be found in
Refs. [18–21] and references therein.

Polyelectrons, in particular those consisting of the same
number of electrons and positrons, possess very high permuta-
tional symmetry. The Hamiltonian of a Psn system is invariant
not only upon permutations of electrons or positrons but also
under charge conjugation (i.e., when all electrons are replaced
with positrons and vice versa). A detailed analysis of the
symmetry of Ps2 is presented in Ref. [26]. In our calculations
we do not explicitly deal with spin-components of the trial
wave function. Instead, we use Young projection operators
for both sets of identical particles, electrons and positrons
[27]. Imposing the charge conjugation symmetry requires an
additional factor in the symmetry projector, 1 ± P14P25P36

(here Pij is a permutation of spatial coordinates of particles
i and j ). For the ground-state calculations with zero angular
momentum, we chose the plus sign as such a configuration
is more energetically favorable. Thus, the total symmetry
projector applied to the basis functions is

(1 ± P14P25P36)(1 − P46)(1 − P13)(1 + P45)(1 + P12). (2)

A common way to study the stability of weakly bound
Coulomb few-body systems is to increase or decrease masses
or charges of particles so that the system under study becomes
relatively strongly bound. Then, assuming that the energy of
that bound state is a continuous function of the masses and
charges, they are changed slowly to approach their actual
values. At the same time the wave function of the system is
evolved gradually by adjusting the shape of the basis functions
and increasing, if necessary, the length of the variational
expansion. This provides a practical possibility to maintain
the boundedness of the trial wave function as long as the
Hamiltonian allows for a bound state. Variational calculations
with ECGs are particularly suitable in such kind of calculations
as they allow very efficient optimization of the basis at each
step.

Because of the symmetry of the Ps3 system, scaling of
all charges and (or) masses simultaneously does not change
the degree of boundedness. Effectively, it only redefines the
units of distance and energy. On the other hand, increasing or
decreasing the charge or mass of only one sort of particles (i.e.,
either positrons or electrons) does not leave the Hamiltonian
invariant with respect to charge conjugation. Therefore, in

order to preserve the symmetry properties of the system and
at the same time be able to vary the boundedness, we split
the potential energy operator into two parts, the attractive
interaction and the repulsive interaction. The Hamiltonian then
looks as follows (atomic units are assumed throughout):

H = −1

2

N∑

i=1

∇2
Ri

+ βVattr + γVrep, (3)

where

Vattr = − 1

R14
− 1

R15
− 1

R16
− 1

R24
− 1

R25

− 1

R29
− 1

R34
− 1

R35
− 1

R36
, (4)

Vrep = 1

R12
+ 1

R13
+ 1

R23
+ 1

R45
+ 1

R46
+ 1

R56
. (5)

Parameters β and γ define the strength of attraction and
repulsion between particles. In the real Ps3 system, they both
are equal to unity. In our calculations we vary γ , starting
with a value γ < 1 that corresponds to a well-bound state
and increasing it gradually up to γ = 1. Parameter β is kept
equal to unity at all times. The bound state of the system could
certainly be maintained in the opposite way, i.e., by keeping
γ = 1 and β > 1.

III. RESULTS

Table I shows the binding energies of both Ps2 and Ps3

as functions of the repulsion strength. The calculations for
the smaller system, Ps2, are converged to beyond the seventh
significant figure regardless of the value of γ when only 500
basis functions are used. Calculations of Ps3 are considerably
more time-consuming due to higher number of degrees of
freedom in that system. In addition, they require the use of
larger basis sets. We increased the size of the basis gradually
from 1000 for a well-bound state with γ = 0.8 to 2000 in
the case of γ = 1.0 in order to provide sufficient accuracy
for higher values of γ and make sure that the relative error
in the binding energy remains small enough. As one can
see from the trends of both the total and binding energies
in Table I, Ps2 remains a relatively well-bound system for
any γ , while the energy of Ps3 approaches the dissociation
threshold when γ → 1. The dissociation of Ps3 in the limit
γ = 1 is seen even better if we look at the trends of the
expectation values 〈r++〉 and 〈r+−〉, which stand for the
average distance between like particles and oppositely charge
particles respectively. One can clearly see that the system gets
broken into two subsystems upon γ approaching unity. In fact,
recalling that this six-particle Ps3 system has 15 interparticle
distances total (6 corresponding to the like particles and 9 to the
oppositely charge particles) it is easy to deduce that the ratio
of 〈r++〉/〈r+−〉 ≈ 1.5 in the limit of two separated fragments
corresponds to the dissociation channel Ps3 → Ps2+ Ps.

Assuming that the Ps3 → Ps2+ Ps dissociation channel
indeed takes place, it should be possible to estimate the
binding energy as a function of γ with a quick calculation.
In this case we effectively have a two-body system with
masses 4 and 2 (in atomic units), which interact via Coulomb
potential V (R) = 4(β − γ )/R, where R is the separation

054501-2



BRIEF REPORTS PHYSICAL REVIEW A 87, 054501 (2013)

TABLE I. The total and binding energies and average interparticle distances of Ps2 and Ps3 molecules as a function of the repulsion strength
parameter γ . All quantities are in a.u. The values in parentheses are estimates of the remaining uncertainty due to finite size of the basis used
in calculations.

Ps2 Ps3

γ Basis Etot Ebind 〈r++〉 〈r+−〉 Basis Etot Ebind 〈r++〉 〈r+−〉
0.8 500 −0.6221182(0) 0.1221182(0) 4.264(0) 3.478(0) 1000 −0.9731482(50) 0.1010299(50) 6.014(0) 4.906(0)
0.9 500 −0.5651750(0) 0.0651750(0) 4.831(0) 3.823(0) 1000 −0.8537477(50) 0.0385727(50) 7.644(2) 6.041(2)
0.95 500 −0.5394609(0) 0.0394609(0) 5.267(0) 4.073(0) 1000 −0.8040834(50) 0.0146224(50) 9.667(10) 7.412(10)
0.98 500 −0.5250743(0) 0.0250743(0) 5.654(0) 4.285(0) 1000 −0.7789411(50) 0.0038668(50) 13.61(10) 10.06(10)
0.985 500 −0.5227641(0) 0.0227641(0) 5.736(0) 4.329(0) 1000 −0.7752591(50) 0.0024951(50) 15.41(10) 11.26(10)
0.99 500 −0.5204813(0) 0.0204813(0) 5.826(0) 4.377(0) 1000 −0.7717932(50) 0.0013119(50) 18.76(10) 13.49(10)
0.993 500 −0.5191254(0) 0.0191254(0) 5.883(0) 4.408(0) 1000 −0.7698473(30) 0.0007219(30) 22.93(20) 16.27(20)
0.995 500 −0.5182274(0) 0.0182274(0) 5.924(0) 4.429(0) 1500 −0.7686265(30) 0.0003991(30) 28.68(50) 20.11(50)
0.997 500 −0.5173342(0) 0.0173342(0) 5.966(0) 4.452(0) 1500 −0.7674843(20) 0.0001501(20) 43.10(1.) 29.72(1.)
0.998 500 −0.5168895(0) 0.0168895(0) 5.988(0) 4.463(0) 1500 −0.7669540(10) 0.0000645(10) 63.19(5.) 43.11(5.)
0.999 500 −0.5164460(0) 0.0164460(0) 6.010(0) 4.475(0) 1500 −0.7664593(8) 0.0000133(8) 125.2(10.) 84.45(10.)
0.9995 500 −0.5162247(0) 0.0162247(0) 6.022(0) 4.481(0) 1500 −0.7662264(6) 0.0000017(6) 207.1(20.) 139.0(20.)
1.0 500 −0.5160038(0) 0.0160038(0) 6.033(0) 4.487(0) 2000 −0.7660021(15) Unbound 487.2(∞) 325.8(∞)

between two fragments, Ps2 and Ps. The binding energy is
then Ebind ∼ 32

3 (β − γ )2. For γ = 0.999 and γ = 0.995 this
gives 0.0000107 a.u. and 0.0000027 a.u. respectively, which
is in good agreement with the calculated binding energies of
0.0000133 a.u. and 0.0000017 a.u. listed in Table I.

While the results shown in Table I cannot be be considered a
rigorous mathematical proof of instability of Ps3, they provide
solid evidence that this system has no bound state with L = 0.

One could also think about the possibility of states with
L > 0. Dipositronium, Ps2, is known to have a stable excited
state with L = 1 [23,28]. However, the existence of that state
is stipulated by a reduced dissociation threshold. Ps2 (L = 1)
cannot decay into two Ps atoms in the ground state due
to the parity conservation. The situation is different in Ps3,
however, because this system dissociates into two nonidentical
fragments. Nevertheless, we have performed Ps3 (L = 1)
calculations for some values of γ and did not observe an
indication that there could be a bound state.

We have also investigated the possibility of resonance states
around the Ps2+ Ps threshold using the complex scaling
method [21]. No narrow resonances were found in the
spectrum. This concerns both the case of pure Coulomb
interaction (β,γ = 1) and the case when the repulsion strength
is weakened (β = 1, γ < 1).

In summary, we have performed a systematic study of the
stability of tripositronium by varying the interaction strength
between particles. The results of our variational calculations,
in which we expanded the wave function of a six-particle
Ps3 system in terms of explicitly correlated Gaussian basis
functions, indicate that tripositronium is not stable against
dissociation.

ACKNOWLEDGMENT

S.B. and O.V.P acknowledge financial support from the
National Science Foundation, Grant CHE-1300118.

[1] J. A. Wheeler, Ann. N.Y. Acad. Sci. 48, 219 (1946).
[2] S. Mohorovičić, Astron. Nachr. 253, 93 (1934).
[3] A. E. Ruark, Phys. Rev. 68, 278 (1945).
[4] M. Deutsch, Phys. Rev. 82, 455 (1951).
[5] A. P. Mills, Phys. Rev. Lett. 46, 717 (1981).
[6] A. P. Mills, Phys. Rev. Lett. 50, 671 (1983).
[7] E. A. Hylleraas and A. Ore, Phys. Rev. 71, 493 (1947).
[8] D. B. Cassidy and A. P. Mills, Nature (London) 449, 195 (2007).
[9] D. B. Cassidy, S. H. M. Deng, and A. P. Mills, Phys. Rev. A 76,

062511 (2007).
[10] C. Surko, Nature (London) 449, 153 (2007).
[11] A. P. J. Mills, D. B. Cassidy, and R. G. Greaves, Mater. Sci.

Forum 445, 424 (2004).
[12] J. Mitroy, M. W. J. Bromley, and G. G. Ryzhikh, J. Phys. B 35,

R81 (2002).
[13] G. F. Gribakin, J. A. Young, and C. M. Surko, Rev. Mod. Phys.

82, 2557 (2010).

[14] M. W. J. Bromley and J. Mitroy, Phys. Rev. A 66, 062504
(2002).

[15] K. Strasburger and H. Chojnacki, in Explicitly Corre-
lated Wave Functions in Chemistry and Physics, edited by
J. Rychlewski, Progress in Theoretical Chemistry and Physics,
Vol. 13 (Springer, Amsterdam, the Netherlands, 2003),
pp. 439–463.

[16] M. V. Pak, A. Chakraborty, and S. Hammes-Schiffer, J. Phys.
Chem. A 113, 4004 (2009).

[17] V. A. Dzuba, V. V. Flambaum, G. F. Gribakin, and C. Harabati,
Phys. Rev. A 86, 032503 (2012).

[18] K. Varga and Y. Suzuki, Phys. Rev. C 52, 2885 (1995).
[19] Y. Suzuki and K. Varga, Stochastic Variational Approach to

Quantum-Mechanical Few-Body Problems, Lecture Notes in
Physics (Springer, Berlin, 1998).

[20] S. Bubin, M. Pavanello, W.-C. Tung, K. L. Sharkey, and
L. Adamowicz, Chem. Rev. 113, 36 (2013).

054501-3

http://dx.doi.org/10.1111/j.1749-6632.1946.tb31764.x
http://dx.doi.org/10.1002/asna.19342530402
http://dx.doi.org/10.1103/PhysRev.68.278
http://dx.doi.org/10.1103/PhysRev.82.455
http://dx.doi.org/10.1103/PhysRevLett.46.717
http://dx.doi.org/10.1103/PhysRevLett.50.671
http://dx.doi.org/10.1103/PhysRev.71.493
http://dx.doi.org/10.1038/nature06094
http://dx.doi.org/10.1103/PhysRevA.76.062511
http://dx.doi.org/10.1103/PhysRevA.76.062511
http://dx.doi.org/10.1038/449153a
http://dx.doi.org/10.4028/www.scientific.net/MSF.445-446.424
http://dx.doi.org/10.4028/www.scientific.net/MSF.445-446.424
http://dx.doi.org/10.1088/0953-4075/35/13/201
http://dx.doi.org/10.1088/0953-4075/35/13/201
http://dx.doi.org/10.1103/RevModPhys.82.2557
http://dx.doi.org/10.1103/RevModPhys.82.2557
http://dx.doi.org/10.1103/PhysRevA.66.062504
http://dx.doi.org/10.1103/PhysRevA.66.062504
http://dx.doi.org/10.1021/jp810410y
http://dx.doi.org/10.1021/jp810410y
http://dx.doi.org/10.1103/PhysRevA.86.032503
http://dx.doi.org/10.1103/PhysRevC.52.2885
http://dx.doi.org/10.1021/cr200419d


BRIEF REPORTS PHYSICAL REVIEW A 87, 054501 (2013)

[21] J. Mitroy, S. Bubin, W. Horiuchi, Y. Suzuki, L. Adamowicz,
W. Cencek, K. Szalewicz, J. Komasa, D. Blume, and K. Varga,
Rev. Mod. Phys. 85, 693 (2013).

[22] G. G. Ryzhikh and J. Mitroy, Phys. Rev. Lett. 79, 4124 (1997).
[23] K. Varga, J. Usukura, and Y. Suzuki, Phys. Rev. Lett. 80, 1876

(1998).
[24] K. Varga, Phys. Rev. Lett. 83, 5471 (1999).

[25] E. Armour, J.-M. Richard, and K. Varga, Phys. Rep. 413, 1
(2005).

[26] D. M. Schrader, Phys. Rev. Lett. 92, 043401 (2004).
[27] M. Hamermesh, Group Theory and Its Application to Physical

Problems (Addison-Wesley, Reading, MA, 1962).
[28] J. Usukura, K. Varga, and Y. Suzuki, Phys. Rev. A 58, 1918

(1998).

054501-4

http://dx.doi.org/10.1103/RevModPhys.85.693
http://dx.doi.org/10.1103/PhysRevLett.79.4124
http://dx.doi.org/10.1103/PhysRevLett.80.1876
http://dx.doi.org/10.1103/PhysRevLett.80.1876
http://dx.doi.org/10.1103/PhysRevLett.83.5471
http://dx.doi.org/10.1016/j.physrep.2005.02.003
http://dx.doi.org/10.1016/j.physrep.2005.02.003
http://dx.doi.org/10.1103/PhysRevLett.92.043401
http://dx.doi.org/10.1103/PhysRevA.58.1918
http://dx.doi.org/10.1103/PhysRevA.58.1918