# Frequency-Temporal Attention Network for Singing Melody Extraction

###### Abstract

Musical audio is generally composed of three physical properties: frequency,
time and magnitude. Interestingly, human auditory periphery also provides
neural codes for each of these dimensions to perceive music. Inspired by these
intrinsic characteristics, a frequency-temporal attention network is proposed
to mimic human auditory for singing melody extraction. In particular, the
proposed model contains frequency-temporal attention modules and a selective
fusion module corresponding to these three physical properties. The frequency
attention module is used to select the same activation frequency bands as did
in cochlear and the temporal attention module is responsible for analyzing
temporal patterns. Finally, the selective fusion module is suggested to
recalibrate magnitudes and fuse the raw information for prediction. In
addition, we propose to use another branch to simultaneously predict the
presence of singing voice melody. The experimental results show that the
proposed model outperforms existing state-of-the-art methods 111Codes are
available in https://github.com/yushuai/FTANet-melodic.

Index Terms—  frequency-temporal attention network, singing melody extraction,
music information retrieval

##  1 Introduction

Singing melody extraction estimates the fundamental frequency (F0) of the
singing melody, which is one challenging and critical task in music
information retrieval (MIR) [1]. Recently it has become an active research
topic with a lot of downstream applications of melody-based AI music, such as
cover song identification [2], query-by-humming [3] and voice separation [4].

With the advances of deep learning technique, several neural network based
methods have been proposed to learn the mapping between the audio and melody
[5, 6, 7, 8]. Bittner et al. [6] proposed a fully convolutional neural network
to learn the pitch salience and achieved promising results. Hsieh et al. [8]
proposed an encoder-decoder architecture and a way to use the bottleneck layer
of the network to estimate the presence of a melody line. Unfortunately, other
type of network without bottleneck layer cannot enjoy the advantage of it. On
the other hand, researchers also attempt to combine the deep learning
technique with human auditory [9, 10, 11]. For example, Gao et al. [11]
proposed to integrate a multi-dilation semantic segmentation model with multi-
resolution representations to respectively simulate top-down and bottom-up
processes of human perception for singing melody extraction. However, prior
works treat the three dimensions (i.e., frequency, time and magnitude) equally
which is inconsistent with human audition. In addition, directly concatenating
the spectral and temporal information may lower the performance of the singing
melody extraction.

Human auditory periphery provides neural codes for frequency, time and
magnitude [12]. When the sound enters the cochlea, different frequencies
within the sound selectively stimulate different regions of the cochlea. Then
human can extract the pitch via temporal patterns generated by unresolved
harmonics in the auditory periphery [13]. In this paper, we focus on designing
a frequency-temporal attention based model for singing melody extraction to
simulate such mechanism.

Based on the Psycho-acoustic research, the human auditory system will select
the stimulated frequency bands in the cochlea. Accordingly, auto-correlation
is performed to capture temporal correlations in the auditory cortex. To mimic
this mechanism, we introduce frequency attention to assign different weights
in the spectrogram along the frequency axis, which corresponds to selecting
the stimulated frequency bands in the cochlea. To simulate the process of
auto-correlation, a temporal attention module is proposed to capture temporal
relations between the adjacent frames, which can capture more complex
nonlinear temporal relationship.

Neuro-physiological research show that some high-frequency signals cannot be
perceived by the spectral mode, but temporal models can weakly perceive them
[12]. We argue that features from spectral and temporal models cannot be
simply concatenated, as the concatenated features may bring noises for melody
extraction which hinders the further improvement of this task. Based on
extending the idea of selective kernel networks used in computer vision [14],
we suggest a selective fusion module to dynamically select the features from
frequency and temporal attention networks. Accordingly, we fuse the features
for thereafter singing melody extraction.

As shown in Fig.1, in diagram (b), the related frequency bands are enhanced
and the unrelated frequency bands are suppressed. Diagram (c) shows the melody
line is more salient via the temporal attention. In diagram (d), the selective
fusion module dynamically selects and fuses features from (b) and (c) to
generate the predicted melody line. We hypothesize that we can learn a
frequency-temporal attention based deep architecture for singing melody
extraction to simulate these potential intrinsic characteristics. Three
technical contributions are made: i) we propose a novel frequency-temporal
attention network to mimic the human auditory assigning different weights in
the time and frequency axis. To the best of our knowledge, there is no such
works for singing melody extraction in the literature. ii) A selective fusion
module is proposed to dynamically assign weight to the spectral features and
temporal features. Then we fuse them for melody extraction. iii) We propose to
use another branch to directly predict the presence of melody.

(a) Input Spectrogram.

(b) Frequency Attention Output.

(c) Temporal Attention Output.

(d) Selective Fusion Output.

Fig. 1: Procedures of singing melody extraction of the 128-256 ms of
‘daisy4.wav’ in the ADC2004 dataset performed by the proposed model.

##  2 Proposed model

The overall architecture of our frequency-temporal attention network is shown
in Fig. 2. It has two branches: the top branch contains stacked convolution
layers performing fast downsampling to obtain high-level semantic
representation for singing melody detection; and the bottom branch contains
the proposed frequency-temporal attention and selective fusion modules. The
Frequency attention and temporal attention, selective fusion module, singing
melody detection in the proposed deep architecture are respectively addressed
in the Section.

###  2.1 Model Input

We choose Combined Frequency and Periodicity (CFP) representation [15] as the
input of our model due to its effectiveness and popularity. The CFP
representation has three parts: the power-scaled spectrogram, generalized
cepstrum (GC) [16, 17] and generalized cepstrum of spectrum (GCoS) [18]. We
use 44100 Hz sampling rate, 2,048-sample window size, and 256-sample hop size
for computing the short-time-Fourier-transformation (STFT).

###  2.2 Frequency-temporal Attention Module

The design of frequency-temporal attention module represents the first
contribution of this work. The detailed architecture is shown in Fig. 2.
Unlike previous works for speech-related tasks, we do not employ classic
attention models [19, 20]. Instead we use 1-D convolution to select relevant
regions (frequency/time) for this task.

Formally, given the input feature map $\mathbf{S}$ $\in\mathbb{R}^{F\times
T\times C^{\prime}}$, an average pooling is firstly applied to the spectrogram
for calculating the distribution of magnitudes along the time axis. Unlike 2-D
average pooling, we use row average pooling to achieve this. The frequency
descriptor $f\in\mathbb{R}^{F\times C^{\prime}}$ can be calculated:

$f=\frac{1}{T}\sum_{i<=T}s_{ij}$ (1)

where $s_{ij}$ is the element in the $i$-th row and $j$-th column in
$\mathbf{S}$.

Fig. 2: The overall architecture of the proposed model. The top branch is
proposed to predict the presence of a melody. The bottom branch consists of
the proposed frequency-temporal attention module and selective fusion module
for F0 estimation. ‘k’ and ‘c’ denotes for the kernel size and the number of
channels, respectively.  Fig. 3: The detailed architecture of the proposed
frequency-temporal attention module. $E_{f}$ and $E_{t}$ are the output of
frequency attention and temporal attention, respectively.

Then 1-D convolution is used to dynamically select task-related frequency for
singing melody extraction. Since 1-D convolution is good at learning the
relationship along the frequency bins or time axis, we choose to employ 1-D
convolution to perform frequency-temporal attention. For a frequency
descriptor $f$, the process of 1-D convolution can be written as:

$V^{l}=f*k^{l}$ (2)

where $k^{l}$ is the kernel of 1-D convolution at $l$-th layer and
$V^{l}\in\mathbb{R}^{F\times C}$ is the newly generated feature map, * stands
for the convolution operation. Finally, we apply a softmax layer to the output
of the 1-D convolution layer $V$ to obtain the frequency attention map
$\mathbf{A_{f}}\in\mathbb{R}^{F\times C}$:

$\mathbf{A_{f}}=softmax(V)$ (3)

Similarly, we can also obtain the temporal attention map
$\mathbf{A_{t}}\in\mathbb{R}^{T\times C}$ through the same process mentioned
above.

Meanwhile, we feed feature map $\mathbf{S}$ into 2-D convolution layers with
the kernel size of $(3\times 3)$ and $(5\times 5)$ to generate two new feature
maps $\\{\mathbf{S_{f}},\mathbf{S_{t}}\\}$$\in\mathbb{R}^{F\times T\times C}$.
Then we perform matrix multiplications between
$\\{\mathbf{S_{f}},\mathbf{S_{t}}\\}$ and
$\\{\mathbf{A_{f}},\mathbf{A_{t}}\\}$ to obtain the final output
$\mathbf{E}=\\{\mathbf{E_{f}},\mathbf{E_{t}}\\}$:

$\begin{split}\mathbf{E_{f}}&=\mathbf{S_{f}}\otimes
broadcast(\mathbf{A_{f}}),\\\ \mathbf{E_{t}}&=\mathbf{S_{t}}\otimes
broadcast(\mathbf{A_{t}}).\end{split}$ (4)

where $broadcast$ is the process of making matrices with different shapes have
compatible shapes for element-wise multiplication.

Fig. 4: The detailed architecture of the proposed selective fusion module.

###  2.3 Selective Fusion Module

The design of selective fusion module (SFM) represents the second contribution
of this work. Inspired by the idea of selective kernel networks used in
computer vision [14], we devise this module to dynamically select spectral and
temporal features and fuse them. The detailed architecture is shown in Fig. 3.
This module takes three inputs: the feature map
$\mathbf{S^{\prime}}\in\mathbb{R}^{F\times T\times C}$ generated from
$\mathbf{S}$ via a $(1\times 1)$ convolution, the frequency attention result
$\mathbf{E_{f}}$ and the temporal attention result $\mathbf{E_{t}}$. Firstly,
an element-wise addition is performed to fuse the three inputs into a new
feature map $\mathbf{\Gamma}$ and then a global average pooling (GAP) is
performed to obtain global descriptor $g\in\mathbb{R}^{C}$.

$g=\frac{1}{F\times T}\sum_{i<=F,j<=T}\Gamma_{ij}$ (5)

After a fully connected (FC) layer for non-linear transformation, three FC
layers are used to learn the importance of each channel of the feature maps.
The softmax layer is applied to obtain the attention map. After obtaining the
attention maps, matrix multiplication is performed between the three inputs
and the attention maps to obtain the weighted feature maps. Finally, we fused
the three weighted feature maps by element-wise addition operation. The fuse
feature map contains rich information selected from the spectral and temporal
modes.

###  2.4 Melody Detection Branch

The design of melody detection branch (MDB) represents the third contribution
of this work. The motivation of this branch is to design a general scheme that
does not depend on special structure and can fast predict the presence of a
melody. In this branch, we employ four stacked convolution layers to directly
perform downsampling on the spectrograms. Specifically, we have the first
three convolution layers with the kernel size $(4\times 1)$ and stride size
$(4\times 1)$ and the last convolution layer with the kernel size $(5\times
1)$ and stride size $(5\times 1)$. As a result, the output of this branch is
$m\in\mathbb{R}^{1\times T}$ stride size of $(5\times 1)$. Following [8], the
output of MDB is concatenated with the salience map to make an $(F+1)\times T$
matrix. By including this branch, the voicing recall (VR) rate can be
improved.

##  3 Experiments

Method | ADC2004 (vocal)
---|---
OA | RPA | RCA | VR | VFA
w/o F-att. | 83.9 | 82.4 | 84.2 | 84.1 | 8.8
w/o T-att. | 83.5 | 82.3 | 84.7 | 85.2 | 8.4
w/o SFM | 84.5 | 83.4 | 85.1 | 85.2 | 8.2
w/o MDB | 84.7 | 84.0 | 85.6 | 87.1 | 11.3
Proposed | 85.9 | 85.2 | 86.7 | 87.5 | 11.1
(a) ADC2004

Method | MIREX 05 (vocal)
---|---
OA | RPA | RCA | VR | VFA
w/o F-att. | 81.7 | 76.5 | 77.7 | 81.0 | 9.2
w/o T-att. | 78.3 | 78.0 | 79.0 | 84.0 | 20.2
w/o SFM | 79.3 | 74.9 | 76.0 | 78.8 | 13.2
w/o MDB | 81.0 | 79.6 | 80.6 | 84.6 | 16.6
Proposed | 84.0 | 80.0 | 80.8 | 85.2 | 9.7
(b) MIREX 05

Table 1: Results of Ablation Study on ADC2004 and MIREX 05 dataset. The values
in the table are percentile. “w/o F-att” and “w/o T-att” denote for without
frequency/temporal attention respectively. “w/o SFM” stands for without
selective fusion module. “w/o MDB” stands for without melody detection branch.

Method | ADC 2004 (vocal)
---|---
OA | RPA | RCA | VR | VFA
MCDNN [5] | 69.8 | 70.2 | 74.3 | 79.0 | 30.7
SegNet [8] | 81.6 | 82.7 | 84.9 | 87.4 | 18.6
MD+MR [11] | 82.8 | 82.4 | 84.6 | 85.3 | 12.7
Proposed | 85.9 | 85.2 | 86.7 | 87.5 | 11.1
(a) ADC 2004

Method | MIREX 05 (vocal)
---|---
OA | RPA | RCA | VR | VFA
MCDNN [5] | 75.6 | 69.4 | 71.1 | 74.9 | 14.7
SegNet [10] | 78.6 | 78.4 | 79.7 | 85.7 | 21.5
MD+MR [11] | 80.7 | 78.6 | 79.5 | 84.3 | 16.9
Proposed | 84.0 | 80.0 | 80.8 | 85.2 | 9.7
(b) MIREX 05

Method | Medley DB (vocal)
---|---
OA | RPA | RCA | VR | VFA
MCDNN [5] | 60.0 | 52.2 | 56.3 | 59.5 | 19.6
SegNet [8] | 65.8 | 59.3 | 64.8 | 72.1 | 27.9
MD+MR [11] | 62.3 | 52.2 | 57.4 | 62.9 | 21.5
Proposed | 66.3 | 58.0 | 63.5 | 69.3 | 23.6
(c) MedleyDB

