Developing Uncertainty Quantification Strategies in Electromagnetic Problems Involving Highly Resonant Cavities

The computational simulation credibility process involves assembling and documenting evidence that can be used to ascertain and communicate the believability of predictions that are produced from computational simulations. The development of capabilities and methods for gathering credibility evidence is a core portfolio to build during the development of a new computational simulation functionality. The collection of credibility evidence often has a specific use case in mind and maps activities to related requirements for that use case. An important aspect of the computational simulation credibility process is uncertainty quantification (UQ). In this paper, we focus on discussing uncertainty quantification strategies for the electromagnetic (EM) modeling of high-quality factor resonant cavities. In general, these cavities exhibit sharp resonance peaks [1–8], whose features are largely dependent on the geometrical parameters of the cavity. We will first show this strong dependence for a sample slotted cylindrical cavity in Sec. 2, supporting our analysis of uncertainty quantification presented in this paper.

Electromagnetic modeling of these cavities can, in general, be performed through numerical simulation or analytical formulations. While the former is generally slow and may hinder a full uncertainty quantification analysis where a large number of samples may be required, the latter provides a much faster pathway while still capturing the physical phenomena, enabling quick turnaround for sensitivity analysis and down selection of important parameters, as will be shown in Sec. 3. Finally, an uncertainty quantification analysis for the most important parameters identified in Sec. 3 is presented in Sec. 4. These analyses are pivotal to develop credibility evidence in support of predictions made using computational simulations for highly resonant cavities in electromagnetics.

In Eq. (1), $Ecavity(r)$ is the field inside the cavity at location $r$⁠, and $E0$ is the external incident electric field strength. Note this quantity also largely depends on the frequency of the excitation. In general, shielding effectiveness is plotted in units of dB, which is obtained by computing $20 log10(SE)$ in Eq. (1).

An example of enclosure is represented by the structure shown in Fig. 1(a): an aluminum-alloy cylinder with thickness d and metal conductivity $σ$⁠, with interior height h and interior radius a. These parameters (in addition to the slot) define the resonant frequencies of transverse magnetic (TM) and transverse electric modes [6,8–10] supported by the cavity. A port of entry is introduced as an azimuthal slot on one side of the cylinder located midway along the cylinder height, with a width w and a length $ℓ$ (in what follows, the projected length $ℓp$ used to setup the geometry in the simulation is computed from the formula $ℓ=2a sin−1[ℓp/(2a)]$⁠). Note that the slot acts as a magnetic current source drive for the interior cavity modes from the exterior fields at the frequencies analyzed in this paper below the slot resonance. To probe TM modes, we excite the cavity with an external plane wave source propagating along the −x-direction and with electric field polarized along the z-direction as shown in Fig. 1(a).

The shielding effectiveness properties of the cavity in Fig. 1(a) are shown in Fig. 1(b) around the resonant frequency of the TM₀₁₀ mode at about 1.129 GHz using the analytical unmatched formulation reported in Appendix  B. One can observe very sharp resonance peaks with SE peak values well above 0 dB, signature of very strong fields within the cavity. Furthermore, one can note the strong dependence of SE to the geometrical parameters of the slot: a mere change of 20 mils in slot width caused a change of 8.2 dB in peak SE. This brief analysis shows that small changes in cavity parameters may result in large changes in SE, thus justifying our analysis of uncertainty quantification presented in this paper.

Sensitivity analysis is a tool to identify important inputs and to characterize their relationship with the output [11]. For the use case investigated in this paper, there were two goals of the sensitivity analysis:

1. Reduce the dimensionality of the input space by identifying whether any uncertain inputs have little effect on SE. These inputs can then be held constant in subsequent analyses.

2. Understand which uncertain inputs have the strongest impact on SE, which may then prompt additional studies of those inputs.

An important feature of sensitivity analysis owes to the fact that probability distributions of each parameter are not required to rank parameter sensitivity.

