**3.3 Simulation of a Brillouin neutron scattering experiment**

Let us imagine to perform a BNS measurement. The instrument will acquire scattering intensity for a certain time, splitting such an acquisition in separate runs of the same duration. We can focus now on the spectrum corresponding to a constant *Q* cut of the *S Q*ð Þ , *E* we are measuring. To visualize this, we can generate simulated experimental data as they were actually measured. **Table 1** provides the parameters


**Table 1.**

*Absolutes parameter values for the model from which the simulated datasets were drawn.*

*Bayesian Inference as a Tool to Optimize Spectral Acquisition in Scattering Experiments DOI: http://dx.doi.org/10.5772/intechopen.103850*

setting used in the simulation study. These otherwise arbitrary parameter values are chosen to reproduce a typical spectrum collected in a neutron or X-ray scattering measurement on an amorphous system. More specifically, we have tuned the parameters so to obtain barely resolved inelastic contributions to the spectrum, which is a typical problem faced in routine measurements. Indeed, because of the limited instrument resolution and finite counting statistics, whenever the excitation features are not sufficiently sharp, that is, adequately separated and broadened by small damping, they can blend to each other and partially merge into the dominant central peak, which makes their detection challenging. In these cases, the limited acquisition time of a single measurement become an even more constraining factor. Therefore, the chosen parameter set is suitable to mimic a typical scenario encountered in scattering studies on disordered systems. In **Figure 2**, we show one of these datasets randomly generated using the following simple model:

$$
\tilde{S}(Q,E) = R(E) \otimes S(Q,E) \tag{6}
$$

where <sup>⊗</sup> indicates the convolution product, *R E*ð Þ¼ <sup>1</sup>ffiffiffiffi <sup>2</sup>*<sup>π</sup>* <sup>p</sup> *<sup>ζ</sup>* exp � *<sup>E</sup>*<sup>2</sup> 2*ζ*<sup>2</sup> � � is the instrument resolution function, assumed to have a zero-centered Gaussian profile with *ζ*<sup>2</sup> variance, which in the present case gives a *FWHM* = 3.1 meV. Once again, we suppose to use the BRISP spectrometer with an incident energy *E*<sup>0</sup> = 80 meV, as achieved by using the (004) reflection from a Pyrolytic Graphite monochromator [15]. The dynamic structure factor is here approximated as:

$$S(Q,E) = A\_\epsilon(Q)\delta(E) + [n(E) + 1]\frac{E}{k\_B T} \left\{ \sum\_{k=1}^2 \frac{2}{\pi} A\_k(Q) DHO\_k(Q,E) \right\} \tag{7}$$

where *δ*ð Þ *E* is the Dirac Delta function describing the elastic response of the system modulated by an intensity factor *Ae*ð Þ *Q* , *DHOk* are *k* inelastic contributions to the spectrum described by Damped Harmonic Oscillator (DHO) functions [21],

#### **Figure 2.**

*Generated spectrum as drawn from the model in Eqs. (6) and (7) at a Q value of 5 nm*�<sup>1</sup>*. This spectrum simulates the data as they could appear after a very short acquisition run. The plotted quantity is in fact the scattered intensity, to which the dynamic structure factor is proportional.*

*n E*ð Þ¼ *eE=kBT* � <sup>1</sup> �<sup>1</sup> is the Bose factor expressing the detailed balance condition, *kB* is the Boltzmann constant, and we have chosen the temperature *T*=1337 *K*. Finally, the simulated experimental data points are corrupted by an additive random fluctuation *ε*ð Þ *Q*, *E* :

$$
\mathcal{Y}(Q,E) = \tilde{\mathcal{S}}(Q,E) + \varepsilon(Q,E), \tag{8}
$$

with *<sup>ε</sup>*ð Þ� *<sup>Q</sup>*, *<sup>E</sup>* <sup>N</sup> 0, *<sup>σ</sup>*2<sup>~</sup> *S Q*ð Þ , *<sup>E</sup>* , for any *<sup>Q</sup>* and *<sup>E</sup>*, where the symblol � means "distributed according to," N denotes the Gaussian distribution and *σ*<sup>2</sup> is a constant factor to be estimated on the experimental data.

We can then generate as many experimental runs as we wish and sum them as usually done. When a scattering measurement is actually performed in the large-scale facilities above mentioned, the beam time granted to researchers is probably the most critical requisite for a successful experiment. The number of equal duration measurement runs on a given sample provides the total time allotted for that sample. In **Figure 3**, we show data as they result from the sum of 20 runs of identical integration time to qualitatively visualize the improvement in data precision and spectral shape definition that can be achieved by enhancing the counting statistics through a factor 20 increase of the acquisition time.

We will try now to fit our experimental data using a Bayesian Markov Chain Monte Carlo (MCMC) [22] algorithm equipped with a Reversible Jump option (RJ) [23], as explained in detail in Refs. [10, 11]. This algorithm allows to draw values from a distribution which is only known up to a normalization constant and thus to simulate the joint posterior distribution of the parameter vector of the model, Θ, as defined in Eq. 7. The analytical evaluation of the normalization constant is in fact usually really hard if not impossible at all. We again stress that in this simplified model the number *k* of inelastic components contributing to the spectrum is in itself a free model parameter to be estimated conditionally on available data. Notice that the RJ option allows the MCMC algorithm to explore models with different numbers *k* of inelastic components with *k* ¼ 1 … *kmax*, *kmax* being the maximum number of excitations allowed. As a first step, the first measurement run is best fitted by the model and the first-level

**Figure 3.** *Sum of 20 generated spectra as drawn from the model in Eqs. (6) and (7).*
