  Two-Body Dirac equations of constraint dynamics provide a covariant framework
to investigate the problem of highly relativistic quarks in meson bound states.
This formalism eliminates automatically the problems of relative time and
energy, leading to a covariant three dimensional formalism with the same number
of degrees of freedom as appears in the corresponding nonrelativistic problem.
It provides bound state wave equations with the simplicity of the
nonrelativistic Schroedinger equation. Unlike other three-dimensional
truncations of the Bethe-Salpeter equation, this covariant formalism has been
thoroughly tested in nonperturbatives contexts in QED, QCD, and nucleon-nucleon
scattering. Here we continue the important studies of this formalism by
extending a method developed earlier for positronium decay into two photons to
tests on the sixteen component quarkonium wave function solutions obtained in
meson spectroscopy. We examine positronium decay and then the two-gamma
quarkonium decays of eta_c, eta'_c, chi_0c, chi_2c, and pi-zero The results for
the pi-zero, although off the experimental rate by 13%, is much closer than the
usual expectations from a potential model.
