  A complete calculation of phase space factors (PSF) for Majoron emitting
$0\nu\beta^-\beta^-$ decay modes is presented. The calculation makes use of
exact Dirac wave functions with finite nuclear size and electron screening and
includes life-times, single electron spectra, summed electron spectra, and
angular electron correlations. Combining these results with recent interacting
boson nuclear matrix elements (NME) we make half-life predictions for the the
ordinary Majoron decay (spectral index $n$=1). Furthermore, comparing
theoretical predictions with the obtained experimental lower bounds for this
decay mode we are able to set limits on the effective Majoron-neutrino coupling
constant $\langle g_{ee}^M\rangle$.
