Frequency Spectra Analysis of Space and Time Averaged Quantum Stress Tensor Fluctuations

Observing physical effects of large quantum stress tensor fluctuations requires knowledge of the interactions between the probe and the particles of the underlying quantum fields. The quantum stress tensor operators must first be averaged in time alone or space and time to confer meaningful results, the details of which may correspond to the physical measurement process. We build on prior results to characterize the particle frequencies associated with quantum fluctuations of different magnitudes. For the square of time derivatives of the massless scalar field in a spherical cavity, we find that these frequencies are bounded above in a power law behavior. Our findings provide a way identify the largest quantum fluctuation that may be probed in experiments relying on frequency-dependent interactions.


I. INTRODUCTION
As the search for a full theory of quantum gravity continues, extensions to the semiclassical theory of gravity remain a promising endeavor for understanding quantum gravitational effects. In the semiclassical theory of gravity, the stress-energy-momentum tensor is treated as a quantum operator, related to the classical Einstein tensor by taking an expectation value [1]. This approach has seen success in a number of scenarios, ranging from Hawking radiation from black holes [2] and the resulting gravitational backreaction [3] to quantum particle creation [4]. Extensions to the semiclassical theory include the addition of higher order derivatives in the metric tensor [5], with consequences that include the production of gravitons in the early universe [6,7], and the incorporation of fluctuations of the quantum stress tensor around its mean value [8]. The latter approach motivates the study of quantum fluctuations of stress tensor operators, which have been shown to source a number of possible effects, such as geodesic focusing [9] and imprints on power spectra [10,11].
In recent years, focus has shifted to large vacuum fluctuations due to hints that they may be observable in experiments. A generic property of normal ordered operators quadratic in the fields is the divergence of the higher moments, a problem formally addressed by averaging these operators in time alone or space and time. Physically, we may interpret this averaging as encoding the details of experimental measurements, though the correspondence between the two remains under investigation. In two-dimensional conformal field theory with a Gaussian temporal sampling function, the probability P (x) of measuring a fluctuation of magnitude x is a shifted Gamma distribution [12][13][14] bounded below by the optimal quantum inequality bound [15]. In four dimensions, the situation is qualitatively similar, and studies have been conducted for time averaging [16,17] and spacetime averaging [18].
These results, numerically verified in Refs. [19,20], suggest the probability distribution P (x) asymptotically falls more slowly with x than a decaying exponential function. We will be particularly interested in the conclusions of Refs. [17,18], as they assume smooth, compactly supported sampling functions. These are functions that are strictly zero outside a finite, bounded domain, and may better reflect the reality that physical measurements necessarily take place in a finite spacetime region.
These findings suggest the probabilities of measuring large vacuum fluctuations of stress tensor operators may not be as negligible as one might have expected, spurring work into possible physical effects such as focusing of geodesics induced by spacetime curvature fluctuations [21,22] and enhanced barrier penetration rates of charged particles due to radiation pressure fluctuations [23]. Related phenomena in other systems that emerge from similar quadratic operators are also expected to occur: false vacuum decay of a self-interacting scalar field may be dominated by a pathway [24] different from the usual instanton approximation [25], Rydberg atoms may exhibit velocity fluctuations in response to a sequence of short laser pulses [26], and low-temperature light scattering experiments may find large variations in the number of scattered photons [27].
In this paper, we investigate a different but related question. Models proposing experimental effects inevitably assume certain interactions between the probe and the particles of the quantum fields, and these interactions are often dependent on the particle frequencies.
The relationship between the magnitudes of quantum fluctuations and the angular frequencies of the constituent particles remains unclear, and better understanding will help produce more accurate models for future experiments. In Sec. II, we review a numerical approach for diagonalizing bosonic operators that are quadratic in the fields, a category that includes stress tensor operators. In Sec. III, we build on previous results to investigate the relationship between the frequencies of the particle contents and the magnitudes of the fluctuations.
In Sec. IV, we discuss numerical simulations that confirm the analytical calculations and, in doing so, provide estimates of constants that are not well predicted analytically. In Sec. V, we briefly describe a physical application as an example showing how these results may be of interest in experimental contexts. Finally, in Sec. VI, we summarize our findings and consider some future investigations.
Units in which the reduced Planck constant and the speed of light are set to unity are used throughout this paper. annihilation operators as where Here F kk (t, r) and G kk (t, r) are matrix elements that depend on the specific choice of operator, and 1 is the identity operator. Because physical measurements take place in a finite region of spacetime, a more meaningful quantity to consider is the operator T (t, r) averaged in space and time with compactly supported sampling functions, g(r) and f (t), We further assume the spatial sampling function g(r) is non-negative, spherically symmetric, and real, with unit integral over all space. Likewise, the temporal sampling function f (t) is non-negative, even, and real, with unit integral over all time. Although experimental data is lacking, we expect the details of these functions to correlate with factors in the measurement process. The Fourier transforms of these sampling functions are defined aŝ We may write the analogous expression to Eq. (1) for the spacetime-averaged operatorT as As we are concerned with vacuum fluctuations of these quadratic operators, let us consider the Minkowski vacuum state |0 a , defined as the state where a k |0 a = 0 for all k. We immediately see that the vacuum state is not an eigenstate of the operatorT , so quantum fluctuations will be present in the vacuum state.
To perform numerical simulations of these vacuum fluctuations, we need to characterize the eigenstates of the operatorT . We follow a numerical method developed in Refs. [19,28], briefly summarized here. For the cases we consider, the matrix elements ofF andḠ are real, i.e.,Ḡ kk =Ḡ * kk . A Bogoliubov transformation [29] relates the original set of creation and annihilation operators {a † k , a k } to a new set {b † k , b k }. The two sets of operators are related by the linear transformation The A and B matrices contain the real Bogoliubov coefficients that transform The new set of creation and annihilation operators {b † k , b k } acts on their own respective particle number states Note that the physical particles of interest are the a particles of the Minkowski vacuum state |0 a , not the b particles that emerge as a consequence of the Bogoliubov transformation. That is to say, any observable physics arising from vacuum fluctuations comes from interactions with the real a particles. Substituting Eq. (7) into Eq. (6), we may rewriteT in the diagonal formT where 1 is again the identity operator. Here {λ k } and C shift are constants that depend on the matrix elements ofF andḠ. These constants may be derived by mandating {b † k , b k } obey the usual commutation relations analogous to those in Eq. (2), in addition to requiring the diagonalization conditions leading to Eq. (9). Theoretical considerations from quantum inequalities bound the eigenvalues ofT from below, so λ k ≥ 0 for all k and C shift < 0. We interpret the b particle number states |m b as the eigenstates of the operatorT . The eigenvalues x are calculated by actingT on an eigenstate |m b , Extracting the relevant physics in terms of the physical a particles may be accomplished by transforming back to {a † k , a k } via the Bogoliubov transformation, Eq. (7).