Table 2: Results of the proposed and baseline methods on ADC 2004, MIREX 05
and Medley DB dataset. The values in the table are percentile.

###  3.1 Experiment Setup

We randomly choose 60 vocal tracks from the Medley DB dataset. To increase the
amount of the training data, we augment the training dataset by copying some
of the chosen vocal tracks. Accordingly, there are 98 clips in the training
dataset. We select only samples having melody sung by human voice from
ADC2004, MIREX 05222https://labrosa.ee.columbia.edu/projects/melody and
MedleyDB for test sets. As a result, 12 clips in ADC2004, 9 clips in MIREX 05
and 12 clips in MedleyDB are selected. Note that there is no overlap between
the training and testing datasets.

To adapt the pitch ranges required in singing melody extraction, following [8]
we set hyperparameters in computing the CFP for our model. For vocal melody
extraction, the number of frequency bins is set to 320, with 60 bins per
octave, and the frequency range is from 31 Hz ($B0$) to 1250 Hz
(${D\kern-0.45003pt\leavevmode\resizebox{}{6.14998pt}{\raisebox{0.0pt}{$\sharp$}}}6$).
We divide the training clips into fixed-length segments of T = 128 frames,
which is 128 milliseconds. Our model is implemented with Keras
333https://keras.io. For model update, we choose the binary cross entropy as
the loss function. The Adam optimizer is used with the learning rate of
0.0001.

Following the convention in the literature, we use the following metrics for
performance evaluation: overall accuracy (OA), raw pitch accuracy (RPA), raw
chroma accuracy (RCA), voicing recall (VR) and voicing false alarm (VFA). We
use mir_eval library [21] with the default setting to calculate the metrics.
For each metric other than VFA, the higher score, the higher performance. In
the literature [1], OA is often considered more important than other metrics.

###  3.2 Ablation Study

To investigate how much the proposed frequency-temporal attention module
contributes to the model, we first remove the frequency attention, and only
temporal attention is used to encode the temporal information. As is shown in
Table 1, the performances on both datasets are decreased. When focusing on OA,
the performance is decreased by 2.4% on ADC2004 and by 2.8% on MIREX 05. We
then remove the temporal attention and only keep the frequency attention. The
performances on both datasets are decreased by 2.9% and by 7.3%, respectively,
on the both datasets.

We then investigate the effectiveness of the selective fusion module. When
focusing on OA, the performances of the ablated version are decreased by 2.8%
and 5.9% on ADC2004 and MIREX 05, respectively. The results justify the
assumption that direct concatenation may hinder the further improvement of the
model. Lastly, we investigate the effectiveness of the proposed melody
detection branch. The results w.r.t. OA on both datasets are decreased by 1.4%
and 3.7%, respectively, on the two datasets. The results demonstrate the
effectiveness of the proposed simple design for melody detection. However,
when focusing on VFA, the ablated version achieves a better score than the
proposed model. We analyzed this phenomenon, and found that the shallow CNNs
cannot capture semantic information well and it is prone to predict a non-
melody frame with complex harmonic patterns as a melody one.

(a) Opera_male3 with ours.

(b) Opera_male3 with SegNet.

(c) Opera_male5 with ours.

(d) Opera_male5 with SegNet.

Fig. 5: Visualization of singing melody extraction results on two opera songs
using different models.

###  3.3 Comparison with existing works

The performances on the three datasets are listed in Table 2. Three baseline
methods are compared in Table 2, including MCDNN [5], SegNet [8] and MD+MR
[11]. We carefully tune the hyper-parameters of the three baseline methods to
ensure that they reach their peak performances on our training dataset. The
proposed model and the three baseline methods are trained on the same dataset.
Compared with the baseline methods, the proposed method achieves the highest
score in general. The results clearly confirm the effectiveness and robustness
of our proposed model. For comparison with other baselines, when focusing on
OA, the proposed method outperforms MCDNN by 23.1% in ADC2004, by 11.1% in
MIREX 05, and by 10.5% in Medley DB. Since all of the three baselines have
melody-detectors, when focusing on VR and VFA, the proposed method achieves
the highest score in ADC2004 and MIREX 05 datasets, and comparable result in
MedleyDB.

To investigate what types of errors are solved by the proposed model, a case
study is performed on two opera songs: “opera_male3.wav” and “opera_male5.wav”
in the ADC2004 dataset. We choose SegNet [8] to compare with due to its
effectiveness and popularity. As depicted in Fig. 5, we can observe that there
are fewer octave errors in diagram (a) and (c) than in diagram (b) and (d).
Moreover, from 1000-1200 ms in the diagram (d), we can find some errors that
predict the wrong frequency bin near the right one, which are correctly
predicted in diagram (c). Through the visualization of the predicted melody
contour, we can say that the performance gains of the proposed model can be
attributed to solving the octave errors and other errors. However, we can also
observe that there seem to be more the melody detection errors (i.e.,
predicting a non-melody frame as a melody one) than in SegNet [8]. The goal of
the melody detection branch is to fast predict the presence of the melody
which does not depend on special structure. We leave this for future research
topic to design more accurate and fast network for improving the quality of
melody detection.

##  4 Conclusion

In this paper, we propose a novel frequency-temporal attention network to
mimic the human auditory for singing melody extraction, which mainly contains
three novel modules: frequency-temporal attention, selective fusion, and
singing melody detection. Frequency-temporal feature learning and singing
melody detection are simultaneously learned in an end-to-end way. Experimental
results show the proposed model outperforms the existing state-of-the-art
models on three datasets. Designing more accurate and faster method to improve
the performance of singing melody detection will be our future work.

##  5 Acknowledgement

This work was supported by National Key R&D Program of China (2019YFC1711800),
NSFC (61671156).

## References

  * [1] Justin Salamon, Emilia Gómez, Daniel PW Ellis, and Gaël Richard,  “Melody extraction from polyphonic music signals: Approaches, applications, and challenges,”  IEEE Signal Processing Magazine, vol. 31, no. 2, pp. 118–134, 2014\.
  * [2] Joan Serra, Emilia Gómez, and Perfecto Herrera,  “Audio cover song identification and similarity: background, approaches, evaluation, and beyond,”  in Advances in Music Information Retrieval, pp. 307–332. Springer, 2010.
  * [3] Chung-Che Wang and Jyh-Shing Roger Jang,  “Improving query-by-singing/humming by combining melody and lyric information,”  IEEE ACM Trans. Audio Speech Lang. Process (TASLP), vol. 23, no. 4, pp. 798–806, 2015.
  * [4] Yukara Ikemiya, Kazuyoshi Yoshii, and Katsutoshi Itoyama,  “Singing voice analysis and editing based on mutually dependent F0 estimation and source separation,”  in Proc. ICASSP, 2015, pp. 574–578.
  * [5] Sangeun Kum, Changheun Oh, and Juhan Nam,  “Melody extraction on vocal segments using multi-column deep neural networks.,”  in Proc. ISMIR, 2016, pp. 819–825.
  * [6] Rachel M Bittner, Brian McFee, Justin Salamon, Peter Li, and Juan Pablo Bello,  “Deep salience representations for f0 estimation in polyphonic music,”  in Proc. ISMIR, 2017, pp. 63–70.
  * [7] Dogac Basaran, Slim Essid, and Geoffroy Peeters,  “Main melody extraction with source-filter NMF and CRNN,”  in Proc. ISMIR, 2018, pp. 82–89.
  * [8] Tsung-Han Hsieh, Li Su, and Yi-Hsuan Yang,  “A streamlined encoder/decoder architecture for melody extraction,”  in Proc. ICASSP, 2019, pp. 156–160.
  * [9] Hsin Chou, Ming-Tso Chen, and Tai-Shih Chi,  “A hybrid neural network based on the duplex model of pitch perception for singing melody extraction,”  in Proc. ICASSP, 2018, pp. 381–385.
  * [10] Ming-Tso Chen, Bo-Jun Li, and Tai-Shih Chi,  “CNN based two-stage multi-resolution end-to-end model for singing melody extraction,”  in Proc. ICASSP, 2019, pp. 1005–1009.
  * [11] Ping Gao, Cheng-You You, and Tai-Shih Chi,  “A multi-dilation and multi-resolution fully convolutional network for singing melody extraction,”  in Proc. ICASSP, 2020, pp. 551–555.
  * [12] William A Yost,  “Pitch perception,”  Attention, Perception, & Psychophysics, vol. 71, no. 8, pp. 1701–1715, 2009.
  * [13] J. F. Schouten, R. J. Ritsma, and B. Lopes Cardozo,  “Pitch of the residue,”  The Journal of the Acoustical Society of America (JASA), vol. 34, no. 9B, pp. 1418–1424, 1962.
  * [14] Xiang Li, Wenhai Wang, Xiaolin Hu, and Jian Yang,  “Selective kernel networks,”  in Proc. CVPR, 2019, pp. 510–519.
  * [15] Li Su and Yi-Hsuan Yang,  “Combining spectral and temporal representations for multipitch estimation of polyphonic music,”  IEEE/ACM Transactions on Audio, Speech and Language Processing (TASLP), vol. 23, no. 10, pp. 1600–1612, 2015.
  * [16] Takao Kobayashi and Satoshi Imai,  “Spectral analysis using generalised cepstrum,”  IEEE Transactions on Acoustics, Speech, and Signal Processing (TASLP), vol. 32, no. 6, pp. 1235–1238, 1984.
  * [17] Keiichi Tokuda, Takao Kobayashi, Takashi Masuko, and Satoshi Imai,  “Mel-generalized cepstral analysis-a unified approach to speech spectral estimation,”  in Proc. ICSLP, 1994.
  * [18] Li Su,  “Between homomorphic signal processing and deep neural networks: Constructing deep algorithms for polyphonic music transcription,”  in Proc. APSIPA ASC, 2017, pp. 884–891.
  * [19] Dzmitry Bahdanau, Kyunghyun Cho, and Yoshua Bengio,  “Neural machine translation by jointly learning to align and translate,”  in Proc. ICLR, 2015.
  * [20] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Lukasz Kaiser, and Illia Polosukhin,  “Attention is all you need,”  in Proc. NeurIPS, 2017, pp. 5998–6008.
  * [21] Colin Raffel, Brian McFee, Eric J Humphrey, Justin Salamon, Oriol Nieto, Dawen Liang, Daniel PW Ellis, and C Colin Raffel,  “mir_eval: A transparent implementation of common mir metrics,”  in Proc. ISMIR, 2014.


# Average Linear and Angular Momentum and Power of Random Fields Near a
Perfectly
Conducting Boundary

Luk R. Arnaut and Gabriele Gradoni

###### Abstract

