Convergence issues in derivatives of Monte Carlo null-collision integral formulations: A solution

Abstract

When a Monte Carlo algorithm is used to evaluate a physical observable A, it is possible to slightly modify the algorithm so that it evaluates simultaneously A and the derivatives ∂ ς A of A with respect to each problem-parameter ς. The principle is the following: Monte Carlo considers A as the expectation of a random variable, this expectation is an integral, this integral can be derivated as function of the problem-parameter to give a new integral, and this new integral can in turn be evaluated using Monte Carlo. The two Monte Carlo computations (of A and ∂ ς A) are simultaneous when they make use of the same random samples, i.e. when the two integrals have the exact same structure. It was proven theoretically that this was always possible, but nothing insures that the two estimators have the same convergence properties: even when a large enough sample-size is used so that A is evaluated very accurately, the evaluation of ∂ ς A using the same sample can remain inaccurate. We discuss here such a pathological example:null-collision algorithms are very successful when dealing with radiative transfer in heterogeneous media, but they are sources of convergence difficulties as soon as sensitivity-evaluations are considered. We analyse theoretically these convergence difficulties and propose an alternative solution.first practical benefit is that the next collision event can be sampled as if the medium was homogenous. Then the choice is made to select a true-collision or a virtual-collision as function of their local respective-amounts and this is how the spatial information is recovered. But several other benefits were recently foreseen in [1] and practically tested in [5,6,7,8,9,10,11,12,13,14,15,16], mainly for radiative-transfer applications. The main idea is that null-collision algorithms transform the non-linearity of Beer-extinction into a linear-recursive problem that Monte Carlo handles without approximation [14]. This was for instance used in [5] to deal with absorption-spectra of molecular gases combining very numerous transitions: the summation over all transitions could be treated by the Monte Carlo algorithm itself, which was previously assumed impossible because this summation was inside the exponential of Beer-extinction. Similarly, the vanishing of the exponential allowed the extension of implicit Monte Carlo algorithms for inversion of absorption and scattering coefficients from intensity measurements [6]. Outside radiative transfer, a very similar idea was used to solve Electromagnetic Maxwell equations for energy propagation in particle-ensembles of statistically-distributed shapes despite of the nonlinearity associated to the square of the electric field [13]. Again similar is the algorithm proposed in [14] solving Boltzmann equation for micro-fluidics applications despite of the nonlinearity of the collision operator.Back to radiative-transfer applications, the ideas suggested in [1] have motivated significant developments in the computer-graphics community for the cinema industry. Here the benef...

Three viewpoints on null-collision Monte Carlo algorithms

“…Roger et al [21] have shown that this is always possible, regardless of the parameter type, and even when the parameter affects the integration domain. However, despite being thoroughly general, this proposition can cause practical difficulties in some contexts, in particular as far as numerical convergence is concerned [26].…”

“…where ∂ ζ k e is the derivative of k e with respect to ζ (see Lataillade et al [4]). It can be easily seen that − ∂ζ kê k−ke tends towards infinity whenk tends towards k e (for more details, see Tregan et al [26]), that is, when the probability of null collisions decreases. For this reason, the variance of the sensitivity estimator becomes infinite and convergence is impossible.…”

“…It also illustrates a solution that could easily be found by changing the viewpoint on null-collisions, which is described in thorough detail in Tregan et al 2020 [26]. In the first case, the null-collision probability needs to be derived, leading to the problematic − ∂ζ kê k−ke mentioned above; in the second case, only p(x) needs to be derived, in whichk does not appear, and the convergence difficulty associated with − 1 k−ke vanishes (see Figure 8).…”

“…• Tregan et al [26] propounded a solution to avoid a convergence issue due to the use of null-collision algorithms to compute sensitivity estimates of a radiative quantity. This was only possible by shifting from the intuitive to the probabilistic point of view on the null-collision method, thereby modifying the integral formulation associated with the first null-collision algorithm into a new formulation (and hence, algorithm) yielding better statistical properties.…”

Monte-Carlo and sensitivity transport models for domain deformation

“…To carry out the computation, the models are stated in an integral formulation, which is then considered as an expectation, which in turn is estimated by sampling, yielding an unbiased estimate of the expectation as well as its confidence interval. Considering sensitivities, a well-known advantage of the Monte-Carlo method is its ability to estimate such expectation and its derivatives by using the very same sampling, avoiding additional computation time [1,11,14,[17][18][19]26] . This advantage has however some limitations when sensitivity to geometrical parameters is considered because the integral formulation of the model raises formalization and implementation difficulties [25] .…”