III. DOMINANT MODE CONTRIBUTIONS TO FLUCTUATIONS OF DIFFER-ENT MAGNITUDES
A useful class of operators to consider is the normal ordered squares of time derivatives of the massless scalar field ϕ(t, r), Here τ is a sampling timescale discussed below. The factor of τ 2p +2 ensures T 2p +1 (t, r) remains dimensionless. We introduce a more convenient label p, related to the number of time derivatives p via Thus, p = 3 corresponds to T 3 (t, r) = : τ 4φ2 (t, r):, p = 5 corresponds to T 5 (t, r) = : τ 6φ2 (t, r):, and so on. These operators are of interest because they are related to physical quantities that may be probed experimentally. For example, the p = 3 case can be used to infer the behavior of electromagnetic energy density fluctuations [16], electromagnetic momentum flux fluctuations [23], and even fluid density fluctuations [27].
The high moments of these operators were analyzed in Ref. [17] for time averaging alone and in Ref. [18] for spacetime averaging, which we summarize here. We consider a class of compactly supported sampling functions with Fourier transforms that asymptotically The constants are assumed to be constrained by C f , C g , β, η > 0, 0 < λ ≤ α < 1, and η 1/λ < τ β 1/α . Here τ and are related to how quickly the sampling functions switch on and off near the bounds of their compact supports. In this paper, τ and are on the order of the sampling times and lengths, respectively, and may be viewed as such, though this need not be true in general. The factor of k λ−2 in the expression forĝ(k) is not necessary but is present in our numerical constructions later.
For T p (t, r) as defined in Eq. (11), after spacetime averaging in rectangular coordinates, we find the form in Eq. (6) with Let us define the dth moment of the spacetime-averagedT p as One may show that Here C 1 and C 2 are constants that depend on p and the sampling functions, and k = |k|.
Note that the two limits are distinguished by the absence or presence of the effects of space averaging. The lower moments are independent of the space averaging, but the higher moments depend uponĝ(k). The effect of space averaging is to reduce the rate of growth of the moment as d increases. The value d * depends on the ratio of τ β 1/α to η 1/λ in some power law behavior: the smaller this ratio, the sooner the effects of space averaging emerge and the earlier the transition to the spacetime-averaged limit. We may intuitively understand this behavior by noting that in the limit of no space averaging, → 0, we reduce to d * → ∞. The expression in the worldline limit is derived in Ref. [17] while the expression in the spacetime-averaged limit can be easily generalized from the p = 3 case in Ref. [18] by noting from Eq. (14) that higher p merely includes more factors of angular frequency.
These moments may be related to those of a probability distribution P (x), from which we may infer [17,18] where c = α/p in the worldline limit and c = α/(p − 2) in the spacetime-averaged limit. We will not be concerned with the remaining constants, which are predicted in the worldline limit but are less clear in the spacetime-averaged limit. Similarly to Eq. (16), the transition from the worldline limit to the spacetime-averaged limit occurs around x ≈ x * , depending again on the ratio of τ β 1/α to η 1/λ and may be generalized from Ref. [18] to be We may now see that the transition from worldline to spacetime-averaged behavior is linked both to the ratio τ / and to the relative rates of decay off (ω) and ofĝ(k) functions of their arguments. If λ ≤ α, and β ≈ η, then we see from Eq. (13) thatf (ω) decreases more rapidly with increasing ω than doesĝ(k) with increasing |k| when τ > .
We now wish to estimate the dominant frequency contributions to µ d from Eq. (16) and the dominant eigenvalue contributions from Eq. (17) by finding the peaks of the integrands in these equations. Assuming that d 1, the integrals in Eq. (16) may be approximated as Notice that λ and η do not explicitly appear in the spacetime-averaged limit here despite the presence ofĝ 2 (q + k) in Eq. (16). When d 1, the integrals may be approximated by taking k 1, and in this limit we have k α k λ for α > λ. We may now estimate ω d , the dominant frequency contribution to the dth moment, as the frequency that maximizes the integrands, giving Similarly, we may estimate x d , the dominant eigenvalue contribution to the dth moment, as the eigenvalue that maximizes the integrand in Eq. (17), finding Here we recall that in the worldline limit, c = α/p, whereas in the spacetime-averaged limit, c = α/(p − 2). We combine Eqs. (21) and (22) to find As the value of a is not well predicted, we are not interested in the proportionality constant, which, at any rate, may be inferred numerically. The crucial behavior here is the power law behavior. The transition from the worldline to spacetime-averaged limit is expected to be given by Eq. (19), though we do not have enough data to verify this prediction.
Reference [20] suggests the exponent in Eq. (19) is well predicted but the proportionality constant is not, an observation that likely holds in our case.
Note that the analytical estimates in this section do not yet have a clear physical interpretation. Equation (23)

