  Recently, QED bremsstrahlung in $B$-meson decays into pair of scalars (\pi's
and/or K's) is of interest. If experimental acceptance must be taken into
account, PHOTOS Monte Carlo is often used in experimental simulations. We will
use scalar QED to benchmark PHOTOS, even though this theory is of limited use
for complex objects. We present the analytical form of the kernel used in the
older versions of PHOTOS, and the new, exact (scalar QED) one. Matrix element
and phase-space Jacobians are separated in the final weight and future
extensions based on measurable electromagnetic form-factors are thus possible.
The massive phase-space is controlled in the the program with no
approximations. Thanks to the iterative solution all leading and next to
leading logarithmic terms are properly reproduced by the Monte Carlo
simulation. Simultaneously, full differential distributions over complete
multiple body phase-space are provided. An agreement of better than 0.01% with
independent calculations of scalar QED is demonstrated.
