Uncertainty Modeling and Analysis of the European X-ray Free Electron Laser Cavities Manufacturing Process

TL;DR

This paper presents a comprehensive uncertainty quantification and sensitivity analysis of the manufacturing process for TESLA-shaped superconducting cavities used in the European X-ray Free Electron Laser, combining measurements and advanced simulation techniques.

Contribution

It introduces a detailed comparison and integration of measurement data with simulation-based uncertainty analysis for accelerator cavities.

Findings

Identified key sensitivities in cavity manufacturing process.

Validated simulation models against measurement data.

Enhanced surrogate modeling techniques for uncertainty quantification.

Abstract

This paper reports on comprehensive efforts on uncertainty quantification and global sensitivity analysis for accelerator cavity design. As a case study object the TESLA shaped superconducting cavities, as produced for the European X-ray Free Electron Laser (EXFEL), are selected. The choice for these cavities is explained by the available measurement data that can be leveraged to substantiate the simulation model. Each step of the manufacturing chain is documented together with the involved uncertainties. Several of these steps are mimicked on the simulation side, e.g. by introducing a random eigenvalue problem. The uncertainties are then quantified numerically and in particular the sensitivities give valuable insight into the systems behavior. We also compare these findings to purely statistical studies carried out for the manufactured cavities. More advanced, adaptive, surrogate modeling techniques are adopted, which are crucial to incorporate a large number of uncertain parameters. The main contribution is the detailed comparison and fusion of measurement results for the EXFEL cavities on the one hand and simulation based uncertainty studies on the other hand. After introducing the quantities of physical interest for accelerator cavities and the Maxwell eigenvalue problem, the details on the manufacturing of the EXFEL cavities and measurements are reported. This is followed by uncertainty modeling with quantification studies.

Tables

Table 1. Table 1 : TM 010 spectra for the European XFEL cavities, data taken from desy_db .

Mode	Frequency, MHz
	Average		St Dev
	RI	EZ	RI	EZ
0	1275.832	1274.734	0.226	0.268
1	1277.951	1276.926	0.212	0.245
2	1281.194	1280.294	0.186	0.208
3	1285.178	1284.460	0.157	0.168
4	1289.423	1288.921	0.124	0.129
5	1293.442	1293.161	0.089	0.090
6	1296.762	1296.649	0.064	0.066
7	1298.944	1298.946	0.052	0.054
8	1299.702	1299.716	0.050	0.052

Table 2. Table 2 : Mean values and standard deviations of the first passband in the tuned configuration.

Mode	Mean [MHz]	Std. dev. [MHz]
$0$	$1,276.45$	$0.15$
$1$	$1,278.50$	$0.13$
$2$	$1,281.64$	$0.11$
$3$	$1,285.60$	$0.08$
$4$	$1,289.85$	$0.06$
$5$	$1,293.84$	$0.04$
$6$	$1,297.11$	$0.02$
$7$	$1,299.25$	$0.00$
$8$	$1,300.00$	$0.00$

Table 3. Table 3 : Sample mean and sample standard deviation of cell-to-cell coupling coefficient for N cav = 826 subscript 𝑁 cav 826 N_{\text{cav}}=826 cavities as well as statistics about corresponding estimated iris deformations.

Manufacturer	$𝔼\left[k_{\text{cc},i}\right]$	$\mathrm{Std}\left[k_{\text{cc},i}\right]$	$𝔼\left[\Delta R_{\text{ir},i}\right]$	$\mathrm{Std}\left[\Delta R_{\text{ir},i}\right]$
RI	1.854	0.016	$0.087\mathrm{mm}$	$0.057\mathrm{mm}$
EZ	1.941	0.021	$0.400\mathrm{mm}$	$0.073\mathrm{mm}$
RI+EZ	1.897	0.047	$0.243\mathrm{mm}$	$0.170\mathrm{mm}$

Equations

$$\nabla \times \mathbf{E}$$

$$\nabla \cdot {ϵ}_0\mathbf{E}$$

$$\mathbf{E}\times \mathbf{n}$$

$$\nabla \times \left(\nabla \times \mathbf{E}\right)$$

$$\mathbf{E}\times \mathbf{n}$$

$$k_{\text{cc}}=2\frac{f_{\pi }-f_0}{f_{\pi }+f_0},$$

$$\mathcal{F}=\frac{{\min }_{i=1,\dots ,N_{\text{c}}}|E_{\text{ax,max}}^{(i)}|}{{\max }_{i=1,\dots ,N_{\text{c}}}|E_{\text{ax,max}}^{(i)}|}$$

$$\mathbf{Y}=[\Delta R_{\text{eq}}^{(1)},\dots ,\Delta R_{\text{eq}}^{(9)},\Delta R_{\text{ir}}^{(1)},\dots ,\Delta R_{\text{ir}}^{(10)}].$$

$$\Delta f_{m}=\bigl{|}\mathbb{E}\bigl{[}f_{\text{RI},m}\bigr{]}-\mathbb{E}\bigl{[}f_{\text{EZ},m}\bigr{]}\bigr{|}$$

$$k_{\text{cc}}$$

$$\begin{array}{r}{\mathbf{Y}}_{\text{prod}}(\theta )=[\Delta R_{\text{eq}}^1(\theta ),\dots ,\Delta R_{\text{eq}}^9(\theta ),\\ \Delta R_{\text{ir}}^{(1)}(\theta ),\dots ,R_{\text{ir}}^{(10)}(\theta )]\end{array}$$

$$\left(\nabla \times \mathbf{E},\nabla \times \mathbf{v}\right)=k^2\left(\mathbf{E},\mathbf{v}\right)\forall \mathbf{v},$$

$$\left(\nabla \times {\mathbf{E}}_h,\nabla \times {\mathbf{v}}_j\right)=k_h^2\left({\mathbf{E}}_h,{\mathbf{v}}_j\right)\forall {\mathbf{v}}_j.$$

$${\mathbf{E}}_h=\sum _{j=1}^{N_{\text{dof}}}{e}_j{\mathbf{v}}_j.$$