TUATIONS
In this section we are primarily concerned with the characterization of the frequency modes that contribute to large vacuum fluctuations. We will find for any given quantum fluctuation, a wide range of frequency modes contribute, necessitating the introduction of some measure that encodes this information. In doing so, we will show that Eq. (23) predicts the greatest angular frequency ω d that substantially contributes to a vacuum fluctuation of magnitude x d .
In the numerical simulations, we work in spherical coordinates to take advantage of the assumed spherical symmetry of the spatial sampling function g(r). The spherical coordinate equivalent of Eq. (14) is [20] F nlm,n l m = τ p+1 (ω nl ω n l ) Here n is a positive integer, l is a non-negative integer, and m is an integer satisfying |m| ≤ l. The function j l (kr) is the lth spherical Bessel function, the function P l (y) is the lth Legendre polynomial, δ l,l is the Kronecker delta, R is the radius of the bounding sphere of the quantized system, and where z nl is the nth zero of j l (kr). Reference [20] suggests reliable results may be obtained by setting l = m = 0, giving the simpler form Functionsf (ω) andĝ(k) with the behavior in Eq. (13) may be constructed numerically [17,18]. A particular implementation is detailed in Sec. IVA1 and IVA2 of Ref. [20], which we continue using here. A choice of τ , α, and the sampling time specifiesf (ω) while a choice of , λ, and the sampling length specifiesĝ(r). We consider the specific case α = λ = 0.5 in Eq. (13), which is less susceptible to numerical error. As in Ref. [20], we work in τ = 1 units, and in these units the sampling time and radius are 2 and 0.0028, respectively, while the radius of the bounding sphere, R, is 0.88. In our construction, the parameter is equal to the radius of the sampling volume, i.e., = 0.0028. Note that these parameter choices are suboptimal given the assumption of flat space: the boundary of the quantization sphere is detectable by an observer at the origin, as the light travel time to the boundary and back is shorter than the sampling time. Numerical instabilities prevent a more optimal parameter selection, a topic for future investigation. We may calculate this quantity by evaluating where the Bogoliubov transformation, Eq. (7), has been used. Here m k labels the particle content of the eigenstates |m b with the convention shown in Eq. (8). This result describes the relation between the b particle content, which is encoded in the eigenstates |m b , and the a particle content, which are the physically interesting particles. Figure 1 shows an example of an eigenstate with two b particles, i.e., k m k = 2. There is some indication that the total mean number of a particles may be equal to the total number of b particles, but more analysis is needed to confirm this. The frequency spectrum displays complex structure that eludes straightforward interpretation, though the higher frequency modes have diminishing contributions, a behavior that is generic to the spectra that we will investigate in more detail. While the mean numbers of a particles fall sharply to negligible values at small and high frequencies, we cannot rule out limitations in our simulations, so we leave further discussion of the consequences for future work. Eigenstates with more than two b particles do not display different behavior, so we will consider only the two-particle sector in our simulations. The higher particle sectors are generally less probable, as argued in Ref. [20], and should not significantly affect our conclusions.