The effect of a perfectly conducting planar boundary on the average linear
momentum (LM), angular (momentum (AM), and their power for a time-harmonic
statistically isotropic random field is analyzed. These averages are purely
imaginary and their magnitude decreases in a damped oscillatory manner with
distance from the boundary. At discrete quasi-periodic distances and
frequencies, the average LM and AM attain their free-space value. Implications
for the optimal placement or tuning of power and field sensors are analyzed.
Conservation of the flux of the mean LM and AM with respect to the difference
of the average electric and magnetic energies and the radiation stresses via
the Maxwell stress dyadic is demonstrated. The second-order spatial
derivatives of differential radiation stress can be directly linked to the
electromagnetic energy imbalance. Analytical results are supported by Monte
Carlo simulation results. As an application, performance based estimates for
the working volume of a reverberation chamber are obtained. In the context of
multiphysics compatibility, mechanical self-stirred reverberation is proposed
as an exploitation of electromagnetic stress.

##  I Introduction

Angular momentum (AM) in electromagnetics was first theorized by Poynting [1]
and consists of spin (or intrinsic) and orbital contributions. Spin angular
momentum (SAM) is represented by familiar left- and right-hand circular field
polarizations, as demonstrated experimentally by Beth [2]. Recently, it was
found that light waves with a helicoidal phase front exhibit another source of
AM whose magnitude may be several times larger than that of SAM [3]. This
additional part was named orbital angular momentum (OAM) and depends on the
transverse spatial field variation. OAM can be identified with a phase front
rotating helicoidally around the axial direction of energy propagation
$\underline{k}$ that governs linear momentum (LM).

More recently, OAM has gained renewed interest in radio engineering because of
its potential for increasing the number of degrees of freedom and channel
capacity of MIMO wireless communications systems [4]. Its mechanism is the
generation of orthogonal transmission modes with several different azimuthal
phase patterns (integer multiples of $2\pi$) and associated AM power flux,
while maintaining the overall point-to-point LM power flux (Poynting vector).

In the OAM literature, the focus has been almost universally on unbounded
propagation in free space. Real scenarios, however, involve antennas near a
ground plane, scatterers, multipath reflection, etc. This raises the question
to what extent LM and (O)AM can then be preserved or detected. One motivation
for studying OAM near a boundary is the use of reverberation chambers for
evaluating (O)AM properties. Furthermore, OAM may arise unintentionally, e.g.,
from a particular phase distribution of edge currents around apertures or
circuit loops. This may give rise to azimuthal AM power in the near field in a
manner that may not be detectable by wire or loop sensors for conventional LM
power.

LM and AM may also be instrumental in the conversion between electromagnetic
(EM) energy and other forms of physical energy, in what may be termed
multiphysics compatibility (MPC). In particular, MPC includes coupling
phenomena between EM fields and their induced or external forces on material
bodies via radiation stress and shear, i.e., electromagnetic-mechanical
compatibility (EMMC); cf. Sec. IV-B. Another subdomain of MPC is
electromagnetic-thermodynamic compatibility (EMTC), originating from ohmic
losses that are already incorporated in EM constitutive relations. EMMC or
EMTC may be significant in e.g. vacuum chambers, aerospace, applications
involving significant EM forces on small or light objects with large surface-
to-volume ratios (“smart dust”), nanoscale components including VLSI circuits
[5], etc., and in high-power applications where mechanical effects of EM
forces and heating may lead to deformation and arcing in extreme cases. A
basic tenet of MPC is that excitation in one subdomain (e.g., EM) should
preserve functionality and nominal conditions of operation that remain within
tolerance levels also for other MPC subdomains. Here too reverberation
chambers offer a primary test bed, because of their ability to generate high
field strengths and power densities, inducing extreme mechanical and thermal
MPC effects.

In this paper, we analyze the effect of a rigid infinite planar perfect
electrically conducting (PEC) surface on the LM and AM power of an ideal
isotropic random field incident from a half-space. The analysis extends
earlier results for LM and AM in unbounded free space [6] and for energy
density near a PEC plane [7], [8]. The approach differs from traditional
studies of OAM, in which a helicoidal wavefront is considered from the start.
Instead, we analyze to what extent statistical AM may be induced or reduced as
a result of interaction with a PEC boundary. A previously developed
methodology [9] based on local angular spectral plane-wave expansions near an
impedance boundary is applied and extended. Throughout this paper, general
time-dependent EM quantities are shown in roman type; time-harmonic quantities
have a suppressed $\exp({\rm j}\omega t)$ dependence and are denoted in
italics.

##  II Linear and Angular Momentum,
Power and Energy

To establish notions and notations, the definitions of LM and AM for general
time-dependent and harmonic fields, and their connection to EM energy and
power are briefly reviewed.

For spatiotemporal fields $\underline{{\rm E}}(\underline{r},t)$ and
$\underline{{\rm B}}(\underline{r},t)$, the (real) local instantaneous AM
density with reference to a location $\underline{r}_{0}$ is [10, ch. 6], [11,
ch. 1]

$\displaystyle\underline{{\rm
M}}(\underline{r},\underline{r}_{0},t)\stackrel{{\scriptstyle\Delta}}{{=}}(\underline{r}-\underline{r}_{0})\times\underline{{\rm
P}}=\mu_{0}\epsilon_{0}(\underline{r}-\underline{r}_{0})\times\underline{{\rm
S}}$ (1)

where

$\displaystyle\underline{{\rm
P}}(\underline{r},t)\stackrel{{\scriptstyle\Delta}}{{=}}\underline{{\rm
D}}(\underline{r},t)\times\underline{{\rm B}}(\underline{r},t)=\underline{{\rm
S}}(\underline{r},t)/{{\rm c}}^{2}$ (2)

is the LM density and $\underline{{\rm S}}(\underline{r},t)=\underline{{\rm
E}}(\underline{r},t)\times\underline{{\rm H}}(\underline{r},t)$ is the local
LM power flux density (Poynting vector). With
$\underline{\nabla}=(\partial/\partial r)\underline{1}_{r}$,
$\underline{\nabla}\times\underline{{\rm E}}=-\partial\underline{{\rm
B}}/\partial t$, and $\underline{\nabla}\cdot\underline{{\rm D}}=\rho$, (1)
becomes

$\displaystyle{\underline{{\rm M}}(\underline{r},\underline{r}_{0},t)}$
$\displaystyle=$
$\displaystyle\epsilon_{0}(\underline{r}-\underline{r}_{0})\times\int^{t}\left[\underline{\nabla}(\underline{{\rm
E}}\cdot\underline{{\rm E}})-(\underline{\nabla}\cdot\underline{{\rm
E}})\underline{{\rm E}}\right]{\rm d}t$ (3) $\displaystyle=$
$\displaystyle(\underline{r}-\underline{r}_{0})\times\int^{t}\left(\underline{\nabla}{{\rm
U}}_{em}-\rho\underline{{\rm E}}\right){{\rm d}}t$

where ${{\rm U}}_{em}(\underline{r},t)={{\rm U}}_{e}+{{\rm
U}}_{m}\stackrel{{\scriptstyle\Delta}}{{=}}\underline{{\rm
E}}\cdot\underline{{\rm D}}/2+\underline{{\rm B}}\cdot\underline{{\rm H}}/2$.
The electric and magnetic energy densities ${{\rm U}}_{e}$ and ${{\rm U}}_{m}$
are quadratic functions of the EM field that are of purely electric and
magnetic type, whereas $\underline{{\rm S}}$ is of mixed types.

For time-harmonic fields $\underline{E}(\underline{r},\omega)$ and
$\underline{B}(\underline{r},\omega)$, the complex Poynting vector is
$\underline{S}(\underline{r},\omega)\stackrel{{\scriptstyle\Delta}}{{=}}(\underline{E}\times\underline{H}^{*})/2$,
with Re$(\underline{S})=\overline{\underline{{\rm E}}\times\underline{{\rm
H}}}$ representing the time averaged LM power flux density. The corresponding
(complex) AM density is

$\displaystyle\underline{M}(\underline{r},\underline{r}_{0},\omega)$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle(\underline{r}-\underline{r}_{0})\times\underline{P}=\mu_{0}\epsilon_{0}(\underline{r}-\underline{r}_{0})\times\underline{S}$
(4) $\displaystyle=$
$\displaystyle(\underline{r}-\underline{r}_{0})\times\frac{{\rm
j}\epsilon_{0}}{2\omega}\left[\underline{\nabla}|\underline{E}|^{2}-(\underline{\nabla}\cdot\underline{E})\underline{E}^{*}\right].~{}~{}$

For a plane wave, (3) with $\underline{\nabla}=-{\rm
j}\underline{k}\equiv-{\rm j}(\omega/{\rm c})\underline{1}_{k}$ yields

$\displaystyle\underline{M}(\underline{r},\underline{r}_{0},\underline{k})=(\underline{r}-\underline{r}_{0})\times\left(\frac{U_{em}}{{\rm
c}}\underline{1}_{k}-{\rm j}\frac{\rho}{2\omega}\underline{E}^{*}\right)$ (5)

where
$U_{em}(\underline{r},\underline{k})=U_{e}+U_{m}=\epsilon_{0}\underline{E}\cdot\underline{E}^{*}/4+\mu_{0}\underline{H}\cdot\underline{H}^{*}/4$
[11]. It is further assumed that $\underline{r}_{0}=\underline{0}$ and
$\rho\ll k||\underline{D}||$, so that the second term in (5) can be neglected.
In Cartesian coordinates,

$\displaystyle\underline{P}(\underline{r},\omega)=\frac{\mu_{0}\epsilon_{0}}{2}\left[(E_{y}H^{*}_{z}-E_{z}H^{*}_{y})\underline{1}_{x}+(E_{z}H^{*}_{x}\right.$
$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}\left.-E_{x}H^{*}_{z})\underline{1}_{y}+(E_{x}H^{*}_{y}-E_{y}H^{*}_{x})\underline{1}_{z}\right]$
$\displaystyle\underline{M}(\underline{r},\underline{0},\omega)=\frac{\mu_{0}\epsilon_{0}}{2}\times$
$\displaystyle\left\\{\left[y(E_{x}H^{*}_{y}-E_{y}H^{*}_{x})-z(E_{z}H^{*}_{x}-E_{x}H^{*}_{z})\right]\underline{1}_{x}\right.$
$\displaystyle\left.+\left[z(E_{y}H^{*}_{z}-E_{z}H^{*}_{y})-x(E_{x}H^{*}_{y}-E_{y}H^{*}_{x})\right]\underline{1}_{y}\right.$
$\displaystyle\left.+\left[x(E_{z}H^{*}_{x}-E_{x}H^{*}_{z})-y(E_{y}H^{*}_{z}-E_{z}H^{*}_{y})\right]\underline{1}_{z}\right\\}.~{}~{}$
(7)

The total AM can be decomposed into SAM and OAM contributions [10, ch. 7] in
dual-symmetrized form [12] as

$\displaystyle\underline{M}(\underline{r},\underline{r}_{0},\omega)$
$\displaystyle=$
$\displaystyle\frac{1}{4}\left[\epsilon_{0}(\underline{E}\times\underline{A}^{*})+\mu_{0}(\underline{H}\times\underline{F}^{*})\right]$
(8)
$\displaystyle+(\underline{r}-\underline{r}_{0})\times\frac{1}{2}\left[\epsilon_{0}(\underline{\nabla}\,\underline{A})\cdot\underline{E}^{*}\right.$
$\displaystyle\left.+\mu_{0}(\underline{\nabla}\,\underline{F})\cdot\underline{H}^{*}\right]$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle\underline{M}_{s}(\underline{r},\omega)+\underline{M}_{o}(\underline{r},\underline{r}_{0},\omega)$

where
$(\underline{\nabla}\,\underline{A})\cdot\underline{E}^{*}\equiv\underline{\nabla}(\underline{A}\cdot\underline{E}^{*})$,
etc. This identification requires knowledge of the magnetic and electric
vector potentials $\underline{A}(\underline{r}|\underline{r}^{\prime})$ and
$\underline{F}(\underline{r}|\underline{r}^{\prime})$, and hence the spatial
distributions of electric and magnetic source currents
$\underline{J}(\underline{r}^{\prime})$ and
$\underline{K}(\underline{r}^{\prime})$ need to be specified, respectively. In
source-free regions, $\underline{A}$ and $\underline{F}$ in (8) are replaced
with ${\rm j}\underline{E}/\omega$ and ${\rm j}\underline{H}/\omega$,
respectively. For unbounded plane waves, both $\underline{M}_{s}$ and
$\underline{M}_{o}$ are purely imaginary.

##  III Average LM and AM of Random Fields

###  III-A Arbitrary Distance or Frequency

To calculate (II) and (7) explicitly, we employ the angular spectral plane-
wave expansion of random fields [7], [9], [13]

$\displaystyle\underline{E}(\underline{r})=\frac{1}{\Omega}\int\int_{\Omega}\underline{\cal
E}(\Omega)\exp(-{\rm j}\underline{k}\cdot\underline{r}){\rm d}\Omega$ (9)

with a similar expansion for $\underline{H}(\underline{r})$, leading to

$\displaystyle E_{\alpha}(\underline{r})H^{*}_{\beta}(\underline{r})$
$\displaystyle=$
$\displaystyle\frac{1}{\Omega}\left(\int\int_{\Omega}\right)^{2}(\underline{\cal
E}_{1}(\Omega_{1})\cdot\underline{1}_{\alpha})(\underline{\cal
H}^{*}_{2}(\Omega_{2})\cdot\underline{1}_{\beta})$
$\displaystyle\times\exp[-{\rm
j}(\underline{k}_{1}-\underline{k}^{*}_{2})\cdot\underline{r}]\,\delta(\Omega_{1}{\Delta}\Omega_{2}){\rm
d}\Omega_{1}{\rm d}\Omega_{2}$

for $\alpha,\beta\in\\{x,y,z\\}$. Here, $\Omega=2\pi$ sr is the solid angle of
the half space of incidence above the boundary ($z\geq 0$);
$\Omega_{1}\Delta\Omega_{2}\stackrel{{\scriptstyle\Delta}}{{=}}(\Omega_{1}\cup\Omega_{2})\setminus(\Omega_{1}\cap\Omega_{2})$
for $\Omega_{1}$, and $\Omega_{2}$ as point sets; $\delta(\cdot)$ is
Kronecker’s delta ($\delta(\emptyset)=1$, and 0 otherwise); ${\rm
d}\Omega_{i}=\sin\theta_{i}{\rm d}\theta_{i}{\rm d}\phi_{i}$ with elevation
angle $\theta_{i}$ and azimuth angle $\phi_{i}$ in standard spherical
coordinates ($0<\theta_{i}\leq\pi/2$, $0<\phi_{i}\leq 2\pi$).

For each plane-wave component $\\{\underline{{\cal E}}_{i},\underline{{\cal
H}}_{i},\underline{k}_{i}\\}$, a TE/TM decomposition is performed with respect
to its plane of incidence $ok_{i}z$ defining $\phi_{i}=0$ [7], [9]. The
incident plus reflected field is aggregated across the angular spectrum by
integration across $\phi_{i}$, $\theta_{i}$ and the uniformly distributed
polarization angle $-\psi_{i}$ in the locally transverse plane. For example,
for $\alpha=x$ and $\beta=y$, substituting [9, eqs. (10)–(15)] into (III-A)
yields

$\displaystyle E_{x}(z)H^{*}_{y}(z)=\frac{{\rm
j}}{\pi}\int^{\pi/2}_{0}\cos\theta_{1}\sin(2k_{1}z\cos\theta_{1})\sin\theta_{1}{\rm
d}\theta_{1}$ (11) $\displaystyle\times\int^{2\pi}_{0}\left({\cal
E}_{1\theta}{\cal H}^{*}_{2\phi}\cos^{2}\phi_{1}-{\cal E}_{1\phi}{\cal
H}^{*}_{2\theta}\sin^{2}\phi_{1}\right){\rm d}\phi_{1}$

with the locally transverse $\underline{{\cal E}}_{i}$ and $\underline{{\cal
H}}_{i}$ ($i=1,2$) given by [9]

$\displaystyle{\cal E}_{i\phi}={\cal E}_{0}\cos\psi_{i},~{}{\cal
E}_{i\theta}=-{\cal E}_{0}\sin\psi_{i},~{}{\cal E}_{ik}=0$ (12)
$\displaystyle{\cal H}_{i\phi}=\frac{{\cal
E}_{0}}{\eta_{0}}\sin\psi_{i},~{}{\cal H}_{i\theta}=\frac{{\cal
E}_{0}}{\eta_{0}}\cos\psi_{i},~{}{\cal H}_{ik}=0$ (13)