$$\varrho (y)=\frac{140}{(u-l{)}^7}\{\begin{array}{ll}(y-l{)}^3(u-y{)}^3,& l<y<u,\\ 0,& \text{else},\end{array}$$

$$\sum_{i=1}^{N_{\text{c}}}\Delta L^{(i)}<$3\text{\,}\mathrm{m}\mathrm{m}$$

$

$$f_{\text{obj}}^{(i)}\bigl{(}L^{(i)}\bigr{)}=f_{0}^{(i)}\bigl{(}L^{(i)}\bigr{)}-$1.3\text{\,}\mathrm{GHz}$,$$

$$\mathbb{E}[k_{\text{cc}}]\approx 1.82802,\mathrm{Std}[k_{\text{cc}}]\approx \mathrm{0.01897.}$$

$$\Delta R_{\text{ir},i}\mapsto k_{\text{cc},i},i=1,\dots ,N_{\text{cav}},$$

$$\mathcal{Q}(\mathbf{Y})$$

$$+{\mathcal{Q}}_{1,2,\dots ,M}(Y_1,\dots ,Y_M),$$

$${\mathcal{Q}}_0$$

$${\mathcal{Q}}_i(Y_i)$$

$${\mathcal{Q}}_{i,j}(Y_i,Y_j)$$

$$\mathbb{V}[\mathcal{Q}]$$

$$+\mathbb{V}[{\mathcal{Q}}_{1,2,\dots ,M}(Y_1,\dots ,Y_M)].$$

$$S_i=\frac{\mathbb{V}[{\mathcal{Q}}_i(Y_i)]}{\mathbb{V}[\mathcal{Q}]}.$$

$$S_i\approx \frac{\frac{1}{N}{\sum }_{j=1}^N\mathcal{Q}({\mathbf{a}}_j)\mathcal{Q}({\mathbf{c}}_j)-{\mathbb{E}}_{\mathbf{A}}^{\text{MC}}[\mathcal{Q}{]}^2}{{\mathbb{V}}_{\mathbf{A}}^{\text{MC}}[\mathcal{Q}]},$$

$$s(Y_i)={\int }_{{\Xi }_{\mathcal{Q}}}|{\rho }_{\mathcal{Q}}(q)-{\rho }_{\mathcal{Q}|Y_i}(q)|\mathrm{d}q.$$

$${\delta }_i$$

$$\displaystyle=\frac{1}{2}\int_{\Xi}\rho_{Y_{i}}(y_{i})\bigg{[}\int_{\Xi_{\mathcal{Q}}}|\rho_{\mathcal{Q}}(q)-\rho_{{\mathcal{Q}}|Y_{i}=y_{i}}(q)|\,\mathrm{d}q\bigg{]}\,\mathrm{d}y_{i}.$$

$${\delta }_i=\frac{1}{2}{\int }_{\Xi \times {\Xi }_{\mathcal{Q}}}|{\rho }_{\mathcal{Q}}(q){\rho }_{Y_i}(y_i)-{\rho }_{\mathcal{Q},Y_i=y_i}(q,y_i)|\mathrm{d}q\mathrm{d}{y}_i.$$

$${\delta }_i\approx \frac{1}{2}\frac{1}{n}\sum _{j=1}^n|{\rho }_{\mathcal{Q}}(q^{(j)}){\rho }_{Y_i}(y_i^{(j)})-{\rho }_{\mathcal{Q},Y_i}(q^{(j)},y_i^{(j)})|.$$

Full text

††preprint: APS/123-QED

I INTRODUCTION

Accelerator devices require advanced, simulation based, design approaches due to demanding performance requirements and a considerable level of technical complexity. This is particularly true for superconducting accelerator cavities, which are a key element of many modern accelerator facilities. A typical design process involves 2D as well as 3D numerical solutions of the Maxwell eigenvalue problem, followed by optimization and uncertainty analysis and quantification studies. The latter have been conducted within the accelerator community from the 1970s, Halbach_1976aa ; Weiland_1977 ; vanRienen_1987 . However, these studies have been mainly based on (local) sensitivity analysis which should be applied with care to quantify uncertainties in the cavities’ geometry. Indeed, the eigenmodes and other measures of interest depend strongly on the shape of the cavity and local measures may not yield reliable results.

The topic of simulation based uncertainty quantification has seen tremendous developments in recent years, also in computational electromagnetics, see e.g. Clenet_2013aa . Nowadays, significant computational resources are available and uncertainty studies, taking into account systematically large parameter variations at all steps of the design process, come into reach. In particular, the concept of global sensitivity analysis saltelli2008global has received much attention. The variance-based approach to global sensitivity analysis, measures the contribution of each parameter (or parameter combination) to the variance of a system output quantity. These Sobol sensitivity indices permit not only to analyze the importance of model input parameters, which in turn is useful in guiding modeling efforts, but also to identify important combined high-order parameter variations. Despite their clear interpretation, Sobol coefficients are not readily extendable to settings with correlated inputs, which are useful in our case to model several parts of the manufacturing and assembly process. Such effects in turn can be quantified by Borgonovo indices Borgonovo_2007 which represent another approach to global sensitivity analysis, which has not received much attention in physics and engineering applications so far. We will apply both methods in the present study, however, their derivations will only be recalled in appendix A. Although the concept of global sensitivity analysis is quite well-established, the efficient computation of sensitivity indices for cavity applications is a difficult task, mainly due to the complexity of the underlying eigenvalue problem Georg_2019aa . This is addressed here by introducing surrogate models which emulate the relation between eigenmodes, or other quantities of interest on the model parameter.

In this respect, this paper reports on comprehensive efforts on uncertainty quantification and global sensitivity analysis for the simulation of TESLA shaped, superconducting, cavities. Such cavities have been produced in considerable quantity for the European X-ray Free Electron Laser (EXFEL) and measurement data is available to substantiate the approach. Each step of the manufacturing chain is documented together with the involved uncertainties. Some of these steps are mimicked on the simulation side, e.g. by introducing a random eigenvalue problem. The uncertainties are then quantified numerically and in particular the sensitivities give valuable insight into the systems behavior. We also compare these findings to purely statistical studies carried out during the manufacturing. However, the simulation of all manufacturing steps would require the solution of several random inverse problems and is postponed to future work.

Uncertainty studies in an accelerator physics context have been reported before, see Xiao_2007aa ; Schmidt_2014aa ; Corno_2015ac . In this work, we use more advanced, adaptive, surrogate modeling techniques, which are crucial to incorporate a large number of uncertain parameters. The main contribution is the detailed comparison and fusion of real manufacturing data on the one hand and simulation based uncertainty studies on the other hand. The paper also clearly points out important directions of future research, which would allow to further combine measurements and simulation.

The structure of the paper is given as follows: Section II introduces quantities of physical interest for accelerator cavities and the Maxwell eigenvalue problem. In Section III, details on the manufacturing of the EXFEL cavities and measurements are reported. This is followed by uncertainty modeling and quantification studies in Section IV and concluding remarks.

II Cavities and Maxwell Eigenvalue Problem

Accelerating cavities are devices used to accelerate particles to higher energies. Elliptical cavities are the accepted geometrical shape for particle velocities close to the speed of light ($\beta\approx 1$, where $\beta$ is the ratio of particle velocity to the speed of light in vacuum). The shape of an elliptical cell is defined by two elliptical arcs connected by a tangent straight line as shown in FIG. 1. The fundamental mode of the cavity is the TM010 mode that is typically used as the operating mode of the cavity.

In order to enhance the accelerating efficiency, multi-cell cavities are created by connecting several cells together via their irises (see FIG. 2 for the EXFEL cavity described in more detail in Section III). In a multi-cell cavity of $N_{\text{c}}$ cells, there are $N_{\text{c}}$ modes of each type due to multi-cell coupling. In this work, we are interested in the TESLA cavity shape Aune2000 , which is composed of $N_{\text{c}}=9$ cells.

Let $\Omega(\bm{Y})\subset\mathbb{R}^{3}$ refer to the inner domain of the multi-cell cavity with boundary $\partial\Omega(\bm{Y})$, where $\bm{Y}$ denotes a vector of shape parameters to be specified. The fields in the source-free, time-harmonic case, are given by Maxwell’s equations:

[TABLE]

where $\mathbf{E}$ and $\mathbf{H}$ denote the electric and magnetic field strength, $\varepsilon_{0}$ and $\mu_{0}$ are the permittivity and permeability of vacuum, respectively. The walls are modeled as perfect electric conductor, i.e.,

[TABLE]

One derives the Maxwell eigenvalue problem from (1) by eliminating $\mathbf{H}$. Introduction of the wave number $k=2\pi f\sqrt{\mu_{0}\epsilon_{0}}$ yields

[TABLE]

for $\mathbf{E}\neq 0$ and $\nabla\cdot{\mathbf{E}}=0$. One should be aware that, although not explicitly specified, the field $\mathbf{E}$ and the eigenfrequencies $f$ depend on the shape parameters $\bm{Y}$. For each mode, i.e., a solution of (3) in the passband, there is a phase shift between fields of neighboring cells that can vary from 0 to $\pi$ radians. The $\pi$-mode, with frequency $f_{\pi}$, of the TM010 passband is used in multi-cell cavities for acceleration.

In order to maximize the voltage across the cavity, the length of the middle-cells $L$ is fixed as $L=\beta\lambda/2$, where $\lambda$ refers to the wavelength of the $\pi$-mode, see Padamsee_2015aa . Thus as the particle traverses a cell, the direction of the field reverses such that the particle is exposed to an identically directed electric field along the whole multi-cell cavity.

The modes in a passband have a small frequency difference. If the frequencies of the modes in the fundamental passband are very close to each other, there is a risk of exciting a mode close to $f_{\pi}$ by the RF generator.

The spread of modes in the first passband is reflected in the cell-to-cell coupling coefficient, which is defined as Belomestnykh2005

[TABLE]

where $f_{0}$ refers to the lowest frequency in the passband. The cell-to-cell coupling coefficient is dimensionless and a sensitive quantity in the design phase. A large aperture radius ($R_{\text{ir}}$) typically gives rise to a stronger cell-to-cell coupling.

If the energy of the fundamental mode is evenly distributed in the cells, the accelerating voltage is maximized (Padamsee1998, , p.129). Furthermore, a uniform field distribution allows for higher field magnitudes before reaching the surface electromagnetic (EM) field limit. The field flatness is a central figure of merit that indicates the uniformity of the field distribution of the fundamental mode between the cells. In this paper, the field flatness is defined as

[TABLE]

where $E_{\text{ax,max}}^{(i)}$ refers to the maximum axial electric field in cell $i$.

The resonant frequencies strongly depend on the geometry parameter in each cell, i.e., the equatorial radii $R_{\text{eq}}^{(i)}$, with $i=1,\ldots,N_{\text{c}}$, and the iris radii $R_{\text{ir}}^{(j)}$, with $j=1,\ldots,N_{\text{c}}+1$. In Section IV we consider perturbations of the form $R_{\text{eq}}^{(i)}+\Delta R_{\text{eq}}^{(i)}$ and $R_{\text{ir}}^{(j)}+\Delta R_{\text{ir}}^{(j)}$ for the equatorial radii and the iris radius, respectively. Then, the parameter vector is given as

[TABLE]

These perturbations change the resonant frequency of the respective cell(s) and consequently affect the frequency and the field distribution of the $\pi$-mode in the multi-cell cavity. It has been observed, that for the $\pi$-mode, the change in the field amplitude of each cell is proportional to the frequency change by a factor of $\propto 1/k_{\text{cc}}$ Sekutowicz1999 ; Nagle1967 . Thus, a small cell-to-cell coupling increases the sensitivity of the field profile with respect to geometrical perturbations.

III Cavity Manufacturing

The European X-ray Free Electron Laser XFEL_report facility is constructed to produce X-ray pulses with the properties of laser light and at intensities much brighter than those produced by conventional synchrotron light sources. The superconducting linear accelerator of the EXFEL has a length of almost 2.1 km and brings electrons to an energy of up to $17.5\text{\,}\mathrm{G}\mathrm{e}\mathrm{V}$. This is achieved by using a total of $808$ superconducting cavities installed in the three main linac sections and the injector. The production of $N_{\text{cav}}>808$ cavities Singer2016 , the largest in the history of cavity production, was realized by the two companies Research Instruments GmbH (RI) and Ettore Zanon S.p.A. (EZ). EXFEL uses nine-cell TESLA cavities build from solid niobium with a nominal $f_{\pi}=1,300$ MHz. Each cavity (see FIG. 2) consists of 10 main sub-components, welded together at the equator area. The sub-components consist of $8$ dumb-bells (DB) and $2$ end-groups (EGL, EGS), referring to short and long end-groups, respectively. A different shape of the end-groups’ half-cells provides the desired asymmetry of the Higher Order Mode (HOM) field distributions and increases the efficiency of their extraction.

Geometric deviations of the inner cavity shape, the cavity length, and the spectra of frequencies as well as deviations in HOM field distributions occur due to random inaccuracies during manufacturing. These uncertainties have a strong impact on the quantities of interest described in Section II. Hence, dedicated measures are undertaken during production to ensure acceptable tolerances according to the EXFEL cavity specification, in particular to obtain $\mathcal{F}>90\%$ and to keep the deviation of $f_{\pi}$ below $100\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$. These measures, together with sources of uncertainty, are summarized in FIG. 3 and described in detail in the following.

Step: Production

In this step $8N_{\text{cav}}$ DBs, consisting of $2$ half cells each, $N_{\text{cav}}$ EGSs and $N_{\text{cav}}$ EGLs are produced. 2. 2.

Step: Trimming

The target of this step is to compensate for shape deviations by trimming the components. It is applied to all components (DB, EGS and EGL) and allows to obtain the necessary cavity length and frequencies with required accuracy. 3. 3.

Step: Selection and Sorting

The manufacturer selects two end groups and $8$ dumb-bells. To minimize the influence of shape variations for the eight DBs on the asymmetry of the HOM field distribution, the DBs are sorted during the cavity assembly: a DB with average frequency is installed at the last position (position 8), the remaining DBs are installed in order of decreasing frequency (see FIG. 4).

Step: Welding

All components (EGL, DBs, EGS) are welded to each other, which induces shape deformations. The materials from the different suppliers exhibit different shrinkages at the welding joint. The resulting equator diameters may be slightly different depending on the cavity position during welding, see Sulimov2015 . 5. 5.

Step: Chemical treatment

The chemical treatment removes impurities and spikes, see Singer2016 , Sulimov2016 . The homogeneity of the removed material from the cavity surface depends on many parameters of the process and the facility. The electrochemical polishing treatment is usually less homogeneous and more unstable than the equator welding. 6. 6.

Step: Tuning

The cavity is tuned, i.e., mechanically stretched or compressed, according to the procedure described by Kreps_1996aa which adjusts $f_{\pi}$ with an accuracy of $\pm$50\text{,}\mathrm{k}\mathrm{H}\mathrm{z}$$ and ensures $\mathcal{F}>98\%$ for the field flatness. 7. 7.

Step: Final preparation (for operation)

The procedures applied in this step vary for different manufacturers and are shown in detail in FIG. 2 of Singer2016 . Those procedures include, e.g., final buffered chemical polishing etching, the integration of the cavities into the Helium tanks, a pressure test using water under the pressure of $6$ bar and the cool down to $2\text{\,}\mathrm{K}$. 8. 8.

Step: Cavity

All $N_{\text{cav}}$ cavities are operational and the statistics of the fundamental mode spectra are measured under $2\text{\,}\mathrm{K}$. A description of the measurement procedures is given below.

Quality assurance by mechanical measurements of the inner surface dimensions becomes impossible after cavity welding and polishing. Measurement data can only be obtained by ultrasonic or RF measurements. These methods are used for control of equator welding stabilities Sulimov2015 or homogeneity of the polishing process Sulimov2016 . In the following the fundamental mode spectra are discussed for both manufacturers (RI and EZ) based on RF measurements.

The measurements were carried out at operation temperature, i.e., cryogenic tests at 2 $\mathrm{K}$ were used (FIG. 5). The measurement results are presented in Table 1. Note, that the $\pi-$mode finally operates at 1.3 $\mathrm{G}\mathrm{H}\mathrm{z}$, which is ensured by an additional tuning step (not described here). The results show that the standard deviation of the $\pi-$mode (mode 8) in the relaxed condition after cool down, is about 50 $\mathrm{k}\mathrm{H}\mathrm{z}$. This value increases with decreasing mode number. A similar behavior can be found for the differences

[TABLE]

of the averaged frequencies for each mode $m$ between the two manufacturers RI and EZ, which also increases with decreasing mode number, see Tab. 1. In particular, the difference of average frequencies is about 14 $\mathrm{k}\mathrm{H}\mathrm{z}$ for mode 8 and 1098 $\mathrm{k}\mathrm{H}\mathrm{z}$ for mode 0, respectively. These results indicate a high deviation of the cell-to-cell coupling coefficient between the manufacturers with

[TABLE]

which are given in the form of average value plus/minus one standard deviation as obtained from database, see also Tab. 3. Preliminary studies Sulimov2013 show that only variations of the iris radius can explain such high deviation in the cell-to-cell coupling coefficient. This will be analysed in more detail numerically in Section IV. Finally, the field flatness distribution for the EXFEL cavities is presented in FIG. 6. It can be observed that at least $70\%$ of the cavities possess a field flatness of more than $95\%$.

IV Uncertainty quantification

A stochastic setting is adopted here to model manufacturing and measurement uncertainties and assess their influence on the cavity design, i.e. each component of $\bm{Y}$ becomes a random variable

[TABLE]

where $\theta$ denotes a random outcome.

Uncertainty quantification encompasses various methods for uncertainty propagation, Bayesian inverse problems, optimal experimental design and robust optimization, among others. The reader is referred to ghanem1991 ; xiu2010 ; lemaitre2010 for a detailed background. In this work, uncertainty propagation is of great interest, in particular, propagating distributions of cell deformation parameters to distributions of quantities of interest, such as the cell-to-cell coupling coefficient. This can be achieved by sampling according to the underlying distribution and repetitively solving the cavity eigenvalue problem. Thereby, surrogate modeling is a key step to keep the computational workload manageable. The probability distributions of the cavity geometry parameters are modeled based on descriptions of the manufacturing process, as described in Section III. A more general approach would consist in inferring these input distributions from measured RF data, which would require the solution of an inverse cavity eigenvalue problem. Such a study is out of the scope of the present work and can only be briefly addressed here, applying significant simplifications. In fact, the present study, should be considered as a step towards a more complete treatment of uncertainties in the cavity eigenvalue problem.

To solve Maxwell’s equations we use the variational formulation of (3). It reads: find $k$ and $\mathbf{E}$ such that

[TABLE]

where we assume that $\mathbf{E}$ and $\mathbf{v}$ are from the functions space of square-integrable vector fields with square-integrable curl. For further information on function spaces in the context of Maxwell’s equations, the reader is referred to Monk_2003aa . In order to numerically solve (8) we use a finite-dimensional spaces such that

[TABLE]

with

[TABLE]

where $N_{\text{dof}}$ denotes the number of degrees of freedom.

In this paper, calculations are carried out with the help of Superlans code SLANS . Superlans is a 2D-axisymmetric finite-element-based code used for the calculation of the monopole modes of azimuthally symmetric geometries. For each simulation, the data describing the contour of the cavity is created by Matlab Matlab and saved in a format readable by Superlans. Superlans is then called by Matlab to triangulate the geometry, solve the resulting eigenvalue problem and calculate the relevant secondary parameters. The results are finally read by Matlab for the post-processing. The detailed description of the uncertainty modeling steps, which are the numerical counterpart of the production chain in Section III, are as following:

Step: Production

We generate $\tilde{N}_{\text{cav}}=10^{6}$ random (”virtual”) cavities by creating, in turn, seven independent random mid cells and two end cells per cavity. Since In particular, those virtual cavities are obtained by drawing random numbers for the vector $\bm{Y}_{\text{prod}}$ defined in (7) where we assume independence of the random variables. Note that this assumption can often be justified for mass-production processes, in particular, as the individual dumb-bells are produced independently of each other. Although, the normal distribution seems to be the right choice, we opt for beta distribution in the range $-0.3\text{\,}\mathrm{m}\mathrm{m}0.3\text{\,}\mathrm{m}\mathrm{m}$. Beta distributions can be used to approximate normal distributions but have bounded support (xiu2010, , Appendix B), which is very important in numerical studies to avoid non-physical parameter configurations. FIG. 7 represents such an approximation with probability density function (PDF)

[TABLE]

where $l=$-0.3\text{,}\mathrm{m}\mathrm{m} and $u=$0.3\text{\,}\mathrm{m}\mathrm{m} denote the lower and upper bound, respectively. In this particular case the shape parameters of the beta distribution were chosen such that a normal distribution with a $2\sigma$ interval of $-0.2\text{\,}\mathrm{m}\mathrm{m}0.2\text{\,}\mathrm{m}\mathrm{m}$ is approximated, as illustrated in FIG. 7.

We note that the large number of virtual cavities $\tilde{N}_{\text{cav}}$ does not lead to prohibitive computational cost, as surrogate modeling is employed which will be explained later. 2. 2.

Step: Trimming

In contrast to real manufacturing, in the simulation approach, we build cavities out of elementary cells instead of dumb-bells. In this setting, it is not clear how to explicitly model the trimming, however it is implicitly taken into account by the cavity length limitation (11) which is considered in the next step. 3. 3.

Step: Selection and Sorting

The virtual 9-cell cavity is obtained from the random building blocks of step 1 as follows: We apply a sorting procedure which computes for all middle cells the fundamental resonance frequency (by solving a one-cell eigenvalue problem), places the cell closest to the average frequency at the 8${}^{\textrm{th}}$ position and orders the remaining cells according to decreasing frequency between positions 2 and 7. Additionally, a total length constraint

[TABLE]

is enforced. If the constraint is violated after tuning (step 6), the corresponding virtual cavity is disregarded. In real manufacturing, the constraint is already incorporated by trimming and compensation effects.

Step: Welding

Since we consider cells instead of dumb-bells for the simulation, the welding is modeled by averaging the iris radii of adjacent cells. Note that the impact of eccentric deformations by (possible) miss-alignment of cells has been investigated in Georg_2019aa . It was found that this type of deformation has negligible impact on the fundamental mode spectrum and is hence omitted here.

The steps of selection, sorting and welding influence the probability distribution. Hence, we introduce a new random vector $\bm{Y}_{\text{sort}}(\theta)$ where the density is estimated with kernel density techniques. Kernel density estimation is a non-parametric technique to infer a continuous probability density function from a sample. Here, we employ an Epanechnikov kernel Epanechnikov on the selected (and sorted) sample $\{\bm{Y}^{(m)}_{\text{sort}}\}_{m=1}^{\tilde{N}_{\text{sel}}}$ of size $\tilde{N}_{\text{sel}}=809641$, which complies with the length constraint. The estimated densities are presented in Fig. 8. 5. 5.

Step: Chemical treatment

An appropriate modeling of chemical treatment would require very fine resolutions or even multi-scale analyses which cannot be carried out in the present setting, see also Sulimov2016 for further information. 6. 6.

Step: Tuning

The virtual tuning procedure we apply, considers each cell individually and is therefore unable to incorporate field flatness constraints directly. However, we observe that an acceptable field flatness of at least $96\%$ was obtained for all eigenvalue problems that are eventually solved. In particular, each cell is tuned to $1.3\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ by changing its length $L^{(i)}$. Computationally, this requires the solution of a non-linear root finding problem for the objective function

[TABLE]

where $f_{0}^{(i)}$ denotes the fundamental eigenfrequency of cell $i$ (one-cell eigenvalue problem). The root finding problem is solved using fzero in Matlab. However, we note that different root finding schemes, e.g., Newton’s method or bisection, could be employed as well.

In FIG. 9 the values of the 9-cell cavity accelerating frequency are depicted for different choices of the length and equatorial radius of the first cell. Considering the $1.3\text{\,}\mathrm{GHz}$ contour (in black in the figure), we observe that the tuning process correctly identifies a length such that a value of $1.3\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ is obtained (the magenta dots in the figure lie on the contour line as expected).

Step: Final preparation (for operation)

Due to insufficient data and involved numerical modeling this step cannot be carried out in the setting employed in this work. 8. 8.

Step: Cavity

The Maxwell eigenvalue problem is solved to obtain the first nine modes for each virtual cavity and statistics of the spectra are computed. It shall be noted that, in general, eigenvalue tracking Georg_2019aa should be employed to ensure a consistent matching of the eigenfrequencies. However, it has been observed that, in this case, the fundamental eigenfrequencies do not cross with respect to parameter changes for the considered variations. Hence, mode tracking is not applied here. The flow diagram of simulation steps is given in Fig. 10.

In order to avoid the tremendous computational cost of repeatedly solving a large number of eigenvalue problems in step 3 and step 8, surrogate modeling is employed. In this work a dimension-adaptive interpolation scheme narayan2014 is used to construct accurate polynomial surrogate models. In particular, the algorithm described in (georg2018, , Chapter 3.2) is adapted to address the case of multiple quantities of interest. This algorithm constructs polynomial models adaptively which yields high computational efficiency. Moreover, we control the associated approximation errors by cross-validation to be sufficiently small. For convenience of the reader, we recall the main ideas of the employed algorithm in appendix B and refer for details to georg2018 .

The 4-variate polynomial surrogate models employed in step 3 are computed by solving 50 one-cell eigenvalue problems while 500 evaluations of the 9-cell eigenvalue problem are used to construct the 19-variate polynomial surrogate model for step 8. For all surrogate models, cross-validation indicates an error below $10\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$ for all fundamental resonance frequencies which is smaller than the standard expected deviations.

Now, all statistical quantities of interest can be obtained from the sample computed in step 8. In particular, we employ nbiased statistical estimates for expectation and standard deviation, which are depicted in FIG. 11 and Table 2. In the same way, we estimate expectation and standard deviation of the cell-to-cell coupling coefficient $k_{\text{cc}}$

[TABLE]

In the following, we employ Sobol indices as sensitivity measure, despite the fact that they are defined for independent parameters. We use the OpenTurns implementation of Saltelli’s algorithm saltelli2002 to estimate Sobol indices of the cell-to-cell coupling coefficient for the $10$ input parameters $\bm{Y}$. Readers interested in details on this algorithm are referred to appendix A. The results presented in FIG. 12 (top) are obtained by evaluating $4\cdot 10^{7}$ times the surrogate model. It can be observed that deviations in the middle cells, in particular of the respective iris radii, have significantly larger contributions to the cell-to-cell coupling coefficient. Additionally, we estimate Borgonovo indices, which represent another global sensitivity index family taking correlations in the input data into account. To this end, we employ an approach based on kernel density estimation, see for instance da2015global and references therein, on the sample $\{\bm{Y}_{\text{sort}}^{(m)}\}_{i=1}^{\tilde{N}_{\text{sel}}}$ and the corresponding surrogate model evaluations of $k_{\text{cc}}$. Contrary to Sobol indices which decompose the output variance with regard to different model input contributions, in the Borgonovo approach, the dependency between input and output parameters is analysed, see Appendix A for more details. Hence, different sensitivity magnitudes can be expected in each case. However, also in the case of Borgonovo sensitivity analysis, similar results are obtained. In particular, the variations in the middle cells can be identified to have a larger impact than the corresponding end-cell variations. The surprisingly large contribution $\Delta R_{\text{eq}}^{(8)}$, in comparison to the variance-based analysis, can most probably be explained by the special treatment of cell $8$ in the sorting procedure. Since we place the cell closest to the average frequency at position $8$, a large variation of $\Delta R_{\text{eq}}^{8}$ might implicitly also imply larger variations of the other middle cells. To confirm this interpretation, we additionally repeated this study with a different sorting procedure which omits the special treatment of cell $8$. In this case, the sensitivity of $\Delta R_{\text{eq}}^{(8)}$ was reduced drastically, as expected.

As in Sulimov2013 , we proceed by making the simplifying assumption, that, for the rest of the paper, all iris radii deviations $\Delta R_{\text{ir}}^{(i)},\,i=1,\ldots,N_{\text{c}}+1$ of a cavity are the same and accordingly drop the index $i$. This allows us to carry out a preliminary inverse analysis, i.e., to compute geometric variations from the cell-to-cell coupling coefficient. To this end, we first repeat the computation of sensitivity indices, by rerunning the full simulation workflow, where we now model the overall iris radii deviation $\Delta R_{\text{ir}}$ as one beta distributed random variable in the range of $\pm 0.3\,\text{nm}$. The corresponding Sobol indices, computed using $2.2\cdot 10^{6}$ surrogate model evaluations in Saltelli’s algorithm, as well as the Borgonovo indices, are shown in Fig. 13.

As expected, the cell-to-cell coupling coefficient is heavily influenced by the iris radius while the equatorial radii have significantly less impact. In particular, the Sobol coefficients indicate that more than $95\%$ of the variance can be attributed to changes in the iris radius $\Delta R_{\text{ir}}$.

To conduct a preliminary inverse analysis, based on the previous findings, we neglect the deformations in the equatorial radii in the following. We note that considering all parameters individually in an inverse study would be numerically much more challenging, as it would require approaches which assure that the problem remains well-posed stuart_2010 . Such a study is considered as out of the scope of the present work. In this case, the task consists in estimating the parameter $\Delta R_{\text{ir}}$ from measurements of the cell-to-cell coupling coefficient $k_{\text{cc}}$. FIG. 14 depicts the associated relation which appears to be monotonic and almost linear in the considered range. We collected measurements of the fundamental mode spectra for $N_{\text{cav}}=826$ cavities from the DESY Database desy_db . Parameter estimation is then carried out by numerically inverting

[TABLE]

where it should be noted that the subscript $i$ refers to the cavity number, while in previous parts of the paper the superscript $i$ referred to the position of a local iris radius variations. This is implemented by reformulating the root-finding problem as an optimization problem and applying the scipy implementation of the L-BFGS-B algorithm byrd1995 ; zhu1997 , i.e., a quasi-Newton method. However, we note that this is not a sensitive choice and other root finding schemes, e.g., bisection or Newton’s method could also be employed. The sample statistics are presented in Table 3.

The XFEL specification requires deviations in the iris radii before welding to be below $0.2\text{\,}\mathrm{m}\mathrm{m}$ with respect to the nominal value. It is expected that the mean values and deviations of the iris radii change during the production chain, for example the chemical treatment is considered to have a significant influence. This offset is estimated by our numerical model to be $\mathbb{E}[\Delta R_{\text{ir},i}]=$0.243\text{,}\mathrm{m}\mathrm{m}, cf. last row of Table [3](#S4.T3). Despite this offset, the standard deviation $\operatorname{Std}[R_{\text{ir},i}]$ from both vendors are still within the limits of the specification (0.17\text{,}\mathrm{m}\mathrm{m}$<$0.2\text{,}\mathrm{m}\mathrm{m}$$) after the full production chain. When considering the vendors separately, then both production processes, standard deviations $0.057\text{\,}\mathrm{m}\mathrm{m}$ and $0.073\text{\,}\mathrm{m}\mathrm{m}$, operate approximately at a three-sigma level.

V Conclusion

In this contribution the manufacturing chain of the EXFEL cavity was summarized and translated into a simulation workflow considering uncertainties. To analyze the sensitivities of the uncertain parameters, we propose an efficient adaptive surrogate modeling technique. The numerical study confirms the expert knowledge that the iris radius is the most critical parameter for the cell-to-cell coupling coefficient. Finally, the surrogate model is used to infer the sample mean and sample standard deviation of the iris radius variations from frequency measurements. For both manufacturers the obtained standard deviations are within the specification.

The practically very relevant determination of statistics for all design parameters from given measurements is still ongoing research, as well as the investigation of global sensitivities for HOM. Incorporating the available measurement data for HOMs

Acknowledgements.

The authors would like to acknowledge the support by the DFG (German Research Foundation) in the framework of the Scientific Network SCHM 3127/1,2 ”Uncertainty quantification techniques and stochastic models for superconducting radio frequency cavities” that provided the basis for this collaborative work. The work of J. Corno, N. Georg and S. Schöps is supported by the Excellence Initiative of the German Federal and State Governments and the Graduate School of Computational Engineering at Technische Universitat Darmstadt. N. Georg’s work is also funded by the DFG grant RO4937/1-1. The work of S. Gorgi Zadeh was/is supported by the German Federal Ministry for Research and Education BMBF under contracts 05H15HRRBA and 05H18HRRB1, respectively. J. Heller’s work was supported by the BMBF under the contracts 05K13HR1 and 05K16HRA. J. Schultz and U. Römer acknowledge the funding by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany´s Excellence Strategy – EXC 2163/1- Sustainable and Energy Efficient Aviation – Project-ID 390881007.

Appendix A Global Sensitivity Analysis

Sensitivity analysis is a powerful tool in order to investigate the model behaviour and analyse the effect of changes in the model inputs Borgonovo2016 . More explicitly, sensitivity analysis determines the key drivers for uncertainty in the model output. In the following $\bm{Y}(\theta)\in\Xi\subset\mathbb{R}^{M}$ denotes the vector of $M$ random input variables and $\mathcal{Q}(\bm{Y}(\theta))\in\Xi_{\mathcal{Q}}\subset\mathbb{R}$ represents the model output, named quantity of interest (QoI) hereafter. Note that in this work the inputs $\bm{Y}$ refer to the geometrical parameters specified in (6) while the output $\mathcal{Q}$ is given by the cell-to-cell coupling coefficient $k_{\text{cc}}$.

Computing derivatives of the QoI is a frequently used approach for sensitivity analysis, however, especially in the case of non-linear models this approach only provides locally information at the nominal value where the derivatives are computed. Therefore, global sensitivity analysis (GSA) methods are developed, that explore the effect of changes in the whole input space Saltelli2008 .

Standard indices used for GSA are Sobol indices saltelli2000sensitivity , which require independent input parameters. Contrary to Sobol indices, Borgonovo indices Borgonovo_2007 do not rely on this assumption. Moreover, they are moment independent quantities. Both indices are explained in more detail in the next two subsections.

A.1 Sobol indices

Sobol indices belong to the class of variance-based sensitivity analysis and are based on the idea of decomposing the variance of the QoI. A detailed description and derivation of the Sobol indices is given in Saltelli2008 and a brief summary, following the exposition in Kadyk2019 , is given hereafter. In a first step, the QoI is decomposed by a so called Hoeffding decomposition as

[TABLE]

where

[TABLE]

Here, $\mathbb{E}_{Y\sim\,i}[\mathcal{Q}|Y_{i}]$ denotes the expected value, conditional on the random variable $Y_{i}$. In the presence of independent input random variables this orthogonal function decomposition can be transformed into a decomposition of the variance as

[TABLE]

Equation (18) can be interpreted as follows: The first sum in (18) describes the individual contribution of the input parameter $Y_{i}$ to the overall variance $\mathbb{V}[\mathcal{Q}]$, while the other terms of the sum contain contributions from combined parameter effects.

Based on this decomposition the first order Sobol indices are defined as

[TABLE]

In addition to the first order indices there are total effect indices $S^{T}_{i}$, which describe the contribution of an individual parameter including higher order (interaction) effects with other parameters to the total variance. As a consequence, by comparing the total indices with the first order indices the modeller can detect potential interaction effects among the input parameters, i.e., if $S^{T}_{i}>S_{i}$.

One approach to compute the Sobol indices is a Monte Carlo based sampling approach introduced in Saltelli2008 and often referred to as Saltelli’s algorithm. Basically, two independent samples of the input vector $\bm{Y}$ are generated. Each sample contains $N$ realizations of each input variable and both samples are stored in a $M\times N$ matrix $\mathbf{A}$ and $\mathbf{B}$, respectively. Additionally a matrix $\mathbf{C}_{i}$ is created, by keeping the i-th row of matrix $\mathbf{A}$ and replacing all other rows by the entries of matrix $\mathbf{B}$. The matrix $\mathbf{C}_{i}$ will be used to approximate the conditional expectation of the QoI, by fixing the variation in one parameter. The first order Sobol index can then be approximated by

[TABLE]

where $\mathbb{E}^{\text{MC}}_{\mathbf{A}},\mathbb{V}^{\text{MC}}_{\mathbf{A}}$ denote the Monte Carlo approximation of the mean value and variance based on the sample stored in $\mathbf{A}$, respectively. Also, $\mathbf{a}_{j}$ and $\mathbf{c}_{j}$ denote one column of the respective sample matrix, which is equal to one independently drawn realization. For further details and the approximation of the total order indices, the reader is referred to Saltelli2008 .

A.2 Borgonovo indices

The Sobol indices presented in the last section are global and quantitative indicators, which do not assume linearity of the model. However, Sobol indices belong to the class of variance based methods that ”implicitly assume that this moment (variance) is sufficient to describe the output variability” saltelli2002sensitivity . In addition, Sobol indices require, that the input variables are statistically independent.

To overcome these limitations Borgonovo Borgonovo_2007 introduced a new moment-independent importance measure, which analyses the effect of each input on the entire output distribution. To this end, the difference between the unconditional pdf $\rho_{\mathcal{Q}}$ and the conditional pdf $\rho_{{\mathcal{Q}}|Y_{i}}$ is computed by the shift function

[TABLE]

The unconditional pdf $\rho_{\mathcal{Q}}$ represents the output distribution, when all input parameter are allowed to vary and the conditional pdf $\rho_{{\mathcal{Q}}|Y_{i}}$ represents the case when the uncertainty in one input is eliminated by fixing it to one value.

Taking the expectation of the shift function w.r.t. $Y_{i}$ and multiplying by a normalization constant yields the Borgonovo sensitivity indices

[TABLE]

In fact it is shown in da2015global that the Borgonovo sensitivity indices are a special case of a global sensitivity analysis framework based on a $f$-divergence dissimilarity measure between the unconditional and conditional distribution of ${\mathcal{Q}}$.

In Borgonovo_2007 several properties of the Borgonovo indices are proven. First of all, the indices are normalized and bounded such that $0\leq\delta_{i}\leq 1$. If $\mathcal{Q}$ is independent on $Y_{i}$ then $\delta_{i}=0$. Since no independence of the input variables is required, the Borgonovo indices can also be applied in the presence of dependent inputs.

Following da2015global , the Borgonovo indices can be computed by reformulating (A.2) in terms of marginal and joint probability densities as

[TABLE]

In this work the joint and marginal densities are estimated via kernel-density estimators based on the given data sample with automated band width selection. The integral in (23) is then estimated via a Monte Carlo approach, where $y^{(j)}_{i},q^{(j)}$ denote a sample point of the corresponding input/output distribution:

[TABLE]

The indices computed by (24) are verified via a second approach, which approximates (A.2) by a histogram-based approach and a good agreement of the two approaches is observed. See Borgonovo_2007 for further information on the histogram-based approach.

Appendix B Surrogate Modeling

As discussed in Section IV, the proposed uncertainty quantification study requires the repeated solution of the parametrized eigenvalue problem (3) for different inputs $\bm{Y}$. In order to achieve a high accuracy, a large number of random (virtual) cavities $\tilde{N}_{\text{cav}}$ is desirable. However, in this case, the computational cost of repeatedly solving (3) with the finite element method would become prohibitive. Hence, we employ surrogate modeling to keep the computational effort feasible. The main idea of polynomial surrogate modeling is to approximate the mapping from input parameters $\bm{y}\in\Xi\subset\mathbb{R}^{M}$, which might, for example, represent a realization of the random vector given in (6), to an output quantity $\mathcal{Q}$, which, e.g., might refer to $f_{1},\,\ldots,f_{9}$ (fundamental mode spectrum) or $k_{\text{cc}}$, using global polynomial approximations. In this case spectral convergence can be expected xiu2010 . In particular, we use approximations of the form

[TABLE]

where $\Psi_{i}$ are global multivariate polynomials and $c_{i}$ denote the associated coefficients. In this work, we rely on the stochastic collocation method nobile2008sparse to compute the approximation (25). In particular, we use the algorithm proposed in narayan2014 , with minor modifications to address the case of multiple quantities of interest. In the following, we recall the main ideas but refer to narayan2014 ; loukrezis2019accessing for the specific details.

The method is based on Leja nodes Leja1957 which are nested interpolation nodes defined by an optimization problem, s.t. for example the uni-variate nodes $\{y_{1}^{(i)}\}_{i}$ are given recursively by

[TABLE]

Since Leja nodes allow, by construction, for a nested and granular refinement, they are well suited for adaptive approximations in the multivariate case narayan2014 ; loukrezis2019adaptive . Note that there are also generalizations of (26), i.e., weighted Leja nodes tailored to the probability density of the input narayan2014 ; loukrezis2019approximation . However, those are not considered here, as we aim for high uniform accuracy.

Corresponding multivariate nodes $\{\bm{y}^{(i)}\}_{i}$ could be obtained by a tensor-product construction of univariate nodes, i.e.

[TABLE]

Stochastic collocation then consists in evaluating the FE model for all $\bm{y}^{(i)}$ and enforcing the corresponding collocation conditions for the surrogate model, i.e.

[TABLE]

For a larger number of parameters (e.g. $M>4$) the computational cost associated to the tensor grid (27) would become prohibitive. Hence, the employed algorithm constructs a sparse-grid, cf. bungartz2004sparse , which utilizes only a subset of (27) and neglects points which do not significantly contribute to the accuracy of the approximation. To this end, the adaptive selection of admissible nodes for refinement is based on the respective point-wise error $\epsilon=|\mathcal{Q}-\mathcal{Q}_{N}|$ if a scalar output quantity is considered. The node associated to the largest error is then chosen to refine the approximation. If several quantities of interest $\mathcal{Q}^{(1)},\mathcal{Q}^{(2)},\ldots$ shall be approximated simultaneously based on the same set of nodes, we replace $\epsilon$ by

[TABLE]

where $w=[w_{1},w_{2},\ldots]$ denotes a vector of weights. The adaptive selection of interpolation nodes based on $\tilde{\epsilon}$ is terminated when a given computational budget is reached. For details on the definition of admissible nodes based on downward-closed multi-index sets and the employed multivariate hierarchical Lagrange polynomials, we refer to narayan2014 ; loukrezis2019accessing .

First, the surrogate modeling technique is applied to construct the 4-variate polynomial surrogate models for step 3 of Section IV. In particular, using a computational budget of 50 FE model evaluations, we construct approximations of the fundamental eigenfrequency $f$ of cell $i$ w.r.t. the variations in both adjacent iris radii $\Delta R_{\text{ir}}^{(i)},\Delta R_{\text{ir}}^{(i+1)}$, the equator radii $\Delta R_{\text{eq}}^{(i)}$ and the length $\Delta L^{(i)}$. The accuracy of the approximations can then be quantified by cross-validation in terms of the empirical $L^{\infty}$ norm. In particular, for a set of random parameter realizations $\{\bm{y}^{\text{cv}}\}_{i=1}^{N^{\text{cv}}}$, we compute

[TABLE]

For the end-cells as well as the middle-cells, the error $E^{\text{cv}}$, computed with a uniformly distributed sample of size $N^{\text{cv}}=1000$, is below $5\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$.

Next, we compute the 19-variate surrogate models employed in step 8 using 500 FE model evaluations. In this case, we consider the 19-dimensional input vector (6) and 10 output quantities $f_{1},\ldots,f_{9},k_{\text{cc}}$ with weight vector $w=[1,\ldots,1,0]$. Cross-validation for the sample $\{\bm{Y}_{\text{sort}}^{(m)}\}_{m=1}^{1000}$, see Section IV, indicates an error below $10\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$ for all fundamental eigenfrequencies and of less than $5\cdot 10^{-4}$ for the cell-to-cell coupling coefficient $k_{\text{cc}}$.

The estimates for expectation as well as standard deviation for each output quantity $\mathcal{Q}^{(i)}$ are then obtained from the respective surrogate models $\mathcal{Q}_{N}^{(i)}$ as

[TABLE]