B. Characteristic frequencies of frequency spectra
Because the frequency spectra are too complex to analyze directly, we would like to extract a single characteristic frequency from each spectrum. Here we proceed with the simplest characterization, which is merely an average of all frequencies, weighted by the mean number of a particles at different frequencies. The denominator normalizes the total mean number to unity. For our case with l = m = 0 and n = 1 ∼ N , we get We will find it convenient to look at the frequencies scaled by R/π because 1 ≤ ω c R/π ≤ N , where N is the total number of modes in the simulation. In Fig. 2

mode. The anomalous behavior of p = 3 remains under investigation, and it remains unclear
why there is a region (12 < ln x < 13 in the top panel of Fig. 2) where the envelope is not given by this subset of states. Regardless, in this region the upper and lower bounds on the characteristic frequencies almost coincide, so our numerical results should still suffice for rough estimates.
However, given the complex behavior of the spectra, we may be concerned whether a finite-mode computation can be used in physical applications, which are better described by the N → ∞ limit. We first recall our observation from Fig. 1: generically, higher frequency modes have diminishing contributions. We supplement this observation with a simple argument. As higher frequency modes are included, the characteristic frequency ω c , a weighted average of all frequencies, increases. If ω c is not bounded above, then in the N → ∞ limit, we have ω c → ∞ for all quantum fluctuations, a nonsensical result. We thus expect that for any given fluctuation, only some frequency modes have substantial contributions.
We may also draw this conclusion numerically in the following manner. The aforementioned argument suggests that as N increases, ω c will converge to some stable value. We confirm this behavior in Fig. 3 for both ω c and x. These results suggest that our numerical analysis of finite-mode systems can be applied to any physical system of interest, provided we limit our analysis to the region that has converged. Recalling that we can identify the states that give the envelope, we focus our analysis on the envelope of the data sets. In doing so, we will find in Sec. IV C that ω d in Eq. (23) may be interpreted as setting an upper bound on the characteristic frequencies ω c . Note that the end of the converged region may be visually distinguished by the kink in the slopes at large x in Fig. 2.
At this point, it may be worth pointing out some subtleties to our interpretation of the theoretical prediction in Eq. (23). Note that ω c is a weighted average, and thus frequencies greater than ω c still contribute, so it is more precise to say that ω d sets upper bounds on the frequencies that contribute substantially, with the understanding that higher frequency modes may have smaller, and eventually negligible, contributions. We should also note that for any fluctuation of magnitude x, there appears to be exactly one corresponding eigenstate and thus one or two corresponding b-mode frequencies, in the two-particle sector. However, there can be a wide range of a-mode frequencies, which is the focus of our discussion. Thus, Eq. (23) predicts upper bounds on the a-mode frequencies that contribute significantly to fluctuations with magnitudes of order x. converge, a behavior that generically holds for arbitrary parameters and eigenstates.