where ${\cal E}_{0}\equiv{\cal E}^{\prime}_{0}-{\rm j}{\cal
E}^{\prime\prime}_{0}$ is the complex amplitude of the circular electric field
of each plane-wave component, with $\langle|{\cal
E}_{0}|^{2}\rangle=2\langle{\cal
E}^{\prime(\prime)^{2}}_{0}\rangle=\langle|\underline{E}|^{2}\rangle/4$ and
$\eta_{0}=\sqrt{\mu_{0}/\epsilon_{0}}$. Integration of (11) followed by
ensemble averaging (denoted as $\langle\cdot\rangle$) over ${\cal E}_{0}$ and
$\psi$ yields

$\displaystyle\langle E_{x}H^{*}_{y}\rangle=-\langle
E_{y}H^{*}_{x}\rangle=-{\rm j}\frac{\langle|{\cal
E}_{0}|^{2}\rangle}{\eta_{0}}j_{1}(2kz)$ (14)

where here and for later use

$\displaystyle j_{0}(2kz)$ $\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle{{\rm sinc}(2kz)}$ (15) $\displaystyle j_{1}(2kz)$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle\frac{\sin(2kz)}{(2kz)^{2}}-\frac{\cos(2kz)}{2kz}=-j^{\prime}_{0}(2kz)$
(16) $\displaystyle j_{2}(2kz)$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle\left(\frac{3}{(2kz)^{3}}-\frac{1}{2kz}\right)\sin(2kz)-\frac{3\cos(2kz)}{(2kz)^{2}}~{}~{}~{}~{}$
(17) $\displaystyle j_{3}(2kz)$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle\left(\frac{15}{(2kz)^{4}}-\frac{6}{(2kz)^{2}}\right)\sin(2kz)$
(18) $\displaystyle-\left(\frac{15}{(2kz)^{3}}-\frac{1}{2kz}\right)\cos(2kz)$

are spherical Bessel functions of the first kind and order zero to three,
respectively. With an analogous calculation,

$\displaystyle\langle E_{x}H^{*}_{z}\rangle=\langle
E_{z}H^{*}_{x}\rangle=\langle E_{y}H^{*}_{z}\rangle=\langle
E_{z}H^{*}_{y}\rangle=0$ (19)

because their kernel’s azimuthal dependence is of the form $\sin\phi$ or
$\cos\phi$, as opposed to their square in (11). Combining (7), (14) and (19)
yields the ensemble averaged LM and AM as

$\displaystyle\langle\underline{P}(\underline{r},\omega)\rangle\equiv\frac{\langle\underline{S}\rangle}{{\rm
c}^{2}}$ $\displaystyle=$ $\displaystyle-{\rm j}\frac{\langle|{\cal
E}_{0}|^{2}\rangle}{{\rm c}^{2}\eta_{0}}j_{1}(2kz)~{}\underline{1}_{z}$ (20)
$\displaystyle\langle\underline{M}(\underline{r},\underline{0},\omega)\rangle$
$\displaystyle=$ $\displaystyle-{\rm j}\frac{\langle|{\cal
E}_{0}|^{2}\rangle}{{\rm
c}^{2}\eta_{0}}j_{1}(2kz)\,\left(y\underline{1}_{x}-x\underline{1}_{y}\right)~{}~{}~{}~{}~{}$
(21)

i.e.,
$\langle\underline{M}\rangle=\underline{r}\times\langle\underline{P}\rangle=-r\langle
P\rangle\underline{1}_{\phi}$. The result for
$\langle\underline{M}(\underline{r},\underline{0},\omega)\rangle$ in (21)
extends to general
$\langle\underline{M}(\underline{r},\underline{r}_{0},\omega)\rangle$ by
subtracting
$\underline{r}_{0}\times\langle\underline{P}(\underline{r},\omega)\rangle$
from $\langle\underline{M}(\underline{r},\underline{0},\omega)\rangle$.

The zero real parts of the power flux densities $\langle{S}\rangle={\rm
c}^{2}\langle{P}\rangle$ and ${\rm c}^{2}\langle{M}\rangle$ indicate zero
time-averaged energy flow in normal (i.e., $\underline{1}_{z}$-directed,
longitudinal) and transverse azimuthal directions, respectively. Thus, the
ensemble averaged incident power and the propagating power after reflection
off a PEC boundary cancel. The combined LM and AM flux densities depend on all
three spatial coordinates via $j_{1}(2kz)$ and the radial transverse distance
$||\underline{r}\times\underline{1}_{z}||=\sqrt{x^{2}+y^{2}}$. At any height
$kz$ above the boundary, $\langle\underline{M}\rangle$ is tangential ($\langle
M_{z}\rangle=0$) and purely solenoidal
($\underline{\nabla}\times\langle\underline{M}\rangle\not=\underline{0}$,
$\underline{\nabla}\cdot\langle\underline{M}\rangle=0$). Similar to
deterministic fields in free space, $\langle\underline{M}\rangle$ is
complementary to the irrotational and normal $\langle\underline{P}\rangle$,
i.e., $\underline{\nabla}\times\langle\underline{P}\rangle=\underline{0}$,
$\underline{\nabla}\cdot\langle\underline{P}\rangle\not=0$. In summary, for
random fields near a planar PEC boundary, the combined linear (longitudinal)
and angular (azimuthal) average power fluxes are reactive and characterized by
a 3-D hybrid vector

$\displaystyle\langle\underline{{\mit\Pi}}(\underline{r},\underline{k})\rangle\stackrel{{\scriptstyle\Delta}}{{=}}-{\rm
j}\frac{\langle|{\cal
E}_{0}|^{2}\rangle}{\eta_{0}}j_{1}(2kz)\left(y\underline{1}_{x}-x\underline{1}_{y}+\underline{1}_{z}\right).$
(22)

Its LM and AM components can be measured, e.g., using a modified magic T or a
turnstile junction with L-shaped in-plane sections that contain inductive
(diaphragm) shunt loads. Alternatively, helical and planar chiral power
sensors with reactive loading may be used to measure the longitudinal and
transverse progression of the phase for SAM and OAM, respectively.

The results (20)–(21) are verified using a Monte Carlo (MC) simulation of the
random plane-wave spectrum (9). We used 370 values of $kz$ ranging from 0.01
to 50 in logarithmic steps of 0.01 and define $x=y=1$ m. For each $kz$, a set
of $n_{\theta}\times n_{\phi}\times n_{\psi}=32\times 16\times 16$ uniformly
spaced angles of incidence and polarization were generated across
$(0,\pi/2]\times(0,2\pi]\times(0,2\pi]$ with $n=30$ complex random fields
${\cal E}_{0}$ for each $\underline{k}$. This yields $245\/760$ plane waves
per $kz$. Fig. 1 shows the resulting average LM power density $\langle
S(kz)\rangle$ and, by extension, the transverse components of the AM power
density ${\rm c}^{2}\langle M(kz)\rangle$ normalized by the values of $x$ and
$y$.

