A boundary integral operator method for stochastic ray tracing in billiards was recently proposed in [1]. In particular, a phase-space boundary integral model for propagating uncertain ray or particle flows was described and shown to interpolate between deterministic and random models of the flow propagation. In this work we describe discretisation schemes for this class of boundary integral operators using piecewise constant collocation
in the spatial variable and either the Nyström method or the collocation method in the momentum variable. The simplicity of the spatial basis means that the corresponding spatial integration can be performed analytically. Convergence properties of the discretisation schemes and strategies for numerical implementation are presented and discussed