C. Numerical fits to the envelope
Let us now analyze the numerical data to confirm the theoretical predictions, Eq. (23).
As we are primarily concerned with the power law behavior, we take natural logarithms to find Note that we have neglected a nonzero additive constant that emerges from the proportionality constant in Eq. (23) and our choice to rescale by a factor of R/π. An adequate data set must contain sufficient data in both the worldline and spacetime-averaged limits to perform numerical fits. Given computational constraints, we must make a compromise to get enough data in both regions. Preliminary analysis following the method in Sec. IV B suggests that the convergence rates as functions of N increase as p increases. Consequently, for similar amounts of data points, we must go to higher N for smaller values of p. We further recall that the transition from the worldline to spacetime-averaged limits scales as τ / raised to some power, so smaller values of give us more data in the worldline limit. These two observations underlie our parameter choices for our data sets.
We interpret ω d as the envelope of the characteristic frequencies ω c and x d as the fluctuation magnitudes x. For the three data sets shown in Fig. 2, we perform least squares linear regressions to the envelopes of the converged regions of the worldline and spacetimeaveraged limits, which we may easily determine by eye as the two distinct slopes in each plot. The transition location, Eq. (19), may be analyzed following the method in Sec. IV B of Ref. [20], which we do not consider here due to the lack of data. The linear fits are shown in Fig. 4, and the best fit values and statistical errors are compiled in Table I. We find fairly good consistency between the predictions and numerical results, verifying our interpretation of Eq. (23) as setting upper bounds on the frequencies that substantially contribute to fluctuations of varying magnitudes. Note that full consistency in the spacetime-averaged limit can be obtained by increasing , at the expense of the worldline data. Full consistency in the worldline limit has not been theoretically explored in much depth, and in particular the dependence of the onset of the worldline limit is unknown. Reducing is thus not guaranteed to provide more data in the worldline limit, and a reliable method to procure better data in this limit is a topic for future work. As a final comment, note that in choosing to fit strictly to the converged region, our numerical results are not dependent on N and by extension the higher frequency modes that are omitted from our finite-mode computation.
Our numerical verification of the theoretical prediction, Eq. (23), is thus mode independent and may be applied to systems with higher degrees of freedom.