---
Figure 1:  Amplitude of mean linear or angular power flux density components
$\langle S_{z}\rangle={{\rm c}}^{2}\langle P_{z}\rangle={{\rm c}}^{2}\langle
M_{x}\rangle/y=-{{\rm c}}^{2}\langle M_{y}\rangle/x$ [in units W/m2] as a
function of height $kz$ above a PEC boundary for $\langle{{\cal
E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2. Blue solid: theory [eq.
(20)]; blue dashed: envelopes ${\mit\Xi}^{\pm}_{s}$ for $kz\gg 1/2$ [eq. (33)]
and for $kz\ll\sqrt{10}$ [eq. (42)]; black dashed: MC simulation.

###  III-B Optimum Locations of Field vs. Power Sensors

On and at asymptotically large distances from the boundary
$(kz\rightarrow+\infty)$, $\langle\underline{P}\rangle$ and
$\langle\underline{M}\rangle$ vanish, which confirms the results for unbounded
random fields [6, eqs. (19)-(20)], [13, eq. (51)]. They also vanish for
$j_{1}(2kz)=0$, i.e., at frequencies and distances that are related as

$\displaystyle\tan(2kz_{s})=2kz_{s}$ (23)

viz., at $z_{s}=0+$ or $z_{s}\simeq(2m+3)\pi/(4k)$ when $kz_{s}\gg 1$
($m=0,1,2,\ldots$). Such measurement locations are preferential when aiming to
avoid the influence of a PEC boundary on the average LM and AM using a power
sensor. Note that these locations $z_{s}$ apply strictly to dot sensors and CW
excitation, and vice versa. Inevitably, the optima become blurred for a sensor
of finite size or for nonzero bandwidths owing to local averaging.

These findings for $\langle{P}\rangle$ and $\langle{M}\rangle$ complement
those for the average energy densities $\langle U_{e}\rangle$ and $\langle
U_{m}\rangle$ of the 3-D total electric or magnetic vector field (no
subscript), the 2-D tangential (subscript $t$), and the 1-D Cartesian
components ($x,y,z$) [7], [8], written here in an alternative but equivalent
form as

$\displaystyle\langle U(kz)\rangle$ $\displaystyle=$
$\displaystyle\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle\left[1\mp\frac{j_{0}(2kz)}{3}\pm\frac{2j_{2}(2kz)}{3}\right]~{}~{}~{}~{}$
(24) $\displaystyle\langle U_{t}(kz)\rangle$ $\displaystyle=$
$\displaystyle\frac{2\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle}{3}\left[1\mp
j_{0}(2kz)\pm\frac{j_{2}(2kz)}{2}\right]~{}~{}~{}$ (25) $\displaystyle\langle
U_{(x)(y)}(kz)\rangle$ $\displaystyle=$
$\displaystyle\frac{\epsilon_{0}\langle{|\cal E}_{0}|^{2}\rangle}{3}\left[1\mp
j_{0}(2kz)\pm\frac{j_{2}(2kz)}{2}\right]$ (26) $\displaystyle\langle
U_{z}(kz)\rangle$ $\displaystyle=$
$\displaystyle\frac{\epsilon_{0}\langle{|\cal E}_{0}|^{2}\rangle}{3}\left[1\pm
j_{0}(2kz)\pm j_{2}(2kz)\right]$ (27)

where upper and lower signs apply to electric $(U=U_{e})$ and magnetic
$(U=U_{m})$ densities, respectively. Unlike (23) for $\langle{S}\rangle$, the
asymptotic values of $\langle U\rangle$ for $kz\rightarrow\infty$ are now
reached for $j_{0}(2kz)=2j_{2}(2kz)$, i.e., when

$\displaystyle\tan(2kz_{u})=\frac{2kz_{u}}{1-2(kz_{u})^{2}},~{}~{}~{}(kz_{u}\not=0)$
(28)

viz., at $z_{u}\simeq(m+1)\pi/(2k)$ for $kz_{u}\gg 1$. The solutions of (23)
and (28) are separated by $kz_{u}-kz_{s}\simeq\pi/4$ for $m\gg 1$.

Similarly, the frequencies and locations for reaching the asymptotic 2-D
tangential energy $\langle U_{t}(kz\rightarrow\infty)\rangle$ follow from (25)
as solutions of $j_{0}(2kz)=j_{2}(2kz)/2$, i.e., when

$\displaystyle\tan(2kz_{u,t})=\frac{2kz_{u,t}}{1-(2kz_{u,t})^{2}},~{}~{}~{}(kz_{u,t}\not=0)$
(29)

viz., at $z_{u,t}\simeq(m+1)\pi/(2k)$ when $kz_{u,t}\gg 1$. For the 1-D
tangential components $\langle U_{x,y}(kz\rightarrow\infty)\rangle$, the
solutions follow from (26) as being identical to those for $\langle
U_{t}\rangle$. Finally, for the 1-D normal component $\langle
U_{z}(kz\rightarrow\infty)\rangle$, the locations are those for $\langle
S\rangle$ but exclude the boundary plane: from (27), $j_{0}(2kz)=-j_{2}(2kz)$
is satisfied when

$\displaystyle\tan(2kz_{u,z})=2kz_{u,z},~{}~{}~{}(kz_{u,z}\not=0)$ (30)

viz., at $z_{u,z}\simeq(2m+3)\pi/(4k)$ when $kz_{u,z}\gg 1$.

The practical significance of these different optimum values is that,
depending on whether one measures either the Cartesian or vector field energy
density (or intensity) using a field sensor (wire or loop probe) or the
reactive LM or AM power flux using a power sensor (aperture antenna), these
devices should be placed at different heights above a PEC boundary in order to
eliminate the effect of the boundary on the measurement.

The first few optimum locations for $kz_{s}$, $kz_{u}$, $kz_{u,t}\equiv
kz_{u,(x)(y)}$ and $kz_{u,z}$ are listed in Tbl. I and shown in Fig. 2. In
each case, the optimal distances are spaced by asymptotically $\pi/2$, as
follows from the asymptotic approximation [17]

$\displaystyle j_{\ell}(2kz)\simeq\frac{\sin(2kz-\ell\pi/2)}{2kz}~{}~{}~{}{\rm
for}~{}kz\gg\frac{\ell(\ell+1)}{4}.$ (31)

Excluding the boundary and DC regime ($kz\not=0$), the shortest optimum
distance ($m=0$) is attained for a 3-D isotropic field sensor ($u$), followed
by a 1-D ($u_{x}$ or $u_{y}$) or 2-D ($u_{t}$) tangential field probe. A
normal field sensor ($u_{z}$) and power sensor ($s$) must be placed farthest,
at more than twice the distance for the isotropic field probe. If the boundary
plane is included, however, then power sensors exhibit the shortest (viz.,
zero) optimum distance.

$m$ | 0 | 1 | 2 | 3 | 4 | 5
---|---|---|---|---|---|---
$kz_{s}$ | 0 | 2.247 | 3.863 | 5.452 | 7.033 | 8.610
$kz_{u}$ | 1.041 | 2.970 | 4.603 | 6.202 | 7.790 | 9.371
$kz_{u,(t)(x)(y)}$ | 1.372 | 3.058 | 4.658 | 6.242 | 7.821 | 9.398
$kz_{u,z}$ | 2.247 | 3.863 | 5.452 | 7.033 | 8.610 | 10.19
TABLE I:  First six locations $kz_{s}$ and $kz_{u}$ for optimal placement of
power or field sensors with various orientations.
---
Figure 2:  First three optimum locations $kz_{{\rm opt}}(m=0,1,2)$ for a
reactive power sensor (black) vs. various field sensors (red) near a PEC
boundary. Symbols: actual optima (cf. Tbl. I), dashed lines: asymptotic laws.

Note that for field intensity (energy) sensors, the optimal $kz$ values are
irrespective of the electric or magnetic type of the sensor. On the other
hand, for a combined electric-magnetic sensor or an (un)intentional receptor
that is sensitive to different components in different measures and
orientations, its optimal distance or frequency can be estimated on a case by
case basis by decomposing its orientation into normal and tangential
components, which provide the weights for an appropriate superposition of the
components of its average response based on (20) and (24)–(26).

###  III-C Maximum Deviations of Mean Power and Energy

####  III-C1 Envelopes

The results in Sec. III-B listed optimum locations $z_{{\rm opt}}$ for a
specified frequency $k$, or vice versa. In EMC practice, signals are often
wideband or being swept across frequency. In addition, surface imperfections,
proximity of scatterers or other surfaces, etc., may perturb the field and
affect the optimality of these locations. Given these spatio-spectral
uncertainty or aggregation effects, a more distant optimum ($kz_{{\rm opt}}\gg
1$) where the fluctuations of $\langle S(kz)\rangle$ are weaker may actually
be preferential, in order to reduce the sensitivity on the optimum placement.
The envelopes of $\langle S(kz)\rangle$ then represent the maximum deviation
from the asymptotic value $\langle S(kz\rightarrow+\infty)\rangle$ that may be
expected as a suboptimal alternative.

The asymptotic analytic representation $\tilde{\langle
S\rangle}\stackrel{{\scriptstyle\Delta}}{{=}}\langle S\rangle+{\rm j}{\cal
H}[\langle S\rangle]$ follows with the Hilbert transformation ${\cal
H}[j_{1}(2kz)]\simeq\cos(2kz-\pi/2)/(2kz)$ for (31) and $-{\rm
j}\equiv\exp(-{\rm j}\pi/2)$ as

$\displaystyle\langle\widetilde{S}(kz\gg 1/2)\rangle\simeq\frac{\langle|{\cal
E}_{0}|^{2}\rangle}{2\eta_{0}kz}\exp\left[{\rm j}(2kz-\pi)\right]$ (32)

whose signed magnitude defines the upper and lower envelopes

$\displaystyle{\mit\Xi}^{\pm}_{s}(kz\gg 1/2)\simeq\pm\frac{\langle|{\cal
E}_{0}|^{2}\rangle}{2\eta_{0}kz}.$ (33)

These are indicated in Fig. 1. Their half separation
$\Delta{\mit\Xi}_{s}(kz)\stackrel{{\scriptstyle\Delta}}{{=}}[{\mit\Xi}^{+}_{s}(kz)-{\mit\Xi}^{-}_{s}(kz)]/2$
measures the maximum absolute deviation from the asymptotic mean power
$\langle
S(kz\rightarrow+\infty)\rangle=[{\mit\Xi}^{+}_{s}(kz)+{\mit\Xi}^{-}_{s}(kz)]/2=0$.
For both electric and magnetic energy densities, similar definitions and
calculations111In (36), the envelopes of $\langle U_{z}\rangle$ are of leading
second order $(kz)^{-2}$, because of cancelling first-order terms. For the
other energies, additional second- and higher-order terms in $(kz)^{-1}$
provide corrections when $kz\not\gg 3/2$. yield

$\displaystyle{\mit\Xi}^{\pm}_{u}(kz\gg 3/2)$ $\displaystyle=$
$\displaystyle\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle\left(1\pm\frac{1}{2kz}\right)$ (34)
$\displaystyle{\mit\Xi}^{\pm}_{u,(x)(y)}(kz\gg 3/2)$ $\displaystyle=$
$\displaystyle\frac{\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle}{3}\left(1\pm\frac{3}{4kz}\right)$ (35) $\displaystyle=$
$\displaystyle\frac{1}{2}{\mit\Xi}^{\pm}_{u,t}(kz\gg 3/2)$
$\displaystyle{\mit\Xi}^{\pm}_{u,z}(kz\gg 3/2)$ $\displaystyle=$
$\displaystyle\frac{\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle}{3}\left(1\pm\frac{3}{4(kz)^{2}}\right)$ (36)

where upper and lower signs now correspond to upper and lower envelopes. These
are shown in Fig. 3 together with (24)–(27).

---
(a)
(b)
Figure 3:  Theoretical averages [eqs. (24)–(27)] (light solid), their MC
simulations (dark dashed), and asymptotic envelopes ${\mit\Xi}^{\pm}$ [eqs.
(34)–(36)] (light dotted), for (a) electric and (b) magnetic energy densities
[in units J/m3] as a function of height $kz$ above a PEC boundary for
$\langle{{\cal E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2. Black/grey: $\langle
U\rangle$; blue/cyan: $\langle U_{t}\rangle$; maroon/red: $\langle
U_{z}\rangle$; olive/green: $\langle U_{x}\rangle$.

The associated maximum relative deviations are

$\displaystyle\Delta{\xi}_{u}(kz\gg 3/2)$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$
$\displaystyle\frac{\Delta{\mit\Xi}_{u}(kz\gg
3/2)}{{\mit\Xi}_{u}(kz\rightarrow+\infty)}=\frac{1}{2kz}$ (37)
$\displaystyle\Delta{\xi}_{u,(x)(y)}(kz\gg 3/2)$ $\displaystyle=$
$\displaystyle\Delta{\xi}_{u,t}(kz\gg 3/2)=\frac{3}{4kz}~{}~{}~{}$ (38)
$\displaystyle\Delta{\xi}_{u,z}(kz\gg 3/2)$ $\displaystyle=$
$\displaystyle\frac{3}{4(kz)^{2}}$ (39)

For example, at $kz=10$ we find $\langle
U_{x}(10)\rangle=0.3104\,\epsilon_{0}\langle|{\cal E}_{0}|^{2}\rangle$,
whereas ${\mit\Xi}^{\pm}_{u,x}(10)=[(1\pm 0.075)/3]\,\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle$ results in
$\Delta{\mit\Xi}_{u,x}(10)=0.025\,\epsilon_{0}\langle{|\cal
E}_{0}|^{2}\rangle$ and $\Delta{\xi}_{u,x}(10)=0.075$. For $\langle S\rangle$,
a different normalization is adopted than for $\langle U\rangle$, viz.,

$\displaystyle\Delta{\xi}_{s}(kz\gg
1/2)\stackrel{{\scriptstyle\Delta}}{{=}}\frac{\Delta{\mit\Xi}_{s}(kz\gg
1/2)}{\langle|{\cal E}_{0}|^{2}\rangle/\eta_{0}}=\frac{1}{2kz}$ (40)

in view of ${\mit\Xi}_{s}(kz\rightarrow+\infty)=0$. The motivation for this
particular choice will become clear in Sec. III-C3 and corresponds to a
normalization of $\Delta{\mit\Xi}_{u}$ and $\langle U\rangle$ by
$\epsilon_{0}\langle|{\cal E}_{0}|^{2}\rangle$.

####  III-C2 Maximum Deviation

The above relative deviations can be used with a specified minimum distance of
the sensor to the PEC boundary to estimate the maximum local deviation from
the asymptotic mean value, or vice versa. Tbl. II shows these deviations for
$kz=\pi/2$, $\pi$ and $2\pi$. The entries indicate that a normally directed
field sensor gives a smaller maximum deviation than a power sensor, at
sufficiently large distances. For example, $\Delta{\xi}_{u,z}=7.60\%$ compared
to $\Delta{\xi}_{s}=15.9\%$ at $z=\lambda/2$ with the chosen normalization,
whereas at $z=\lambda/4$ the maximum deviations of the normalized $\langle
U_{z}\rangle$ and $\langle S_{z}\rangle$ are more similar ($30.4\%$ vs.
$31.8\%$, respectively).

$kz$ | $\Delta{\xi}_{s}$ | $\Delta{\xi}_{u}$ | $\Delta{\xi}_{u,(t)(x)(y)}$ | $\Delta{\xi}_{u,z}$
---|---|---|---|---
$\pi/2$ | 0.318 | 0.318 | 0.477 | 0.304
$\pi$ | 0.159 | 0.159 | 0.239 | 0.076
$2\pi$ | 0.080 | 0.080 | 0.119 | 0.019
TABLE II:  Maximum relative deviation for envelopes of $\langle
S_{z}\rangle/(\langle|{\cal E}_{0}|^{2}\rangle/\eta_{0})$ and $\langle
U_{\alpha}\rangle/(\epsilon_{0}\langle|{\cal E}_{0}|^{2}\rangle)$ at quarter-,
half-, and full-wavelength distances from a PEC boundary.

At the other extreme of very low frequencies and/or small distances from the
boundary, using the approximation [17]

$\displaystyle
j_{\ell}(2kz)\simeq\frac{(2kz)^{\ell}}{(2\ell+1)!!}~{}~{}~{}{\rm
for}~{}kz\ll\sqrt{\frac{(2\ell+2)(2\ell+3)}{\ell+1}}$ (41)

the asymptote for quasi-static deviations from $\langle S(kz)\rangle$ is

$\displaystyle{\mit\Xi}^{-}_{s}(kz\ll\sqrt{10})\simeq-\frac{2\langle|{\cal
E}_{0}|^{2}\rangle kz}{3\eta_{0}}.$ (42)

####  III-C3 Working Volume of a Reverberation Chamber

The local deviations from the asymptotic $\langle S\rangle$ or $\langle
U\rangle$ may be used to estimate the working volume (WV) of a reverberation
chamber before measurement. Consider a cubic cavity of side length $L_{0}$.
Its WV is a symmetrically located interior cube of sides $L_{0}-2z_{W}$ for
some $z_{W}$. On the boundary of this WV, a maximum tolerable relative
deviation $\Delta{\xi}$ of the mean energy or power is specified with
reference to its asymptotic (ideal) value when $L_{0}\rightarrow+\infty$,
which yields $kz_{W}$ from (37)–(40). The relative WV thus follows as

$\displaystyle\frac{V_{W}(\Delta{\xi})}{V_{0}}=\left(\frac{kL_{0}-2kz_{W}(\Delta{\xi})}{kL_{0}}\right)^{3}.$
(43)

This dependence of ${V_{W}}/{V_{0}}$ on $\Delta{\xi}$ is shown in Fig. 4 for
average power and energy components, for two cavity size values of
$V_{0}/\lambda^{3}$. It is seen that the total (vector) $\langle U_{e}\rangle$
or $\langle U_{m}\rangle$ and $\langle S\rangle\equiv\langle S_{z}\rangle$
yield the same WV. By contrast, $\langle U_{z}\rangle$ gives a larger WV,
whereas $\langle U_{t}\rangle$ and $\langle U_{(x)(y)}\rangle$ give a smaller
WV. E.g., for $\Delta{\xi}=0.1$, the value of $V_{W}/V_{0}$ based on $\langle
U_{e}\rangle$, $\langle U_{m}\rangle$ or $\langle S\rangle$ is $31.7\%$,
whereas for $\langle U_{(e)(m),z}\rangle$ it is $56.3\%$, whilst for $\langle
U_{(e)(m),x}\rangle$ it is merely $14.3\%$. The larger value of WV for
$\langle U_{z}\rangle$ can be traced to the more rapid decay of its envelope,
viz., according to $(kz)^{-2}$, causing the specified uniformity to be reached
closer to the boundary.

In practice, adjacent and opposite walls create additional standing waves,
whence (43) serves merely as a first estimate, whose accuracy nevertheless
increases with chamber volume $V_{0}$. The estimates and results serve as a
guide to a more precise evaluation to account for the effects of all
boundaries, e.g., based on a full-wave simulation or measurement.

Note that, for a combined electric-magnetic field sensor, $V_{W}/V_{0}=1$ for
any specified $\Delta{\mit\Xi}$ because the dependence of $\langle
U_{e,\alpha}\rangle+\langle U_{m,\alpha}\rangle$ on $kz$ cancels for any
component $\alpha$.

---
Figure 4:  Relative working volume $V_{W}$ as a percentage of cavity volume
$V_{0}$ for specified relative deviation $\Delta{\xi}$ of average power ($S$)
or average electric or magnetic energy ($U$). Blue solid:
$V_{0}=(5\lambda)^{3}$; red dashed: $V_{0}=(50\lambda)^{3}$.

##  IV Conservation of Energy and Momentum

###  IV-A Poynting Theorem for Average Power Flux and Energy

For deterministic fields in a lossless medium, the Poynting theorem states
that

$\displaystyle\underline{\nabla}\cdot\underline{S}+{\rm
j}\omega(U_{m}-U_{e})=0.$ (44)

Application of (III-A) yields

$\displaystyle\frac{\partial}{\partial
z}\left\\{\begin{array}[]{l}E_{x}H^{*}_{y}\\\
E_{y}H^{*}_{x}\end{array}\right\\}$ $\displaystyle=$
$\displaystyle\pm\frac{{\rm
j}2k}{3}\delta(\Omega_{1}\Delta\Omega_{2})\left({\cal E}_{1\theta}{\cal
H}^{*}_{2\phi}-{\cal E}_{1\phi}{\cal H}^{*}_{2\theta}\right)$ (47)
$\displaystyle\times\left[j_{0}(2kz)-2j_{2}(2kz)\right]$ (48)

where the upper and lower signs refer to $E_{x}H^{*}_{y}$ and
$E_{y}H^{*}_{x}$, respectively. Substituting (12)–(13) herein followed by
ensemble averaging yields $\partial\langle S_{z}\rangle/\partial z$, which is
found to coincide with the result obtained by direct differentiation of (20).
Thus,
$\langle\underline{\nabla}\cdot\underline{S}\rangle=\underline{\nabla}\cdot\langle\underline{S}\rangle$.
Moreover, by differentiation of (24) using

$\displaystyle\frac{{{\rm d}}j_{n}(2kz)}{{{\rm
d}}z}=\frac{2k}{2n+1}\left[nj_{n-1}(2kz)-(n+1)j_{n+1}(2kz)\right]$ (49)

it is found that (44) also applies to $\langle S\rangle$, viz.,

$\displaystyle\underline{\nabla}\cdot\langle\underline{S}\rangle$
$\displaystyle=$ $\displaystyle-{\rm j}\omega\frac{2\epsilon_{0}\langle|{\cal
E}_{0}|^{2}\rangle}{3}\left[j_{0}(2kz)-2j_{2}(2kz)\right]$ (50)
$\displaystyle=$ $\displaystyle-{\rm j}\omega(\langle U_{m}\rangle-\langle
U_{e}\rangle).$

The result (50) is confirmed by MC simulation shown in Fig. 5. Thus,
Poynting’s theorem for (physical) deterministic power and energy extends to
their (arithmetic) averages for random fields, and also to the average power
flux (divergence). The extended theorem indicates that, for arbitrary $kz$,
the spatial local flow of the reactive (nonradiating) linear and angular
average power fluxes ${\rm c}^{2}\langle{P}\rangle$ and ${\rm
c}^{2}\langle{M}\rangle$ are associated and in phase with temporal
oscillations of the imbalance between electric and magnetic average energies.
For increasing $kz$, the sign of their difference (i.e., the dominance of
local average magnetic over electric energy density, or vice versa)
alternates, while the magnitude (i.e., strength of the imbalance) decreases.

Figure 5:  Conservation of mean linear power flux and rate of change of energy
density $\underline{\nabla}\cdot\langle\underline{S}\rangle\equiv{\rm
j}\partial{{\rm Im}}(\langle{S_{z}}\rangle)/\partial z={\rm j}\omega(\langle
U_{e}\rangle-\langle U_{m}\rangle)$ [in units W/m3] as a function of height
$kz$ above a PEC boundary for $\langle{{\cal
E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2 at $f=100$ MHz. Solid: MC
simulation; dashed: theory.

###  IV-B Conservation of Average LM and AM

####  IV-B1 Average EM Stress

For stochastic fields, the (random) symmetrized Maxwell stress dyadic (cf.
Appendix)

$\displaystyle\underline{\underline{T}}(\underline{r},k)$ $\displaystyle=$
$\displaystyle\underline{\underline{T}}_{e}+\underline{\underline{T}}_{m}\stackrel{{\scriptstyle\Delta}}{{=}}\frac{\epsilon_{0}}{2}\left(\underline{E}\,\underline{E}^{*}-\frac{1}{2}|\underline{E}|^{2}\right)$
(51)
$\displaystyle{+}\frac{\mu_{0}}{2}\left(\underline{H}^{*}\underline{H}-\frac{1}{2}|\underline{H}|^{2}\underline{\underline{I}}\right)$

characterizes the random radiation pressure
$\underline{\underline{T}}\cdot\underline{1}_{z}$ (LM flux) and shear
$\underline{\underline{T}}\times\underline{1}_{z}$ with reference to the
surface normal $\underline{1}_{z}$. The resultant exerted random EM force
follows by integrating $\underline{\underline{T}}$ across the oriented
boundary surface (or $\underline{\nabla}\cdot\underline{\underline{T}}$ across
the enclosed volume). To evaluate $\langle\underline{\underline{T}}\rangle$,
note that

$\displaystyle\langle E_{\alpha}E^{*}_{\beta}\rangle=\langle
H_{\alpha}H^{*}_{\beta}\rangle=0,~{}~{}~{}{\forall}\alpha\not=\beta\in\\{x,y,z\\}$
(52)

because the azimuthal dependence of their kernels are of the form $\sin\phi$,
$\cos\phi$ or $\sin(2\phi)$, or because $\langle{\cal E}_{1\phi}{\cal
E}^{*}_{2\theta}\rangle\propto\langle\sin\psi_{1}\cos\psi_{2}\rangle=0$. From
(51) and (52), with the aid of [8, eqs. (8)–(10)], it follows that
$\langle\underline{\underline{T}}\rangle$ is diagonal, isotropic, and
homogeneous (i.e., nondispersive), viz.,

$\displaystyle\langle\underline{\underline{T}}(kz)\rangle$ $\displaystyle=$
$\displaystyle-\frac{\langle
U_{em,z}\rangle}{2}\underline{\underline{I}}_{t}+\frac{\langle
U_{em,z}\rangle-\langle U_{em,t}\rangle}{2}\underline{1}_{z}\underline{1}_{z}$
(53) $\displaystyle=$ $\displaystyle-\frac{\epsilon_{0}\langle|{\cal
E}_{0}|^{2}\rangle}{3}\underline{\underline{I}}$

where
$\underline{\underline{I}}_{t}\stackrel{{\scriptstyle\Delta}}{{=}}\underline{1}_{x}\underline{1}_{x}+\underline{1}_{y}\underline{1}_{y}$
is the transverse unit dyadic, with

$\displaystyle\langle U_{em,z}\rangle$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$ $\displaystyle\langle
U_{e,z}\rangle+\langle U_{m,z}\rangle=\frac{2}{3}\langle\epsilon_{0}|{\cal
E}_{0}|^{2}\rangle,$ $\displaystyle\langle U_{em,t}\rangle$
$\displaystyle\stackrel{{\scriptstyle\Delta}}{{=}}$ $\displaystyle\langle
U_{e,t}\rangle+\langle U_{m,t}\rangle=\frac{4}{3}\langle\epsilon_{0}|{\cal
E}_{0}|^{2}\rangle,~{}\forall kz\geq 0~{}~{}~{}~{}$ (54)

and $\langle U_{(e)(m),t}\rangle\stackrel{{\scriptstyle\Delta}}{{=}}\langle
U_{(e)(m),x}\rangle+\langle U_{(e)(m),y}\rangle$. The negative sign in (53)
indicates that the average stress constitutes a pressure, as opposed to a
tension. While all $\langle U_{(e)(m),(\alpha)(t)}\rangle$ are individually
dispersive, the sum of each matching pair $\langle U_{e,\alpha}\rangle+\langle
U_{m,\alpha}\rangle$ – and hence $\langle U_{em,t}\rangle$, $\langle
U_{em,z}\rangle$ and $\langle\underline{\underline{T}}\rangle$ – are real and
dispersionless with respect to $kz$. Thus,

$\displaystyle\langle\underline{\underline{T}}\rangle=\langle
T\rangle\underline{\underline{I}}=\left\langle\underline{\underline{\overline{{\rm
T}}}}\right\rangle$ (55)

and

$\displaystyle\underline{\nabla}\cdot\langle\underline{\underline{T}}\rangle=\underline{0},~{}~{}~{}\frac{\partial\langle{T}_{\alpha\alpha}\rangle}{\partial
z}={0},~{}~{}~{}\forall\alpha\in\\{x,y,z\\}.$ (56)

Fig. 6 shows MC results for $\langle\underline{\underline{T}}(kz)\rangle$. The
diagonal elements are easily verified to correspond to their theoretical
constant real values (53). The residual off-diagonal $\langle
T_{\alpha\beta}\rangle$ originate from finite-precision errors in the
numerical quadrature of (52), yet demonstrate that
$\langle\underline{\underline{T}}\rangle$ is Hermitean.

Unlike the overall $\langle\underline{\underline{T}}\rangle$, the average
individual electric $\langle\underline{\underline{T}}_{e}\rangle$ and magnetic
$\langle\underline{\underline{T}}_{m}\rangle$ stress dyadics, i.e.,

$\displaystyle\langle\underline{\underline{T}}_{(e)(m)}(kz)\rangle$
$\displaystyle=$ $\displaystyle-\frac{\langle
U_{(e)(m),z}\rangle}{2}\underline{\underline{I}}_{t}$ (57)
$\displaystyle+\frac{\langle U_{(e)(m),z}\rangle-\langle
U_{(e)(m),t}\rangle}{2}\underline{1}_{z}\underline{1}_{z}~{}~{}~{}$

are dispersive, as follows from (26)–(27). These are also shown in Fig. 6.
From Fig. 3 and (57), it follows that the average normal electric stress
$\langle T_{e,zz}\rangle$ changes from tension ($>0$) to pressure ($<0$) at
$kz\simeq 1.14$. By contrast, the average normal magnetic stress always occurs
as a pressure ($\langle T_{m,zz}\rangle<0$).

---
Figure 6:  MC simulated $\langle\underline{\underline{T}}\rangle$ and
$\langle\underline{\underline{T}}_{(e)(m)}\rangle$ [in units J/m3 or N/m2] as
a function of height $kz$ above a PEC boundary for $\langle{{\cal
E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2. Blue: Re[$\langle
T_{\alpha\beta}\rangle$]; red: Im[$\langle T_{\alpha\beta}\rangle$]; green:
Re[$\langle T_{e,\alpha\alpha}\rangle$]; cyan: Re[$\langle
T_{m,\alpha\alpha}\rangle$]; black: Re[$\langle
T_{e,\alpha\alpha}\rangle-\langle T_{m,\alpha\alpha}\rangle$].

As a practical EMMC application of $\underline{\underline{T}}$ within the
realm of mode-stirred reverberation, consider the concept of electromechanical
self-stirring. In this scenario, LM or AM arises from the EM stress caused by
a source field impinging onto suspended or free-flowing small scatterers,
causing their motion or morphing if inertia is sufficiently small. In turn,
this affects the cavity field distribution, hence
$\underline{\underline{T}}(\underline{r})$ via (51), and therefore
$\underline{P}(\underline{r})$,
$\underline{M}(\underline{r},\underline{r}_{0})$ via (58)–(59), etc. The
result is akin to that of chaff in radar, except that the mechanism for its
dynamics here is purely EM and equally feasible in vacuum. The efficiency and
control of self-stirring is governed by the field strength and mechanical
properties of each scatterer. Contrary to conventional mechanical mode
stirring, self-stirring is most efficient for electrically small scatterers,
as it relies on a net nonzero integrated $\underline{\underline{T}}$ across
their surface. As is well known, mechanical pressure exerted onto the walls of
an overmoded microwave cavity easily results in substantial changes to mode
degeneracy, coupling and spectra, even for geometric distortions smaller than
$\lambda/100$ [14, Fig. 1]. A fortiori, corresponding effects of EM stress are
relevant to macro- or mesoscale structures with small inertia. If left
uncontrolled, self-stirring produces noise additional to thermal noise caused
by ohmic dissipation.

####  IV-B2 Conservation of LM and AM

For deterministic fields in vacuum, the conservation of LM density states that
[10], [11]

$\displaystyle\underline{\nabla}\cdot\underline{\underline{T}}-{\rm
j}\omega\underline{P}=\underline{0}.$ (58)

This follows from $\underline{P}$ in (4) upon adding
$(\underline{\nabla}\cdot\underline{H})\underline{H}^{*}\equiv\underline{0}$,
dual-symmetrizing $\epsilon_{0}\underline{\nabla}|\underline{E}|^{2}$, and
recalling that
$\underline{\nabla}|\underline{E}|^{2}=\underline{\nabla}\cdot(|\underline{E}|^{2}\underline{\underline{I}})$
and
$(\underline{\nabla}\cdot\underline{E})\underline{E}^{*}=\underline{\nabla}\cdot(\underline{E}\,\underline{E}^{*})$.
The corresponding conservation of AM follows by pre-multiplying (58) by
$\underline{r}\times$ and using the dyadic identities
$\underline{r}\times(\underline{\nabla}\cdot\underline{\underline{T}})=-\underline{r}\cdot(\underline{\nabla}\times\underline{\underline{T}})=-\underline{\nabla}\cdot\left(\underline{\underline{T}}\times\underline{r}\right)$
to yield

$\displaystyle\underline{\nabla}\cdot\left(\underline{\underline{T}}\times\underline{r}\right)+{\rm
j}\omega\underline{M}=\underline{0}.$ (59)

For arbitrary $kz$ from the PEC boundary, (59) is in fact a manifestation of
Noether’s theorem with generator $\underline{1}_{z}\times\underline{r}\equiv
r\underline{1}_{\phi}$, which holds because of rotational symmetry around
$\underline{1}_{z}$, such that $M_{z}(kz)=0$ – as already found in (21) – and
hence
$[\underline{\nabla}\cdot\left(\underline{\underline{T}}(kz)\times\underline{r}\right)]_{z}=0$
follows from (59).

Since $\underline{P}$ and $\underline{M}$ are proportional to linear and
angular EM power flux, eqs. (58) and (59) provide a conduit between EM and
mechanical effects. On account of Gauss’s theorem, the rate of change of the
LM inside a finite volume can thus be observed from the net flux of stress
$\underline{\nabla}\cdot\underline{\underline{T}}$ through its boundary where
EM forces act. For random fields, it appears that the average flux
$\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle$ and the
fluxes of the average
$\underline{\nabla}\cdot\langle\underline{\underline{T}}_{(e)(m)}\rangle$
enable alternative physical interpretations of $\langle\underline{P}\rangle$
in (20), as will be shown next.

##### Average Flux

To extend and evaluate (58)–(59) for random fields near a PEC plane, note that
spatial variations of $E_{\alpha}E^{*}_{\beta}$ and $H^{*}_{\alpha}H_{\beta}$
are in the normal direction
($\underline{\nabla}=\underline{1}_{z}(\partial/\partial z)$) for all
$\alpha,\beta\in\\{x,y,z\\}$, whence

$\displaystyle\underline{\nabla}\cdot\underline{\underline{T}}=\sum_{\alpha=z}\sum_{\beta=x,y,z}\frac{\partial{T}_{\alpha\beta}}{\partial\alpha}\underline{1}_{\beta}.~{}~{}~{}$
(60)

Correspondingly, the dyadic skew product
$\underline{\underline{T}}\times\underline{r}$ may be calculated termwise,
i.e.,
$(\underline{E}\,\underline{E}^{*})\times\underline{r}=\underline{E}(\underline{E}^{*}\times\underline{r})$,
etc., followed by differentiation, resulting in

$\displaystyle\underline{r}\times(\underline{\nabla}\cdot\underline{\underline{T}})=\sum_{\alpha=z}\sum_{\beta=x,y,z}(\underline{r}\times\underline{1}_{\beta})\frac{\partial{T}_{\alpha\beta}}{\partial\alpha}.~{}~{}~{}~{}~{}~{}$
(61)

It can be easily shown from (III-A) that $\langle T_{zx}\rangle=\langle
T_{zy}\rangle=0=\langle T_{(e)(m),zx}\rangle=\langle T_{(e)(m),zy}\rangle$ and
that

$\displaystyle\left\langle E_{\alpha}\frac{\partial E^{*}_{\beta}}{\partial
z}\right\rangle=\left\langle H_{\alpha}\frac{\partial H^{*}_{\beta}}{\partial
z}\right\rangle=0,~{}{\forall}\alpha\not=\beta\in\\{x,y,z\\}~{}~{}~{}~{}~{}~{}$
(62)

for the same reasons as those given for (52). Therefore, $\langle{\partial
T_{zx}}/{\partial z}\rangle=\langle{\partial T_{zy}}/{\partial z}\rangle=0$.
Ensemble averaging of (60) and (61) yields

$\displaystyle\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle$
$\displaystyle=$
$\displaystyle\langle\left(\underline{\nabla}\cdot\underline{\underline{T}}\right)_{z}\rangle=\left\langle\frac{\partial
T_{zz}}{\partial z}\right\rangle\underline{1}_{z}$ (63)
$\displaystyle\langle\underline{r}\times(\underline{\nabla}\cdot\underline{\underline{T}})\rangle$
$\displaystyle=$ $\displaystyle y\left\langle\frac{\partial T_{zz}}{\partial
z}\right\rangle\underline{1}_{x}-x\left\langle\frac{\partial T_{zz}}{\partial
z}\right\rangle\underline{1}_{y}.$ (64)

Thus, whereas $\underline{\nabla}\cdot{\underline{\underline{T}}}$ and
$\underline{r}\times(\underline{\nabla}\cdot{\underline{\underline{T}}})$
depend on the $z$-derivative of both the fluctuating EM stress ($T_{zz}$) and
shear $(T_{zx},T_{zy}$) components as generated by each plane wave component,
only the dependence on the EM stress survives after averaging.

Figs. 7 and 8 show
$\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle$ and
$\langle\underline{\nabla}\cdot(\underline{\underline{T}}\times\underline{r})\rangle$
as a function of $kz$, respectively, for a MC simulation based on
$\underline{\underline{T}}$ calculated from (III-A) and after finite
differencing of $\nabla_{z}$ for discrete $kz$. It is seen that conservation
of LM and AM also holds between their statistical averages and the average
flux of EM stresses and their moments (or torque, in case of SAM),
respectively, i.e.,

$\displaystyle\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle-{\rm
j}\omega\langle\underline{P}\rangle$ $\displaystyle=$
$\displaystyle\underline{0}$ (65)
$\displaystyle\langle\underline{\nabla}\cdot\left(\underline{\underline{T}}\times\underline{r}\right)\rangle+{\rm
j}\omega\langle\underline{M}\rangle$ $\displaystyle=$
$\displaystyle\underline{0}.$ (66)

While (65) and (66) follow of course trivially from (58) and (59) from a
purely mathematical perspective, the demonstration of their validity here
through an independent MC calculation of $\underline{\underline{T}}$,
$\langle\underline{\underline{T}}\rangle$,
$\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle$ and
$\langle\underline{\nabla}\cdot(\underline{\underline{T}}\times\underline{r})\rangle$
starting from the plane-wave expansion serves to validate (58) and (59) as a
starting point based on first principles.

---
Figure 7:  Shear ($x,y$) and normal ($z$) components of average spatial
gradient $\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle$ [in
units J/m4 or N/m3] as a function of height $kz$ above a PEC boundary for
$\langle{{\cal E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2 at $f=100$ MHz. Blue solid:
MC real part; red solid: MC imaginary part; blue dashed, theoretical ${\rm
j}\omega\langle{P}_{z}\rangle$.
---
Figure 8:  Angular shear ($x,y$) and stress ($z$) components of average
spatial gradient
$\langle\underline{\nabla}\cdot(\underline{\underline{T}}\times\underline{r})\rangle$
[in units J/m3 or N/m2] as a function of height $kz$ above a PEC boundary for
$\langle{{\cal E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2 at $f=100$ MHz. Blue solid:
MC real part; red solid: MC imaginary part; blue dashed: theoretical ${\rm
j}\omega\langle{M}_{(x)(y)}\rangle$.

##### Flux of Average Electric and Magnetic Stresses

From (20), (56) and (58), it follows that in general

$\displaystyle\underline{\nabla}\cdot\langle\underline{\underline{T}}\rangle-{\rm
j}\omega\langle\underline{P}\rangle$ $\displaystyle\not=$
$\displaystyle\underline{0}$ (67)
$\displaystyle\underline{\nabla}\cdot\langle\underline{\underline{T}}\times\underline{r}\rangle+{\rm
j}\omega\langle\underline{M}\rangle$ $\displaystyle\not=$
$\displaystyle\underline{0}$ (68)

except at the optimum locations $kz_{s}$ where $\langle\underline{P}\rangle$
or $\langle\underline{M}\rangle$ vanish (cf. Sec. III-B). Generally, spatial
differentiation and statistical averaging operations commute (in particular
$\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle=\underline{\nabla}\cdot\langle\underline{\underline{T}}\rangle$)
provided that the probability distribution (of $\underline{\underline{T}}$) is
spatially homogeneous [15, Sec. 7.2]. This condition has previously been shown
to be violated for $U_{e}$ near an EM boundary in its normal direction [16,
Figs. 5-7]. Hence (67) and (68) are consistent with this result. On the other
hand, in Sec. IV-A it was found that
$\langle\underline{\nabla}\cdot\underline{S}\rangle=\underline{\nabla}\cdot\langle\underline{S}\rangle$
even though $S$ and $\langle S\rangle$ (as well as $U_{m}-U_{e}$ and $\langle
U_{m}\rangle-\langle U_{e}\rangle$) are dispersive near the boundary. This
result is again not inconsistent with commutativity, because the boundary zone
field is statistically inhomogeneous. Note, however, that the average total
energy density is nondispersive, viz., $\langle U_{em}\rangle\equiv\langle
U_{e}\rangle+\langle U_{m}\rangle=2\epsilon_{0}\langle|{\cal
E}_{0}|^{2}\rangle$.

Expressions for
$\underline{\nabla}\cdot\langle\underline{\underline{T}}_{(e)(m)}\rangle$ and
$\underline{\nabla}\cdot\langle\underline{\underline{T}}_{(e)(m)}\times\underline{r}\rangle$
are readily obtained as follows. Application of (III-A) to the calculation of
$\partial E_{\alpha}E^{*}_{\beta}/\partial{z}$ and $\partial
H_{\alpha}H^{*}_{\beta}/\partial{z}$ produces zero when $\alpha\not=\beta$,
and

$\displaystyle\frac{\partial}{\partial
z}\left\\{\begin{array}[]{l}{|E_{(x)(y)}|^{2}}\\\
{|H_{(x)(y)}|^{2}}\end{array}\right\\}=\pm\frac{2k}{5}\delta(\Omega_{1}\Delta\Omega_{2})\left[\left\\{\begin{array}[]{l}{\cal
E}_{1\theta}{\cal E}^{*}_{2\theta}\\\ {\cal H}_{1\theta}{\cal
H}^{*}_{2\theta}\end{array}\right\\}\right.$ (73)
$\displaystyle\left.\times\left[3j_{1}(2kz)-2j_{3}(2kz)\right]+\left\\{\begin{array}[]{l}{{\cal
E}_{1\phi}{\cal E}^{*}_{2\phi}}\\\ {{\cal H}_{1\phi}{\cal
H}^{*}_{2\phi}}\end{array}\right\\}5j_{1}(2kz)\right]$ (76)
$\displaystyle\frac{\partial}{\partial
z}\left\\{\begin{array}[]{l}{|E_{z}|^{2}}\\\
{|H_{z}|^{2}}\end{array}\right\\}$ $\displaystyle=$
$\displaystyle\mp\frac{8k}{5}\delta(\Omega_{1}\Delta\Omega_{2})\left\\{\begin{array}[]{l}{{\cal
E}_{1\theta}{\cal E}^{*}_{2\theta}}\\\ {{\cal H}_{1\theta}{\cal
H}^{*}_{2\theta}}\end{array}\right\\}$ (81)
$\displaystyle\times\left[j_{1}(2kz)+j_{3}(2kz)\right]$ (82)

otherwise. Substituting (12)–(13) followed by ensemble averaging leads to

$\displaystyle\frac{\partial\langle U_{(x)(y)}\rangle}{\partial z}$
$\displaystyle=$ $\displaystyle\pm\frac{k\epsilon_{0}\langle|{\cal
E}_{0}|^{2}\rangle}{5}\left[4j_{1}(2kz)-j_{3}(2kz)\right]$ (83)
$\displaystyle\frac{\partial\langle U_{z}\rangle}{\partial z}$
$\displaystyle=$ $\displaystyle\mp\frac{2k\epsilon_{0}\langle|{\cal
E}_{0}|^{2}\rangle}{5}\left[j_{1}(2kz)+j_{3}(2kz)\right]$ (84)
$\displaystyle\frac{\partial\langle\underline{\underline{T}}_{(e)(m)}\rangle}{\partial
z}$ $\displaystyle=$ $\displaystyle\pm\frac{k\epsilon_{0}\langle|{\cal
E}_{0}|^{2}\rangle}{5}\left\\{\left[j_{1}(2kz)+j_{3}(2kz)\right]\underline{\underline{I}}_{t}\right.$
(85)
$\displaystyle\left.-5j_{1}(2kz)\underline{1}_{z}\underline{1}_{z}\right\\}$

where upper and lower signs again refer to electric ($U_{e}$,
$\underline{\underline{T}}_{e}$) and magnetic ($U_{m}$,
$\underline{\underline{T}}_{m}$) quantities, respectively. Eqs. (83)–(85) also
follow from (26)–(27) with the aid of (49). Eq. (85) confirms (56) for the sum
$\langle\underline{\underline{T}}\rangle\equiv\langle\underline{\underline{T}}_{e}\rangle+\langle\underline{\underline{T}}_{m}\rangle$
and enables (67)–(68) to be refined to

$\displaystyle\underline{\nabla}\cdot\langle\underline{\underline{T}}_{(e)(m)}\rangle$
$\displaystyle=$ $\displaystyle\mp{\rm j}\omega\langle\underline{P}\rangle=\mp
k\epsilon_{0}\langle|{\cal E}_{0}|^{2}\rangle
j_{1}(2kz)\underline{1}_{z}~{}~{}~{}~{}~{}$ (86)
$\displaystyle\underline{\nabla}\cdot\langle\underline{\underline{T}}_{(e)(m)}\times\underline{r}\rangle$
$\displaystyle=$ $\displaystyle\pm{\rm j}\omega\langle\underline{M}\rangle.$
(87)

Furthermore, (86) and (87) reveal a conservation law between the flux of the
differential average EM stress
$\Delta\langle\underline{\underline{T}}\rangle=\langle\underline{\underline{T}}_{e}\rangle-\langle\underline{\underline{T}}_{m}\rangle$,
and the rate of oscillation of the average LM or AM, viz.,

$\displaystyle\underline{\nabla}\cdot(\langle\underline{\underline{T}}_{e}\rangle-\langle\underline{\underline{T}}_{m}\rangle)+{\rm
j}2\omega\langle\underline{P}\rangle$ $\displaystyle=$
$\displaystyle\underline{0}$ (88)
$\displaystyle\underline{\nabla}\cdot\left[(\langle\underline{\underline{T}}_{e}\rangle-\langle\underline{\underline{T}}_{m}\rangle)\times\underline{r}\right]-{\rm
j}2\omega\langle\underline{M}\rangle$ $\displaystyle=$
$\displaystyle\underline{0}.$ (89)

Eq. (88) also follows more directly by substituting (25) and (27) into (57),
which yields $\Delta{T}_{zz}=j_{0}(2kz)$, and applying (49). To illustrate and
validate these results, Figs. 6 and 9222The transverse diagonal elements
($\alpha=x,y$; cf. top two plots in Fig. 9) are only relevant in
configurations for which ${\nabla}_{x}$ and/or ${\nabla}_{y}$ are also
nonzero, e.g., when one or more additional adjacent boundaries are
perpendicular to $\underline{1}_{x}$ or $\underline{1}_{y}$.  show
$\langle\Delta
T_{\alpha\alpha}(kz)\rangle\stackrel{{\scriptstyle\Delta}}{{=}}\langle{T}_{e,\alpha\alpha}(kz)\rangle-\langle{T}_{m,\alpha\alpha}(kz)\rangle$
and their divergence, respectively. The latter Figure confirms (88) for
$\alpha=z$ and (85) for arbitrary $\alpha$.

---
Figure 9:  Spatial gradients of average differential EM stress components
$\underline{\nabla}\cdot\langle\Delta\underline{\underline{T}}\rangle$ [in
units J/m4 or N/m3] as a function of height $kz$ above a PEC boundary and
$\langle{{\cal E}^{\prime}_{0}}^{2}\rangle=\langle{{\cal
E}^{\prime\prime}_{0}}^{2}\rangle=10^{-4}$ (V/m)2 at $f=100$ MHz. Blue solid:
MC real part, red solid: MC imaginary part, black dashed: theoretical
$\partial\langle{T}_{e,\alpha\alpha}{-}{T}_{m,\alpha\alpha}\rangle/\partial z$
for $\alpha=x,y,z$ [eq. (85)].

###  IV-C Relation Between EM Stress and Energy

For deterministic fields, substitution of (58) into (44) enables the second-
order spatial derivatives (curvatures) of the radiation stress components to
be related to the EM energy imbalance $U_{m}-U_{e}$ (or the second-order time
derivative of the total energy, for general time-dependent fields) without
recourse to the LM or power, viz.,

$\displaystyle{\rm
j}\omega\underline{\nabla}\cdot\underline{P}=k^{2}(U_{m}-U_{e})=\underline{\nabla}\cdot(\underline{\nabla}\cdot\underline{\underline{T}})$
(90)

where
$\underline{\nabla}\cdot(\underline{\nabla}\cdot\underline{\underline{T}})=\partial^{2}T_{zz}/\partial
z^{2}$ in the present configuration. For an averaged random field, the
corresponding identity follows from (50) and (65) as

$\displaystyle\underline{\nabla}\cdot\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle-k^{2}(\langle
U_{m}\rangle-\langle U_{e}\rangle)=0.$ (91)

Hence
$\langle\underline{\nabla}\cdot\underline{\nabla}\cdot\underline{\underline{T}}\rangle=\underline{\nabla}\cdot\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle$
but
$\underline{\nabla}\cdot\langle\underline{\nabla}\cdot\underline{\underline{T}}\rangle\not=\underline{\nabla}\cdot\underline{\nabla}\cdot\langle\underline{\underline{T}}\rangle=0$
in view of the foregoing analysis, except asymptotically for
$kz\rightarrow+\infty$ or at discrete locations and frequencies $kz_{s}$ at
which the dispersion vanishes. In terms of the individual electric and
magnetic average stresses and energies, (91) with (88) degenerates into
separate equations, viz.,

$\displaystyle\underline{\nabla}\cdot\underline{\nabla}\cdot\langle\underline{\underline{T}}_{(e)(m)}\rangle\pm
k^{2}\langle U_{(e)(m)}\rangle=0$ (92)

where the upper and lower signs refer to
$(\langle\underline{\underline{T}}_{e}\rangle,\langle U_{e}\rangle)$ and
$(\langle\underline{\underline{T}}_{m}\rangle,\langle U_{m}\rangle)$,
respectively.

##  V Conclusion

The average LM (20) and AM (21) (and the corresponding linear and angular
power flux densities) of a random field near a PEC boundary exhibit a
dependence on the electric distance $kz$ that contains similarities as well as
differences compared to the dependence of the average electric or magnetic
energy densities (24)–(27). Similar to the energy, the strength of the average
momentum and power decays with increasing distance in a damped oscillatory
manner and vanishes far away from the boundary. By contrast, the average
momentum and power flux vanish on the boundary itself, their rate of decay
with $kz$ is different, and their asymptotic deep values are attained at
values of $kz$ in (23) that are interleaved with (i.e., spaced by
approximately $\pi/4$ from) those for the total energy density (28). As an
application to reverberation chambers, it was shown that by specifying a
maximum tolerable deviation (nonuniformity) of the average boundary power or
energy from its asymptotic free-space value, a performance based metric for
the size of the working volume of a chamber can be defined and calculated via
(43).

The Poynting theorem for conservation between energy imbalance and power flux
in the absence of ohmic losses was found to remain satisfied for their
(arithmetic) statistical averages in the case of random fields near a PEC
plane, extending its validity beyond deterministic fields to (50).
Conservation of the average LM and AM can be expressed either for the average
flux of the full EM stresses (65) and its moment (66), or in terms of the flux
of electric or magnetic stress (86) or their difference (88), and their
moments (87) and (89), respectively. The average EM energy imbalance and
stress are directly related as (91) or individually as (92).

As a final comment, LM and AM of random fields were found to exhibit a partial
and statistical behaviour, rather than being an all-or-nothing property in the
case of deterministic fields. This is similar to other wave characteristics of
inhomogeneous random fields, in particular for the degree of polarization
[18], [19]. This is not surprising in view of the fact that polarization
content is already comprised as SAM within AM (8).

## References

  * [1] J. H. Poynting, “The wave motion of a revolving shaft, and a suggestion as to the angular momentum in a beam of circularly polarised light,” Proc. R. Soc. Lond., Ser. A, vol. 82, pp. 560–567, 1909.
  * [2] R. Beth, “Mechanical detection and measurement of the angular momentum of light,” Phys. Rev., vol. 50, pp. 115–125, 1936.
  * [3] L. Allen, M. W. Beijersbergen, R. J. C. Spreeuw, and J. P. Woerdman, “Orbital angular-momentum of light and the transformation of Laguerre–Gaussian laser modes,” Phys. Rev. A, vol. 45, pp. 8185–8189, 1992.
  * [4] B. Thidé, F. Tamburini, H. Then, C. G. Someda, and R. A. Ravanelli, “The physics of angular momentum radio,” arXiv:1410.4268v3, 2015\.
  * [5] G. Ya. Slepyan, A. Boag, V. Mordachev, E. Sinkevich, S. Maksimenko, P. Kuzhir, G. Miano, M. E. Portnoi, and A. Maffucci, “Nanoscale electromagnetic compatibility: quantum coupling and matching in nanocircuits,” IEEE Trans. Electromagn. Compat., vol. 57, no. 6, pp. 1645–1654, Dec. 2015.
  * [6] M. Migliaccio, G. Gradoni, and L. R. Arnaut, “Electromagnetic reverberation: the legacy of Paolo Corona,” IEEE Trans. Electromagn. Compat., vol. 58, no. 3, pp. 643–652, Jun. 2016.
  * [7] J. M. Dunn, “Local, high-frequency analysis of the field in a mode-stirred chamber,” IEEE Trans. Electromagn. Compat., vol. 32, no. 1, pp. 53–58, Feb. 1990.
  * [8] L. R. Arnaut and P. D. West, “Electromagnetic reverberation near a perfectly conducting boundary,” IEEE Trans. Electromagn. Compat., vol. 48, no. 2, pp. 359–371, May 2006.
  * [9] L. R. Arnaut, “Spatial correlation functions of inhomogeneous random electromagnetic fields,” Phys. Rev. E, vol. 73, no. 3, 036604, Mar. 2006.
  * [10] J. D. Jackson, Classical Electrodynamics. 2nd ed. Wiley: New York, NY, 1975.
  * [11] J. A. Kong, Electromagnetic Wave Theory. 2nd ed. Wiley: New York, NY, 1990.
  * [12] M. V. Berry, “Optical currents,” J. Opt. A., vol. 11, 094001, 2009.
  * [13] D. A. Hill, Electromagnetic Theory of Reverberation Chambers. NIST Techn. Note 1506, U.S. Dept. of Commerce, Boulder, CO, Dec. 1998.
  * [14] M. Schröder, “Eigenfrequenzstatistik und Anregungsstatistik in Räumen: Modellversuche mit elektrischen Wellen,” Acustica, vol. 4, pp. 456–468, 1954.
  * [15] S. R. de Groot, The Maxwell Equations. North-Holland: Amsterdam, 1969.
  * [16] L. R. Arnaut, “Probability distribution of random electromagnetic fields in the presence of a semi-infinite isotropic medium,” Radio Sci., vol. 42, RS3001, 2007.
  * [17] G. B. Arfken, H. J. Weber, and F. E., Harris, Mathematical Methods for Physicists. 7th ed. Academic: Waltham, MA, 2013.
  * [18] L. R. Arnaut, “Compound exponential distributions for undermoded reverberation chambers,” IEEE Trans. Electromagn. Compat., vol. 44, no. 3, pp. 442–457, Aug. 2002.
  * [19] M. Migliaccio, J. J. Gil, A. Sorrentino, F. Nunziata, and G. Ferrara, “The polarization purity of the electromagnetic field in a reverberating chamber,” IEEE Trans. Electromagn. Compat., vol. 58, no. 3, pp. 694–700, Jun. 2016.
  * [20] I. V. Lindell, Methods for Electromagnetic Field Analysis. Clarendon: Oxford, U.K., 1992.

## Appendix: Maxwell Stress Dyadic for Time-Harmonic Random Fields

The Maxwell stress dyadic is usually formulated in the time domain as (e.g.,
[10, Sec. 6.8])

$\displaystyle\underline{\underline{{\rm
T}}}(\underline{r},t)\stackrel{{\scriptstyle\Delta}}{{=}}\epsilon_{0}\underline{{\rm
E}}\,\underline{{\rm E}}+\mu_{0}\underline{{\rm H}}\,\underline{{\rm
H}}-\frac{1}{2}\left(\epsilon_{0}\underline{{\rm
E}}^{2}+\mu_{0}\underline{{\rm H}}^{2}\right)\underline{\underline{I}}.$ (93)

With an assumed $\exp({\rm j}\omega t)$ dependence, the corresponding
expression for the complex $\underline{\underline{T}}(\underline{r},\omega)$
for a time harmonic random field ($\underline{E},\underline{H},\omega)$ can be
derived as follows. First, consider excitation by a single plane wave
($\underline{{\cal E}},\underline{{\cal H}},\underline{k})$ from the angular
spectrum of ($\underline{E},\underline{H})$. Being a sum of products of
complex field quantities, the complex Lorentz force can be written as

$\displaystyle\underline{{\cal
F}}(\underline{r},\underline{k})=\frac{1}{2}\left(\varrho\underline{{\cal
E}}^{*}+\underline{{\cal J}}\times\underline{{\cal B}}^{*}\right).$ (94)

Applying Gauss’s and Ampère’s laws to express $\varrho$ and $\underline{{\cal
J}}$ in terms of their source fields gives

$\displaystyle\underline{{\cal
F}}=\frac{1}{2}\left[(\underline{\nabla}\cdot\underline{{\cal
D}})\underline{{\cal E}}^{*}+(\underline{\nabla}\times\underline{{\cal
H}})\times\underline{{\cal B}}^{*}-{\rm j}\omega\underline{{\cal
D}}\times\underline{{\cal B}}^{*}\right].$ (95)

The last term in (95) equals $-{\rm j}\omega\underline{{\cal P}}$. Using
dyadic algebra and the vector identity $\underline{{\cal
H}}^{*}\times(\underline{\nabla}\times\underline{{\cal
H}})=\underline{\nabla}(\underline{{\cal H}}^{*}\cdot\underline{{\cal
H}})-(\underline{\nabla}\cdot\underline{{\cal H}}^{*})\underline{{\cal
H}}=\underline{\nabla}\,|\underline{{\cal
H}}|^{2}-\underline{\nabla}\cdot(\underline{{\cal H}}^{*}\underline{{\cal
H}})$, this yields in vacuum

$\displaystyle\underline{{\cal F}}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\left[\epsilon_{0}\underline{\nabla}\cdot(\underline{{\cal
E}}\,\underline{{\cal E}}^{*})+\mu_{0}\underline{\nabla}\cdot(\underline{{\cal
H}}^{*}\underline{{\cal H}})\right.$ (96)
$\displaystyle\left.-\mu_{0}\underline{\nabla}\cdot(|\underline{{\cal
H}}|^{2}\underline{\underline{I}})\right]-{\rm j}\omega\underline{{\cal
P}}=\underline{0}.$

Expressing ${\rm j}\omega\underline{{\cal P}}$ as
$\underline{\nabla}\cdot\underline{\underline{{\cal T}}}$ in view of (58)
gives

$\displaystyle\underline{\underline{{\cal
T}}}(\underline{r},\underline{k})=\frac{1}{2}\left(\epsilon_{0}\underline{{\cal
E}}\,\underline{{\cal E}}^{*}+\mu_{0}\underline{{\cal H}}^{*}\underline{{\cal
H}}-\mu_{0}|\underline{{\cal H}}|^{2}\underline{\underline{I}}\right).$ (97)

This expression for $\underline{\underline{{\cal T}}}$ is dual-asymmetric
because (94) is as such. Formally applying ${\cal H}={\cal E}/\eta_{0}$
enables the last term in (97) to be dual-symmetrized [20, sec. 4.2] as
$\mu_{0}|\underline{{\cal
H}}|^{2}\underline{\underline{I}}=\left(\epsilon_{0}|\underline{{\cal
E}}|^{2}+\mu_{0}|\underline{{\cal H}}|^{2}\right)\underline{\underline{I}}/2$.
Spherical integration of each dyad ${\cal E}_{\alpha}{\cal E}^{*}_{\beta}$ and
${\cal H}^{*}_{\alpha}{\cal H}_{\beta}$ according to (III-A) finally results
in the stress dyadic for the time-harmonic random field as

$\displaystyle\underline{\underline{T}}(\underline{r},\omega)$
$\displaystyle=$
$\displaystyle\frac{1}{2}\left[\epsilon_{0}\underline{E}\,\underline{E}^{*}+\mu_{0}\underline{H}^{*}\underline{H}\right.$
(98)
$\displaystyle\left.-\frac{1}{2}\left(\epsilon_{0}|\underline{E}|^{2}+\mu_{0}|\underline{H}|^{2}\right)\underline{\underline{I}}\right].$