  Despite its generality and powerful convergence properties, Milstein's method
for functionals of spatially bounded stochastic differential equations is
widely regarded as difficult to implement. This has likely prevented it from
being utilised in applications. In this paper, we design and analyse in detail
one such implementation. The presented method turns out to be on par with
other, popular schemes in terms of computational cost---but with a (nearly)
linear weak convergence rate under the usual smoothness requirements on
coefficients and boundary. Two byproducts of theoretical interest are a new,
non-standard rank-one update formula, and a connection between numerics of
bounded diffusions and Eikonal equations. Three examples are worked out,
confirming the accuracy and robustness of the method.