V. PHYSICAL APPLICATIONS
In most models proposing physical effects arising from large vacuum fluctuations, an interaction between a probe and the particles of quantum field is assumed. The details of these interactions often depend crucially on the frequencies of the quantum particles because the eigenstates of spacetime-averaged quadratic operators are multimode squeezed vacuum states with nonzero particle content, and it is the effects of these particles that are potentially observable. We have discovered that we may characterize the particle frequencies of fluctuations of various magnitudes, which may aid the development of better experimental proposals.
Let us consider a brief example how one may take advantage of the results in prior sections. The case p = 3 is particularly interesting due to its broad applicability, and here we explore which photons arising from energy density fluctuations may be trapped in a metal cavity. Let u be the energy density of the massless scalar field averaged in some region of spacetime. The moments of u are expected to behave similarly to those ofφ 2 , 1 so in our notation the dimensionless measures of the fluctuation magnitudes, x, are related to the actual energy density measurements via Recall that the parameter τ may be considered as an estimate of the sampling duration, so we now reintroduce τ to account for variable measurement times. In Sec. IV C, we found numerical estimates for the exponents in Eq. (23) from the slopes of the data in Fig. 4, assuming τ = 1. The proportionality constant, which is not well predicted theoretically, may also be numerically estimated from the vertical intercepts of the data in those figures.
Accounting for the scaling factor R/π, we find Here we recall that ω d predicts the upper bound on the set of characteristic frequencies ω c that contribute to a fluctuation of magnitude x d ≈ x. We may estimate the transition between the two limits, x * , by the intersection of the worldline and spacetime-averaged limits, from which we deduce This rough estimate of the transition between the two limits gives ln(x * ) ≈ 12, consistent with the data for p = 3 in Fig. 2. We may rewrite Eq. (32) in terms of the ratio u d /u * , finding Here we follow Eq. (31) and define u * = x * τ −4 , and we have used the numerical estimate of x * given by Eq. (34). Here we wish to point out a subtlety regarding our interpretation of τ . In general, our numerical procedure for construction of the temporal sampling function does not require the switch-on time τ and the sampling duration, say t 0 , to be approximately equal. However, the specific functions used in this paper have τ ≈ t 0 . A set of two timescale functions where t 0 τ is possible was discussed in Sec. IID in Ref. [26], where the dependence of bothf (ω) and P (x) upon the ratio t 0 /τ was considered. However, no numerical studies using these more general functions have yet been performed.
Note that the value of t 0 /τ is independent of the choice of units for τ itself. The latter is simply a choice of scale which does not alter the physical description. The meaningful quantities are the dimensionless variables, such as ωτ of x = τ 4 u. We argue that our focus on a subset of temporal sampling functions should not be a source of worry and may in fact be necessary in certain cases. The correspondence between the measurement process and the sampling functions is not always clear, though we expect some relationship between the two.
We now return to the question of reflection in a metal cavity. Such reflection requires the angular frequencies of the photons to be less than the plasma frequency of the metal, ω p .
From Eq. (35), we may then find a constraint on the sampling time given a fluctuation of magnitude u. Requiring ω c < ω p , we find For fixed u/u * , the constraint ω c < ω p thus implies that an increased τ results in a decreased ω c . That is to say, a measurement over longer timescales will observe contributions from smaller frequencies, consistent with our intuition that higher frequency modes are suppressed in these cases. We end our discussion with a specific example. For aluminum, the plasma frequency is ω p ≈ 15 eV [30]. Recalling that in units where the reduced Planck constant is set to unity, 1 eV −1 = 0.66 fs, we find This result may be used to estimate the timescales in which the most of the contributing particle frequencies are small enough to allow reflection from boundaries of an aluminum cavity. Note that because ω d in Eq. (35)is an upper bound on the characteristic frequencies, timescales shorter than the bounds in Eq. (37) may still allow reduced levels of reflection.
Practical application of these results will require an estimate of u/u * and may depend on the particular context of the experiment. Here we merely note that because u and x differ only by a scaling constant, Eq. (31), the ratios u/u * and x/x * are equivalent, so numerical simulations of the sort explored in this paper may be used directly for these estimates.
Although our numerical calculations were performed assuming quantization in a spherical cavity, the analytic arguments given in Sec. III suggest that ω d is independent of the cavity geometry.
A different physical application of the characteristic frequencies involves the effects of vacuum radiation pressure fluctuations on Rydberg atoms [26]. These fluctuations involve an operator with p = 7, as shown in Ref. [23]. The model proposed in Ref. [26] involves the excitation of an atom to a highly excited Rydberg state by a laser pulse, and its subsequent deexcitation by a second pulse. The combined effect can be viewed as a measurement of the p = 7 operator with a two-timescale temporal sampling function, resulting in velocity fluctuations of the atom. Upper bounds on the characteristic frequency ω c associated with a fluctuation of magnitude x is predicted by Eq. (30) to scale as x 1/7 in the worldline limit.
Ideally, we would like to have ω c small compared to the energy level separations of the atom, for which a low upper bound is a sufficient but not necessary condition. Whether this is the case will require a more detailed analysis of involving a two-timescale sampling function, which will be a topic for future research.