The matched power balance analytical model briefly summarized in Appendix  A was used for sensitivity analysis as it is a computationally efficient representation of the full-wave eiger simulation method, a higher fidelity method of moments code developed at Sandia National Laboratories [12,13]. For this analysis, six uncertain inputs were considered with ranges specified in Table 1. These ranges are not meant to capture physical uncertainties, but rather to highlight how these parameters affect shielding effectiveness; a proper determination of their uncertainties will be performed in future work.

The quantity $Tj$ in Eq. (4) gives the proportion of the variance in SE that can be attributed to the $jth$ input, along with its interactions with other inputs. If $Tj$ is much larger than $Sj$⁠, it indicates that there are significant interactions with the $jth$ input that contribute to the variance in SE. Alternatively, when the difference between these quantities is small, the interactions with the $jth$ input are negligible.

The expectations and variances in Eqs. (3) and (4) are estimated using high-dimensional integrals that often require a Monte Carlo integration approach to solve [11]. The Monte Carlo approach requires many runs of the model to estimate the integrals, and therefore, a surrogate model (i.e., a computationally efficient approximation) is often used in place of the model. In this case, a polynomial chaos expansion (PCE) model was used as a surrogate to the power balance model to estimate the indices [18]. The PCE approach employs bases of multivariate orthogonal polynomials to capture the functional relationship between a response and input random variables. Once the expansion has been constructed, Eqs. (3) and (4) may be evaluated in closed form, yielding estimates of the variance-based indices. dakota implements the generalized PCE scheme, in which the particular polynomials used are based on the distributions of the input variables.

In this instance, the PCE models were constructed via regression. dakota features several compressed sensing techniques for performing regression; here, we employed orthogonal matching pursuit. Additional explanation and implementation details have been reported in Ref. [19]. A training dataset was created from the results of 448 runs of the matched power balance model (see Appendix  A). The 448 points in parameter space were selected using Latin hypercube sampling. Each run produced a prediction of the shielding effectiveness at 51 equally spaced frequencies between 1 and 3 GHz, and a separate PCE model was constructed at each frequency. To reduce overfitting, tenfold cross validation (CV) was used to separately select the order of each expansion. dakota uses the mean squared error for cross validation of PCE models. Total orders between 1 and 4 were explored, and for all frequencies, fourth-order PCEs produced the lowest CV scores, which were all between $3×10−4$ and $7×10−4$⁠.

The resulting variance-based indices over the whole frequency range considered are shown in Fig. 2. There it can be seen that slot width w and slot length $ℓ$ were the highest contributors to the uncertainty in the SE across all frequencies, with slot width contributing 70–80% of the variation in SE. The cavity parameters (i.e., height h, radius a, and conductivity $σ$⁠) and slot depth d contributed a negligible amount, though for low frequencies (i.e., less than 1.5 GHz), slot depth contributed slightly more than the other three parameters. The total-order indices $Tj$ shown as dotted lines in Fig. 2 are similar to the first-order indices $Sj$⁠, especially at lower frequencies, implying that there are likely not significant interactions present. These results were confirmed using the unmatched formulation model in Appendix  B at select resonant modes.

While the results in Fig. 2 allow us to rank the parameters affecting SE, to quantify the SE change from each individual parameter in a more direct way, we plot SE versus frequency while varying each input parameter in Fig. 3. In particular, the plots in Fig. 3 were generated by uniformly sampling each input across its range while holding the other inputs constant at their midpoints. The yellow line, i.e. the second line, represents the mean SE, while the black lines, i.e., the first and third lines, represent the minimum and maximum SE. The width between the black lines gives an estimate of the effect of each input on the SE. This plot confirms the results derived from the Sobol indices in Fig. 2, and both plots help build a more robust credibility evidence.

Since the uncertainty in the cavity parameters was found to contribute a negligible amount to the uncertainty in SE across frequency, it was decided to hold the cavity parameters fixed at nominal values. Namely, the cavity height, radius, and conductivity were fixed at 24 in., 4 in., and $2.6×107 S/m$⁠, respectively, for the remainder of the analyses. This process showed that the analytical power balance models could be used to down select uncertain input parameters for future analyses.