VI. CONCLUSION
Phenomena associated with a full theory of quantum gravity may be explored through extensions to the semiclassical theory of gravity. Many such phenomena have been proposed over the years, including geodesic focusing [9,21,22], impacts on power spectra [10,11], increased barrier penetration rates of charged particles [23], alternative false vacuum decay pathways for self-coupled scalar fields [24], velocity fluctuations of Rydberg atoms [26], and large variations in scattered photons in low-temperature light scattering experiments [27].
Many of these effects are not gravitational and may be interesting in their own right.
Thus far, experiments searching for these effects remain lacking. To provide future proposals with the machinery for more accurate calculations, we investigated the particle frequencies driving large vacuum quantum fluctuations of the normal ordered square of time derivatives of the massless scalar field. We focused on the case where these fluctuations are observed in a finite duration of time and finite region of space, which we speculate is more relevant for experimental purposes. We discovered that these quantum fluctuations are composed of a variety of particles at different angular frequencies, and the corresponding frequency spectra display complex behavior that remains to be understood. We extracted a characteristic frequency from each spectra to assign to each quantum state, assumed to be representative of the particle frequencies contributing to the fluctuation associated with that state. In doing so, we found that the characteristic frequencies and fluctuation magnitudes are related by a power law behavior with two distinct limits, the worldline and spacetimeaveraged limits. Furthermore, the numerical results obtained here are in good agreement with the analytic predictions, as illustrated in Table I. These results allow us to describe the particle frequencies that substantially contribute to different fluctuations. We described a simple physical application with photons arising from energy density fluctuations reflecting off the boundary of an aluminum cavity. We showed that the measurement duration for which a majority photons are reflected is bounded below. Such constraints may be used to infer the viability of experiments probing large quantum fluctuations.
Further investigations in this direction may be of interest, some of which we briefly outline here. The qualitatively different behavior of the characteristic frequencies for the p = 3 case remains to be addressed. Perhaps the effects of other frequency modes, such as those where l = 0 and m = 0, are more crucial. Missing in our analysis is also any notion of probabilities: there is no guarantee that the upper bound on the characteristic frequencies is approximately equal to the most likely characteristic frequency, and the latter may also be of theoretical interest. Analysis of the most likely frequency can be done by revising the method in Secs. IV B and IV C to incorporate the probabilities P (x) of measuring a fluctuation of magnitude x. The connection to the probability distribution may run deeper than we considered in this paper. The transition from the worldline to spacetime-averaged limit, Eq. (19), is expected to be identical for both the probability distribution, Eq. (18), and the envelope of the characteristic frequencies, Eq. (23). This equivalence is not guaranteed, and any discrepancies between the probability and frequency data could provide promising avenues to further understand the physics behind large quantum fluctuations.

VII. ACKNOWLEDGMENTS
We would like to thank Christopher Fewster for reading and commenting upon the manuscript. This work was supported by the National Science Foundation under Grant