Once the dimensionality of the uncertain input space was reduced as shown in Sec. 3, the next step was to characterize the uncertainty in SE. The goals of this second analysis were again twofold:

1. Estimate the range (i.e., minimum and maximum) of SE based on the range of the remaining three uncertain inputs. This is referred to as an interval uncertainty analysis [20].

2. Compare the results of an analytical code to a full-wave, higher fidelity code to assess whether the analytical code could be used for a full probabilistic uncertainty analysis.

It is important to briefly discuss the distinction between an interval and a probabilistic uncertainty analysis. In the former case, uncertain inputs are bounded by minimum and maximum values, and the objective is to estimate the minimum and maximum SE, as shown in the notional example in Fig. 4. The benefits of this type of analysis are that probability distributions do not need to be defined for the uncertain inputs (a nontrivial task) and that in situations where the input/output relationship is simple (e.g., roughly linear, as is the case in our example), a relatively small number of model runs are required to accurately estimate the range of SE. The primary disadvantage is that statements about the likelihood of different values within that range cannot be made. For example, it would not be possible to estimate the probability that the SE is above, say, 20 dB using this approach. In the latter case, probability distributions must be defined for the uncertain inputs, and these probability distributions are typically sampled in a Monte Carlo procedure and then propagated through the model. The result is a probability distribution on the QoI from which probabilistic statements can be made (see Fig. 5).

The first step was to compare an analytical model against a higher fidelity model for an interval analysis as this can be done with a relatively small number of model runs. We employ here the unmatched formulation reported in Appendix  B that provides SE spectra around resonant modes. The results from this analytical model are compared to those of eiger for the TM₀₁₀ resonant mode.

The interval uncertainty analysis was performed considering slot depth d, slot width w, and projected slot length $ℓp$ with the minimum and maximum values that are defined in Table 2 (cavity radius, height, and metal conductivity were kept fixed at 4 in., 24 in., and $2.6×107 S/m$⁠, respectively). A $43$ full factorial design [21] was used, where the exponent “3” represents the number of uncertain inputs and the base “4” represents the number of levels for each input. The four levels were chosen to be equally spaced within the range of each input as detailed in Table 2. This resulted in 64 model runs comprising of every combination of the input parameters.

This $43$ full factorial design was chosen as the relationship between the inputs and outputs was expected to be roughly linear, but it was of interest to confirm this by assessing whether quadratic or cubic effects were present. It also allowed for the estimation of the range of SE with a relatively limited number of model runs. This design was run for both the unmatched formulation in Appendix  B and eiger, and the results are shown in Fig. 6. Immediately obvious is the good agreement between the unmatched formulation and eiger across the entire input space; this indicates that the unmatched formulation is a good choice for a more computationally expensive probabilistic uncertainty analysis. Additionally, these results are consistent with those seen on the matched bound sensitivity analysis described in Sec. 3. In particular, slot width and length have the most effect on the peak SE, and there do not appear to be significant interactions among different input parameters. This is easily identifiable by the fact that the lines in Fig. 6 are roughly parallel.

Because of the good agreement between the unmatched formulation and eiger results shown in Fig. 6, we decided that the former could be used for a probabilistic uncertainty analysis. Currently, there is not a strong physical basis for the choice of input distributions for the uncertain inputs. However, it was of interest to exercise dakota for uncertainty propagation and to understand how different distributional choices might affect the QoI. Therefore, a probabilistic uncertainty analysis was performed under both uniform and normal distributional assumptions, as defined in Table 3. The mean and standard deviation for the normal distributions were derived from the lower and upper bounds of the uniform distributions; the means are the midpoints, and the standard deviations are 1/6 of the range. This standard deviation was chosen such that approximately 99.7% of the distribution falls within the bounds specified for each input. This choice was somewhat arbitrary—alternatives could be increasing the number of standard deviations within the range or making the variance of the normal distribution equal to that of the uniform. In addition to assessing how different input distributions affect the QoI, we also wished to examine the consequences of reducing the initial set of six input parameters listed in Table 1 to the three most influential ones (the slot dimensions) identified by our sensitivity analysis and listed in Table 2.

To achieve well-converged results, we once again employed PCE surrogates, which can be sampled very inexpensively. Using a procedure similar to the one described in Sec. 3, we created a total of four training sets and models, one for each combination of distribution assumption and variable set (all six input parameters versus three). dakota was used to perform Latin hypercube sampling on unmatched formulation model to obtain training data. For the two six-input cases, the training sets included 112 samples, and for the three-input cases, 64. (For the three-input cases, the three parameters that were not varied were fixed at their means/midpoints.) PCE orders up to five were considered; Table 4 lists the selected order and their cross-validation scores.

The PCE models for the four cases considered were sampled 10⁵ times, and the resulting histograms of SE are shown in Fig. 7. One can see that the two SE distributions (full sets of parameters versus down selected parameters) are virtually the same (besides small differences around the right tail of the distributions), confirming that the slot parameters are the most important parameters for SE determination.

The inset of Fig. 7 illustrates some of the consequences of input distribution selection. While the two distributions are centered at similar values of SE, the distribution is noticeably wider in the uniform case. This observation is consistent with the statistics that are reported in Table 5. For example, the standard deviation of SE for the uniform case is more than 75% greater (5.5 dB versus 3.1 dB) than for the normal case. These results show the importance of selecting input distributions appropriately, particularly when low-probability events in the tails of the distributions matter.

The development of methods and capabilities for use in uncertainty quantification for EM problems is an important aspect of gathering credibility evidence to support predictions derived from computational simulation analyses. The work presented in this paper gives an exemplar for uncertainty quantification and sensitivity analysis in this application space. In particular, we performed a sensitivity analysis for a slotted cylindrical cavity and ranked the contribution to the QoI from each of the input parameters. We observed no interaction among different inputs and determined that the slot parameters are the most important for the evaluation of SE. We then compared the results from the analytical, unmatched bound code to eiger to assess whether the unmatched bound code could accurately bound the SE. Finally, we performed an uncertainty quantification analysis of SE assuming two distributions for the uncertain input parameters and observed a strong dependence on the distribution of the QoI from the distribution of the inputs. This exemplar will be used as a prototype for the application of uncertainty and sensitivity analysis methods to more complex electromagnetic computational simulation analyses in the future.

We acknowledge Dr. Aaron J. Pung, Sandia National Laboratories, for the construction of Fig. 1(a). This paper describes objective technical results and analysis. Any subjective views or opinions that might be expressed in the paper do not necessarily represent the views of the U.S. Department of Energy or the U.S. Government.

Laboratory Directed Research and Development program at Sandia National Laboratories. Sandia National Laboratories is a multimission laboratory managed and operated by National Technology and Engineering Solutions of Sandia, LLC, a wholly owned subsidiary of Honeywell International, Inc., for the U.S. Department of Energy's National Nuclear Security Administration (Contract No. DE-NA-0003525).

a =

interior cavity radius

A =

interior cavity wall surface area

$Am,p,n$ =

coefficient to account for coupling through the slot

CV =

cross validation

d =

cavity thickness and slot depth

$Ez$ =

interior electric field for fundamental TM modes

$E0$ =

external incident electric field strength

$Ecavity(r)$ =

field inside the cavity at location r

EM =

electromagnetic

h =

interior cavity height

$H0$ =

incident magnetic field

$〈|Hj|2〉A$ =

mean squared magnetic field component on the wall

$ℓ$ =

slot length

$ℓp$ =

projected slot length

$Prec$ =

received power of the aperture with backing cavity

$Pwall$ =

absorption in the cavity walls

PCE =

polynomial chaos expansion

QoI =

quantity of interest

r =

location point inside the cavity

$RS$ =

surface resistance

$Sj$ =

first-order sensitivity index

SE =

shielding effectiveness

$SEmax$ =

approximate extreme value of the mean shielding effectiveness

$Tj$ =

total-order sensitivity index

TM =

transverse magnetic

UQ =

uncertainty quantification

w =

slot width

$μ0$ =

absolute permeability of free space

$σ$ =

metal conductivity

$ω$ =

angular frequency