# Unsupervised Deep Learning Image Verification Method

Enoch Solomon  Department of Computer Science
Virginia State University
Richmond, Virginia
<EMAIL_ADDRESS>Abraham Woubie  Silo AI
Helsinki, Finland
<EMAIL_ADDRESS>Eyael Solomon Emiru  Department of Information
Engineering and Computer Science
University of Trento, Italy
<EMAIL_ADDRESS>
###### Abstract

Although deep learning are commonly employed for image recognition, usually
huge amount of labeled training data is required, which may not always be
readily available. This leads to a noticeable performance disparity when
compared to state-of-the-art unsupervised face verification techniques. In
this work, we propose a method to narrow this gap by leveraging an autoencoder
to convert the face image vector into a novel representation. Notably, the
autoencoder is trained to reconstruct neighboring face image vectors rather
than the original input image vectors. These neighbor face image vectors are
chosen through an unsupervised process based on the highest cosine scores with
the training face image vectors. The proposed method achieves a relative
improvement of 56% in terms of EER over the baseline system on Labeled Faces
in the Wild (LFW) dataset. This has successfully narrowed down the performance
gap between cosine and PLDA scoring systems.

###### Index Terms:

autoencoder, biometrics, deep learning, face recognition, image recognition

††publicationid: pubid: 979-8-3503-9341-5/23/$31.00 ©2023 IEEE

##  I Introduction

Face recognition is a system used for identifying or verifying a person from
an image or video [1, 2, 4, 5, 6, 7, 8, 3]. Deep learning techniques, as a
front-end approach, have demonstrated proficiency in acquiring intricate
features [9], including bottleneck features (BNF) [10]. These features are
subsequently integrated into traditional frameworks or employed in the
extraction of face image vectors for the face recognition process. Deep
learning methods prove effective in compressing face images into embeddings,
as exemplified in [11] in addition to different speech and security
applications [12, 15, 13, 14, 16, 17, 18, 19, 20, 21, 22].

Several works have been proposed to enhance face image verification
performance through unsupervised learning techniques [23, 24, 25]. For
instance, [26] introduces a vector representation of face images using RBM
adaptation. On the backend, [27] proposes various imposter selection
algorithms aimed at narrowing the performance gap between supervised and
unsupervised face recognition methods. They apply DBN adaptation as a backend
for face images. These algorithms necessitate training a separate model for
each target person, incurring computational costs. However, the results still
fall short of the supervised approach, which benefits from actual face image
labels.

In this work, we focus on mitigating the necessity for labeled face image data
in face verification. Our objective is to narrow the performance disparity
between supervised and unsupervised techniques, all without relying on face
image labels. Unlike traditional DNN-based classifiers, the training of
autoencoders is an unsupervised process that operates independent of labeled
data. Our contribution lies in the novel training approach for the
autoencoder, specifically designed to account for session variability among
face images in the absence of labeled data. Rather than reconstructing the
same training face image, we train the autoencoder to reconstruct neighboring
face images. Following this training, we extract face image vector embeddings
for testing purposes. During experimental trials, we apply cosine scoring to
evaluate the face image embedding vectors, which have exhibited a remarkable
discriminatory capability. The experimental findings demonstrate a substantial
56% enhancement in performance over the baseline system that employs the
cosine scoring technique, underscoring the efficacy of our proposed
autoencoder training methodology.

The rest of the paper is as follows: Section 2 provides an in-depth
description of the proposed autoencoder training method and the process of
neighbor face image selection. Section 3 describes the experimental
configuration and outlines the dataset employed. Section 4 describes the
experimental results. Finally, Section 5 gives the conclusion of the proposed
work.

Figure 1: Proposed approach for training the autoencoder.

##  II Proposed Method

In this work, we introduce an unsupervised deep learning method for face
verification. The training procedure of the proposed method is outlined as
follows: From a pool of unlabeled large-scale face image data, we identify the
k nearest neighbor images that closely resemble the specified face image,
which are then employed for training. Diverging from the conventional deep
neural network training, our method computes the loss function for each input
face image k times.

PLDA scoring technique for face image vectors necessitates labels [28], which
can be challenging to obtain in practical scenarios. Conversely, the cosine
scoring technique operates without the need for labels, albeit at the expense
of performance degradation. To bridge the achievement gap between these two
scoring methods in the absence of labeled training data, we introduce a novel
framework for autoencoder training. This approach is entirely unsupervised,
distinguishing it from conventional DNN classifiers and PLDA. Additionally, we
advocate training the autoencoder to reconstruct similar face images, rather
than exclusively focusing on the same training face image. Through this
method, the autoencoder, trained in an unsupervised manner, effectively
addresses session variability among face images without relying on labels.

###  II-A Autoencoder training

As shown in Fig 1, an autoencoder comprises an encoder and a decoder. The
encoder functions to compress the input face image vector $x$ into a lower-
dimensional space, while the decoder reconstructs it back to its original
form. Conventional training involves minimizing the Mean Square Error (MSE)
between the input $X$ and the reconstructed $X^{\wedge}$.

In this work, we proposed the training of an autoencoder to minimize the loss
function $MSE\left(X^{\wedge},N\right)$, illustrated in Figure 1, where $N$
represents a face image vector similar to $X$ and $X^{\wedge}=$ decoder
$(\operatorname{encoder}(X))$. We put forth a method for the automatic
selection of similar face image vectors. For each training input face image
vector $X$, multiple similar face image vectors can be taken into
consideration. Thus, each images in the dataset could have different number of
similar images. In Section 4, we will conduct a comparative analysis between
the results obtained from our proposed training approach and the conventional
training method, which involves reconstructing the same training face image
vectors.

###  II-B Selection of similar face image vectors

All the training face image vectors are assessed against one another using the
cosine scoring technique. For each face image vector, a group of similar face
image vectors is chosen. A simple method for selecting these similar face
image vectors involves setting a threshold on the cosine scores between the
training face image vectors. Face image vectors with scores surpassing this
threshold are designated as similar face image vectors. It’s worth noting that
this approach may result in a varying number of similar face image vectors for
each training face image vector.

Another approach is to adopt a constant value, denoted as k, and choose k
similar face image vectors for each training face image vector. This method
ensures a consistent number of similar face images for every training face
image vector, resulting in a balanced training approach.

Figure 2 provides a visual representation of the process for selecting similar
face image vectors when a constant $k$ number of neighbors is considered. Let
$X_{i}$ be a training face image vector, where $i=(1,\ldots,n)$ and $n$
represents the total number of training face image vectors. Initially, we
compute cosine scores for all training face image vectors relative to one
another. Subsequently, we pick the top $k$ face image vectors with the highest
scores to serve as similar face image vectors for each $X_{i}$. These selected
similar face image vectors are denoted as $N_{ij}$, where $N_{ij}$ is the
$j^{\text{th }}$ similar of $i^{\text{th }}$ training face image vectors.

Thus, we have a total of $n\times(k-1)$ training data, which automatically
generated on the fly which address the issue of having large training data. We
run extensive experiment to determine the the values of threshold and $k$.
Thus, the autoencoder is able to learn information about the intrinsic
variability of face image vectors, without necessarily using actual face image
labels.

Consequently, we accumulate a total of $n\times(k-1)$ samples for the
autoencoder training process. The specific values for the threshold and $k$
are established through experimental investigation, and a detailed discussion
on these parameters will be provided in Section 4. Through this approach, the
autoencoder effectively assimilates information regarding the session
variability inherent in face image vectors, all without the dependency on
actual face image labels.

###  II-C Autoencoder vector extraction

After the autoencoder is trained with the chosen neighbor face image vectors,
we proceed to convert the testing face image vectors into a fresh face image
vector representation, utilizing the autoencoder, as illustrated in Figure 1.
At this stage, we extract the target face image vectors from the output of the
autoencoder. In the experimental trials, these extracted vectors have
demonstrated an enhanced discriminatory capability for face image vectors, all
achieved without reliance on face image labels. Leveraging these extracted
vectors, we carry out the experimental trials using the cosine scoring
technique.

Figure 2: Proposed method of selecting the neighbor face image.

##  III Experimental Setup and Dataset

###  III-A Experimental Setup

For the training of our model, we utilized the Keras deep learning library
[29]. The architecture of the autoencoder employed in this work comprises
three hidden layers. This design maintains symmetry between the encoder and
decoder sections. Specifically, both the first and third hidden layers contain
an identical number of neurons, as illustrated in Fig. 2. Each of the encoder
input and decoder output layers consists of 112 by 112 neurons. In the second
layer of both the encoder and decoder, there are 800 neurons. The encoder’s
output layer is composed of 300 neurons, determining the size of the resulting
embedding vector. Within the DNN, the dimension of the face image embedding
layer was fixed at 400, while the classification layer encompassed 1000
neurons.

The training of the autoencoder was conducted over 400 epochs or until the
error stops decreasing and showed no further decrease. In all layers of the
autoencoder, we applied the Rectified Linear Unit (ReLU) activation function,
with the exception of the final layer, which utilized a linear activation
function. The learning rate was set to 0.03, incorporating a decay rate of
0.0002. We implemented a logarithmically decaying learning rate, spanning from
$10^{-2}$ to $10^{-8}$. The batch size was configured at 100.

###  III-B Datasets

The CelebA dataset, introduced by Liu et al. [9], encompasses a vast
collection of over 200,000 images featuring 10,177 celebrities. This dataset
incorporates diverse elements such as pose variations and background clutter,
providing a comprehensive representation of real-world scenarios. To
facilitate model training and evaluation, the dataset is partitioned into
distinct sets for training, validation, and testing purposes. In this study,
we leveraged the training and test segments of the CelebA dataset as our
unlabeled dataset for training the autoencoder. In contrast, supervised
training was conducted using the validation portion of the CelebA dataset.

The Labeled Faces in the Wild dataset (LFW) [30] comprises a collection of
13,233 images featuring 5,749 distinct individuals. For testing purposes, this
database is randomly and uniformly divided into ten subsets. Within each
subset, 300 pairs of matched images (depicting the same person) and 300 pairs
of mismatched images (depicting different persons) are randomly selected. In
essence, this testing protocol employs a total of 3,000 matched pairs and
3,000 mismatched pairs [30].

##  IV Experimental Results

###  IV-A Comparison with different threshold values

TABLE I: The Equal Error Rate (EER) was computed for both the proposed face image vectors and the baseline, assessed across various threshold values when employing cosine scoring. Approach |  threshold |  EER(%)
---|---|---
[1] Baseline | No | 16.84
[2] Proposed | 0.7 | 15.30
[3] Proposed | 0.6 | 15.04
[4] Proposed | 0.5 | 13.82
[5] Proposed | 0.4 | 11.23
[6] Proposed | 0.3 | 9.61
[7] Proposed | 0.2 | 8.13
[8] Proposed | 0.1 | 8.24
Fusion of [1] and [7] | combined | 7.12

In Tables I and II, approach refers to the method where face image vectors are
extracted using a conventionally trained autoencoder. Table I presents a
performance comparison between our proposed method and the baseline, where a
threshold is set to select a specific number of neighbor face image vectors.
All vectors are then evaluated using the cosine scoring technique. We
experimented with different threshold values to fine-tune the system and
achieve the best results. The table clearly illustrates that our proposed
method outperforms the baseline system. As we decrease the threshold value,
the system’s performance improves. The optimal Equal Error Rate (EER) of 8.24%
was achieved with a threshold set to 0.2. This represents a 51% relative
improvement over the baseline method. Consequently, this bridges the
performance gap between cosine and Probabilistic Linear Discriminant Analysis
(PLDA) scoring techniques.

By performing a score-level fusion of the baseline and our proposed method,
using a threshold of 0.2, we achieved a further enhancement in performance, as
evidenced by a reduced Equal Error Rate (EER) of 7.12%. The optimal weights
for the fusion were determined empirically, resulting in values of 0.55 and
0.45 for the baseline and our proposed method, respectively. This fusion
strategy effectively leverages the strengths of both approaches to yield
superior results.

###  IV-B Comparison with different k values

TABLE II: The Equal Error Rate (EER) was computed for both the proposed face image vector and the baseline method across different values of k all evaluated using the cosine scoring technique. Approach |  k |  EER(%)
---|---|---
[1] Baseline | No | 18.54
[2] Proposed | 1 | 15.21
[3] Proposed | 2 | 13.13
[4] Proposed | 3 | 12.71
[5] Proposed | 4 | 10.12
[6] Proposed | 10 | 8.40
[7] Proposed | 11 | 8.03
[8] Proposed | 12 | 8.94
Fusion of [1] and [11] | combined | 7.01

In Table II, we present a performance comparison between our proposed method
and the baseline, utilizing a constant value of k for the selection of similar
face image vectors. We conducted experiments with various values of k to
optimize the system for achieving the best results. The table illustrates that
employing a fixed value of k has demonstrated improvement compared to both the
threshold-based approach and the baseline. Beginning with k set to 1, an Equal
Error Rate (EER) of 15.21% was attained. This signifies a 17% relative
improvement over the baseline system by considering only the nearest neighbor
for each image vector. As we increment the value of k, the system’s
performance continues to improve. The most favorable EER of 8.03% was achieved
when k was set to 11, leading to a remarkable 56% relative improvement over
the baseline method. This achievement effectively bridges the performance gap
between cosine and PLDA scoring techniques. Moreover, this approach enables a
balanced training regimen, ultimately enhancing the system’s performance. It
is worth noting that a further increase in the value of k may incorporate very
distant neighbors, potentially degrading the system’s performance.

By combining the scores from the baseline and our proposed method, both
utilizing k set to 11, we achieved a further improvement in performance,
specifically in terms of the Equal Error Rate (EER). This fusion approach
yielded an EER of 7.01%, which is in close proximity to the results obtained
using Probabilistic Linear Discriminant Analysis (PLDA) with actual face image
labels. The optimal weights for the fusion were fine-tuned through
experimental iteration and were ultimately set to 0.49 and 0.51 for the
baseline and our proposed method, respectively. This fusion strategy
effectively leverages the strengths of both approaches to yield highly
competitive results.

###  IV-C Comparison with score level fusion

TABLE III: The Equal Error Rate (EER) was computed for the score-level fusion between Probabilistic Linear Discriminant Analysis (PLDA) and cosine scoring, using both the proposed face image vectors with k set to 15, and the baseline method. Approach |  Scoring |  EER(%)
---|---|---
[1] Baseline | PLDA | 6.71
[2] Proposed (k = 11) | Cosine | 8.03
Fusion of [1] and [2] | combined | 6.23

As depicted in Table III, conducting a score-level fusion between
Probabilistic Linear Discriminant Analysis (PLDA) and cosine scoring,
utilizing the proposed method with k set to 11 (yielding the best results),
leads to an Equal Error Rate (EER) of 6.23%. This showcases an impressive
relative improvement of nearly 7% over the PLDA approach alone. The optimal
fusion weights were determined to be 0.04 for the baseline and 0.96 for the
proposed method, respectively. This fusion strategy effectively capitalizes on
the strengths of both approaches, resulting in a highly competitive
performance.

###  IV-D Comparison with the state of the art on LFW dataset

TABLE IV: Comparison of the proposed method result with some of state-of-the-art methods on LFW testing dataset. Model |  Training size |  Labeled/ Unlabeled |  Testing size |  Accuracy(%)
---|---|---|---|---
UniformFace | 6.1M | Labeled | 6K | 99.80 [31]
ArcFace | 5.8M | Labeled | 6K | 99.82 [32]
GroupFace | 5.8M | labeled | 6K | 99.85 [33]
CosFace | 5M | Labeled | 6K | 99.73 [34]
Marginal Loss | 4M | Labeled | 6K | 99.48 [35]
CurricularFace | 3.8M | Labeled | 6K | 99.80 [36]
RegularFace | 3.1M | Labeled | 6K | 99.61 [37]
AFRN | 3.1M | Labeled | 6K | 99.85 [38]
Stream Loss | 1.5M | Labeled | 6K | 98.97 [39]
MDCNN | 1M | Labeled | 6K | 99.38 [40]
ULNet | 1M | Labeled | 6K | 99.70 [41]
Ben Face | 0.5M | Labeled | 6K | 99.20 [42]
F2C | 5.8M | Labeled | 6K | 99.83 [43]
PCCycleGAN | 0.5M | Unlabeled | 6K | 99.52 [44]
CAPG GAN | 1M | Unlabeled | 6K | 99.37 [45]
UFace | 200K | Unlabeled | 6K | 99.40 [46]
Proposed | 200K | Unlabeled | 6K | 99.53

Despite the slightly higher accuracy of methods like ArcFace, GroupFace,
Marginal Loss, and CosFace, it’s noteworthy that our proposed method achieves
competitive results with a significantly smaller training dataset (around 200K
images). In contrast, many state-of-the-art methods rely on much larger
training sets, often in the millions of images. This demonstrates the
efficiency and effectiveness of our approach, particularly in scenarios where
acquiring a massive dataset may not be feasible or practical.

The experimental result of the proposed method become particularly notable
when compared to state-of-the-art approaches. For instance, ArcFace, which
attains 99.82% accuracy, utilized a substantially larger dataset of 5.8
million labeled images. In contrast, our method achieved comparable
performance with significantly fewer training images (around 200K). This
underscores the efficiency and effectiveness of our approach, especially in
scenarios where obtaining a massive dataset is challenging.

##  V Conclusion

The proposed work introduces a novel unsupervised deep learning method for
face verification. The training process involves selecting k images from
unlabeled data that are most similar to a specified face image, and using them
for training. Unlike standard deep neural network training, our proposed
method calculates the loss function for each input image k times. The primary
goal is to bridge the performance gap between cosine and PLDA scoring
techniques when labeled training data is unavailable. This is achieved by
training an autoencoder to reconstruct a set of neighbor face image vectors,
instead of reconstructing the same training face image vectors. These neighbor
face image vectors are selected based on their cosine scores in relation to
the training face image vector. After training, face image vectors are
extracted for test face images as the output of the autoencoder. In
experimental trials, face image vectors are scored using cosine scoring. The
evaluation conducted on the LFW dataset demonstrates that our proposed method
achieves a relative improvement of 56% in terms of Equal Error Rate (EER) over
the baseline system. This successfully narrows the performance gap between
cosine and PLDA scoring systems.

## References

  * [1] Jain, A.; Li, S. Handbook of Face Recognition; Springer: Berlin/Heidelberg, Germany,
  * [2] Solomon, E. Face Anti-Spoofing and Deep Learning Based Unsupervised Image Recognition Systems. (2023)
  * [3] Solomon, E. Face Anti-Spoofing and Deep Learning Based Unsupervised Image Recognition Systems. VCU Theses And Dissertations. (2023), https://scholarscompass.vcu.edu/etd/7482
  * [4] Solomon, E. & Cios, K. FASS: Face Anti-Spoofing System Using Image Quality Features and Deep Learning. Electronics. 12, 2199 (2023)
  * [5] Solomon, E. & Cios, K. HDLHC: Hybrid Face Anti-Spoofing Method Concatenating Deep Learning and Hand-Crafted Features. 2023 IEEE 6th International Conference On Electronic Information And Communication Technology (ICEICT). pp. 470-474 (2023)
  * [6] Solomon, E., Woubie, A. & Emiru, E. Deep Learning Based Face Recognition Method using Siamese Network. ArXiv Preprint ArXiv:2312.14001. (2023)
  * [7] Solomon, E., Woubie, A. & Emiru, E. Autoencoder Based Face Verification System. ArXiv Preprint ArXiv:2312.14301. (2023)
  * [8] Woubie, A., Solomon, E. & Emiru, E. Image Clustering using Restricted Boltzman Machine. ArXiv Preprint ArXiv:2312.13845. (2023)
  * [9] Liu, Z., Luo, P., Wang, X. & Tang, X. Deep learning face attributes in the wild. Proceedings Of The IEEE International Conference On Computer Vision. pp. 3730-3738 (2015)
  * [10] Deng, L., Yu, D. & Others Deep learning: methods and applications. Foundations And Trends® In Signal Processing. 7, 197-387 (2014)
  * [11] Snyder, D., Garcia-Romero, D., Sell, G., Povey, D. & Khudanpur, S. X-vectors: Robust dnn embeddings for speaker recognition. 2018 IEEE International Conference On Acoustics, Speech And Signal Processing (ICASSP). pp. 5329-5333 (2018)
  * [12] Abraham Woubie and Tom Bäckström: Federated learning for privacy-preserving speaker recognition, In IEEE Access, 2021; pp. 149477–149485.
  * [13] Abraham Woubie and Tom Bäckström: Federated Learning for Privacy Preserving On-Device Speaker Recognition, In ISCA Symposium on Security and Privacy in Speech Communication, 2021; pp. 1–5.
  * [14] Abraham Woubie, Pablo Zarazaga and Tom Bäckström: The Use of Audio Fingerprints for Authentication of Speakers on Speech Operated Interfaces, In ISCA Symposium on Security and Privacy in Speech Communication, 2021; pp. 6–9.
  * [15] Abraham Woubie and Tom Bäckström: Voice-quality Features for Deep Neural Network Based Speaker Verification Systems, In 29th European Signal Processing Conference (EUSIPCO), 2021; pp. 176–180.
  * [16] Chekole, E., Chattopadhyay, S., Ochoa, M., Guo, H. & Cheramangalath, U. CIMA: Compiler-Enforced Resilience Against Memory Safety Attacks in Cyber-Physical Systems. Computers & Security. 94 pp. 101832 (2020), https://www.sciencedirect.com/science/article/pii/S0167404820301061
  * [17] Chekole, E., Ochoa, M. & Chattopadhyay, S. SCOPE: Secure Compiling of PLCs in Cyber-Physical Systems. International Journal Of Critical Infrastructure Protection. 33 pp. 100431 (2021), https://www.sciencedirect.com/science/article/pii/S1874548221000238
  * [18] Chekole, E. & Huaqun, G. ICS-SEA: Formally Modeling the Conflicting Design Constraints in ICS. Proceedings Of The Fifth Annual Industrial Control System Security (ICSS) Workshop. pp. 60-69 (2019), https://doi.org/10.1145/3372318.3372325
  * [19] Chekole, E., Castellanos, J., Ochoa, M. & Yau, D. Enforcing Memory Safety in Cyber-Physical Systems. Computer Security. pp. 127-144 (2018)
  * [20] Chekole, E., Chattopadhyay, S., Ochoa, M. & Huaqun, G. Enforcing Full-Stack Memory-Safety in Cyber-Physical Systems. Engineering Secure Software And Systems. pp. 9-26 (2018)
  * [21] Chekole, E., Thulasiraman, R. & Zhou, J. EARIC: Exploiting ADC Registers in IoT and Control Systems. Applied Cryptography And Network Security Workshops. pp. 245-265 (2023)
  * [22] Wong, A., Chekole, E., Ochoa, M. & Zhou, J. On the Security of Containers: Threat Modeling, Attack Analysis, and Mitigation Strategies. Computers & Security. 128 pp. 103140 (2023), https://www.sciencedirect.com/science/article/pii/S0167404823000500
  * [23] Liu, P., Han, S., Meng, Z. & Tong, Y. Facial expression recognition via a boosted deep belief network. Proceedings Of The IEEE Conference On Computer Vision And Pattern Recognition. pp. 1805-1812 (2014)
  * [24] Kurup, A., Ajith, M. & Ramón, M. Semi-supervised facial expression recognition using reduced spatial features and deep belief networks. Neurocomputing. 367 pp. 188-197 (2019)
  * [25] Teh, Y. & Hinton, G. Rate-coded restricted Boltzmann machines for face recognition. Advances In Neural Information Processing Systems. 13 (2000)
  * [26] Alphonse, A., Shankar, K., Jeyasheela Rakkini, M., Ananthakrishnan, S., Athisayamani, S., Robert Singh, A. & Gobi, R. A multi-scale and rotation-invariant phase pattern (MRIPP) and a stack of restricted Boltzmann machine (RBM) with preprocessing for facial expression classification. Journal Of Ambient Intelligence And Humanized Computing. 12 pp. 3447-3463 (2021)
  * [27] Nguyen, H. & Bai, L. Cosine similarity metric learning for face verification. Asian Conference On Computer Vision. pp. 709-720 (2010)
  * [28] Muslihah, I. & Muqorobin, M. Texture characteristic of local binary pattern on face recognition with probabilistic linear discriminant analysis. International Journal Of Computer And Information System (IJCIS). 1, 22-26 (2020)
  * [29] Chollet, F. & Others Keras. https:// keras. io. (2015)
  * [30] Huang, G., Mattar, M., Berg, T. & Learned-Miller, E. Labeled faces in the wild: A database forstudying face recognition in unconstrained environments. (2008)
  * [31] Duan, Y., Lu, J. & Zhou, J. Uniformface: Learning deep equidistributed representation for face recognition. Proceedings Of The IEEE/CVF Conference On Computer Vision And Pattern Recognition. pp. 3415-3424 (2019)
  * [32] Deng, J., Guo, J., Xue, N. & Zafeiriou, S. Arcface: Additive angular margin loss for deep face recognition. Proceedings Of The IEEE/CVF Conference On Computer Vision And Pattern Recognition. pp. 4690-4699 (2019)
  * [33] Kim, Y., Park, W., Roh, M. & Shin, J. Groupface: Learning latent groups and constructing group-based representations for face recognition. Proceedings Of The IEEE/CVF Conference On Computer Vision And Pattern Recognition. pp. 5621-5630 (2020)
  * [34] Wang, H., Wang, Y., Zhou, Z., Ji, X., Gong, D., Zhou, J., Li, Z. & Liu, W. Cosface: Large margin cosine loss for deep face recognition. Proceedings Of The IEEE Conference On Computer Vision And Pattern Recognition. pp. 5265-5274 (2018)
  * [35] Deng, J., Zhou, Y. & Zafeiriou, S. Marginal loss for deep face recognition. Proceedings Of The IEEE Conference On Computer Vision And Pattern Recognition Workshops. pp. 60-68 (2017)
  * [36] Huang, Y., Wang, Y., Tai, Y., Liu, X., Shen, P., Li, S., Li, J. & Huang, F. Curricularface: adaptive curriculum learning loss for deep face recognition. Proceedings Of The IEEE/CVF Conference On Computer Vision And Pattern Recognition. pp. 5901-5910 (2020)
  * [37] Zhao, K., Xu, J. & Cheng, M. Regularface: Deep face recognition via exclusive regularization. Proceedings Of The IEEE/CVF Conference On Computer Vision And Pattern Recognition. pp. 1136-1144 (2019)
  * [38] Kang, B., Kim, Y., Jun, B. & Kim, D. Attentional feature-pair relation networks for accurate face recognition. Proceedings Of The IEEE/CVF International Conference On Computer Vision. pp. 5472-5481 (2019)
  * [39] Rashedi, E., Barati, E., Nokleby, M. & Chen, X. “Stream loss”: ConvNet learning for face verification using unlabeled videos in the wild. Neurocomputing. 329 pp. 311-319 (2019)
  * [40] Huang, X.; Zeng, X.; Wu, Q.; Lu, Y.; Huang, X.; Zheng, H. Face Verification Based on Deep Learning for Person Tracking in Hazardous Goods Factories. Processes 2022, 10, 380.
  * [41] Boragule, A.; Akram, H.; Kim, J.; Jeon, M. Learning to Resolve Uncertainties for Large-Scale Face Recognition. Pattern Recognit. Lett. 2022, _160_ , 58–65.
  * [42] Ben Fredj, H.; Bouguezzi, S.; Souani, C. Face recognition in unconstrained environment with CNN. Vis. Comput. 2021, 37, 217–226.
  * [43] Wang, K.; Wang, S.; Zhang, P.; Zhou, Z.; Zhu, Z.; Wang, X.; Peng, X.; Sun, B.; Li, H.; You, Y. An efficient training approach for very large scale face recognition. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, New Orleans, Louisiana, USA, 21-24 June 2022; pp. 4083–4092.
  * [44] Liu, Y.; Chen, J. Unsupervised face frontalization for pose-invariant face recognition. Image Vis. Comput. 2021, _106_ , 104093.
  * [45] Hu, Y.; Wu, X.; Yu, B.; He, R.; Sun, Z. Pose-guided photorealistic face rotation. In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, Salt Lake City, UT, USA, 18-23 June 2018; pp. 8398–8406.
  * [46] Solomon, E., Woubie, A. & Cios, K. UFace: An Unsupervised Deep Learning Face Verification System. Electronics. 11, 3909 (2022)


# On the effect of localized and delocalized quasiparticles on a nonlinear
conductivity of superconductors

A. Gurevich  Department of Physics and Center for Accelerator Science, Old
Dominion University, Norfolk VA 23510

###### Abstract

Recently Deyo et al. [Phys. Rev. B106, 104502 (2022)] suggested that the
surface resistance $R_{s}(H_{0})$ of a superconductor under strong
electromagnetic field can be affected by field-induced quasiparticle bound
states at the surface. This Comment shows that the existence of such bound
states and their contribution to $R_{s}$ are not substantiated and the
phenomenological models of Ref. corn  give incorrect field and temperature
dependencies of $R_{s}(T,H_{0})$.

###### pacs:

74.25.-q, 74.25.Ha, 74.25.Op, 74.78.Na

Deyo et al. analyzed the surface states in a clean superconductor, using the
Bogoliubov-de-Gennes (BdG) equations for a film of thickness $L=5\lambda$ in a
magnetic field $H(t)=H_{0}\sin\omega t$ applied along the $x$-axis:

$\displaystyle\\!\\!\\!\\!\frac{\hbar^{2}}{2m}\left[k_{x}^{2}+\left(k_{y}-\frac{eA}{\hbar
c}\right)^{2}\\!\\!-\frac{d^{2}}{dz^{2}}\right]\\!u+\Delta v=(E+E_{F})u,$ (1)
$\displaystyle\\!\\!\\!\\!\frac{\hbar^{2}}{2m}\left[k_{x}^{2}+\left(k_{y}+\frac{eA}{\hbar
c}\right)^{2}\\!\\!-\frac{d^{2}}{dz^{2}}\right]\\!v-\Delta u=(E_{F}-E)v,$ (2)
$\displaystyle u(0)=v(0)=0,\qquad u(L)=v(L)=0.$ (3)

Here $[u({\bf r}),v({\bf r})]=[u(z),v(z)]e^{ik_{x}x+ik_{y}y}$,
$A(z,t)=H_{0}\lambda e^{-z/\lambda}\sin\omega t$, $\lambda$ is the London
penetration depth, and the energy $E$ is counted from the Fermi surface.

Figure 1:  The BdG quasiparticle bound states for $k_{x}=k_{F}$, $k_{y}=0$, at
a field of 40 MV/m. For each state the wave functions $u$ and $v$ are plotted
with a vertical shift according to their energy. The situation is somewhat
analogous to a particle in a potential $\Delta(z)$ with a hard wall on the
left. Reproduced from Ref. corn .

Fig. 1 shows the bound states in a deep well in $\Delta(z)$ for quasiparticles
moving parallel to the surface with velocities $v_{x}=v_{F}$ obtained in Ref.
corn . Here the well in $\Delta(z)$ is caused by the Meissner currents but not
the boundary conditions (3) as it only occurs at the film side exposed to the
magnetic field. The strong suppression of $\Delta(z)$ at the surface shown in
Fig. 1 appears inconsistent with the magnitude of current pairbreaking tinkh .
This can be seen from the Ginzburg-Landau (GL) equation,

$\xi^{2}\psi^{\prime\prime}+\bigl{[}1-\left(2\pi\xi\lambda
H_{0}/\phi_{0}\right)^{2}e^{-2z/\lambda}\bigr{]}\psi-\psi^{3},$ (4)

where $\psi(z)=\Delta(z)/\Delta_{0}$, $\Delta_{0}$ is the order parameter at
zero field, $\xi$ is the coherence length, and $\phi_{0}$ is the flux quantum.
The term $\propto H_{0}^{2}$ in Eq. (4) causes a weak suppression in
$\psi(z)=1+\delta\psi(z)$, where the field-induced correction $\delta\psi(z)$
satisfies the linearized GL equation:

$\xi^{2}\delta\psi^{\prime\prime}-2\delta\psi=(H_{0}/\sqrt{2}H_{c})^{2}e^{-2z/\lambda},$
(5)

where $H_{c}=\phi_{0}/2^{3/2}\xi\lambda$ is the GL thermodynamic critical
field. The solution with $\Delta^{\prime}(0)=0$ like in Fig 1 is:

$\frac{\Delta(z)}{\Delta_{0}}=1+\frac{(H_{0}/2H_{c})^{2}}{2(\xi/\lambda)^{2}-1}\biggl{[}e^{-2x/\lambda}-\frac{\xi\sqrt{2}}{\lambda}e^{-\sqrt{2}x/\xi}\biggr{]}.$
(6)

If $\kappa=\lambda/\xi\gg 1$ the well in $\Delta(z)$ at $H_{0}=H_{c1}$ is very
shallow, $\delta\psi\sim-\kappa^{-2}$. Likewise,
$\delta\psi\sim-2^{-5/2}\kappa$ at $\kappa\ll 1$ and $H_{0}=H_{c}$. For a
clean Nb with $\xi=\lambda$, Eq. (6) also gives a very shallow well
$\Delta(0)/\Delta_{0}=0.994-0.974$ in the field region
$H_{0}=(0.25-0.5)H_{c}$, where the decrease of $R_{s}(H_{0})$ with $H_{0}$ in
Nb cavities has been observed raise . For $H_{0}\approx H_{c}^{Nb}=200$ mT,
Eq. (6) gives $\Delta(0)/\Delta_{0}\approx 0.9$, much different from
$\Delta(0)/\Delta_{0}\approx 0.25$ in Fig. 1. At $T\ll T_{c}$ the well in
$\Delta(z)$ in the clean limit is even shallower than the above GL estimates
because a uniform current has no effect on $\Delta$ at $T=0$ and the
superfluid velocities $v_{s}<v_{c}=\Delta/p_{F}$ maki . However, to obtain
this result, the calculation of $\Delta(z)$ at $T\ll T_{c}$ requires not only
summations over $k_{x}$ and $k_{y}$ but also taking into account lots of
($\sim Lk_{F}\sim 10^{3}$ for Nb) rapidly oscillating eigenfunctions
$u_{n}(z)$ and $v_{n}(z)$ with periods $l({\bf
v}_{F})\sim\hbar/mv_{z}\ll\lambda$ and $|E|<2\hbar\omega_{D}$, where $v_{z}$
is the $z$ component of the Fermi velocity ${\bf v}_{F}$. There is no
information in Ref. corn  about how the results shown in Fig. 1 were produced
and what eigenstates were taken into account in the self-consistency equation
for $\Delta(z)$. If only $u(z)$ and $v(z)$ varying slowly over
$\lambda_{F}=2\pi/k_{F}$ were taken into account, while lots of itinerant
states with $v_{z}\sim v_{F}$ were disregarded, as was briefly mentioned in
Ref. corn , the bound states shown in Fig. 1 are not real.

Equations (1)-(2) are ill-suited for the numerical calculations of itinerant
states $u({\bf r})=e^{i{\bf k_{F}r}}\alpha({\bf r})$ and $v({\bf r})=e^{i{\bf
k_{F}r}}\beta({\bf r})$, where $\alpha({\bf r})$ and $\beta({\bf r})$ vary
slowly over $\lambda_{F}$. In this quasi-classical Andreev approximation the
BdG equations become:

$\displaystyle-i\hbar
v_{z}\alpha^{\prime}+v_{y}p_{s}\alpha+\Delta\beta=E_{-}\alpha$ (7)
$\displaystyle i\hbar
v_{z}\beta^{\prime}+v_{y}p_{s}\beta+\Delta^{*}\alpha=E_{+}\beta$ (8)

Here $E_{\pm}=E\pm\hbar(v_{x}s_{x}+v_{y}s_{y})$, ${\bf s}={\bf k}-{\bf
k}_{F}$, and ${\bf p}_{s}=m\hbar{\bf v}_{s}$ defines the Meissner current
density ${\bf J}=en_{s}{\bf v}_{s}$. If $p_{s}$ is independent of $z$, the
solutions $(\alpha,\beta)\propto e^{is_{z}z}$ give the well-known anisotropic
energy spectrum bardeen ; maki ; tinkh :

$E({\bf k})=p_{s}v_{z}\pm\sqrt{\xi_{{\bf k}}^{2}+\Delta^{2}},$ (9)

where $\xi_{\bf k}=\hbar(({\bf k-k}_{F})\cdot{\bf v}_{F}$. As follows from Eq.
(9), the energy gap $E_{g}({\bf v}_{F})=\Delta_{0}-{\bf v}_{F}{\bf p}_{s}$
depends on the direction of ${\bf v}_{F}$ relative to ${\bf J}$, the minimum
gap $E_{g}(t)=\Delta_{0}-|p_{s}(t)|v_{F}$ is smaller than $\Delta_{0}$ bardeen
. Inserting Eq. (9) in the self-consistency equation for $\Delta$ and
integrating over the orientations of ${\bf v}_{F}$ at $T=0$ gives $\Delta$
independent of $p_{s}$ if $p_{s}<\Delta/v_{F}$ maki . This also follows from
the Eilenberger equations eli . If $\alpha(z)$ and $\beta(z)$ vary slowly over
$\lambda_{F}$ in the case of $p_{s}(z)=(e\lambda H_{0}/c)e^{-z/\lambda}$, the
eigenstates of Eqs. (7)-(8) are in a narrow energy belt around the Fermi
surface where the BCS pairing occurs.

The dependencies of the surface energy levels on $(k_{x},k_{y})$ and $H_{0}$
were not investigated in Ref. corn , where Fig. 1 was presented as the only
evidence for such states. Since the well in $\Delta(z)$ calculated self-
consistently is likely much shallower than $\Delta(z)$ shown in Fig. 1, the
existence of these bound states is questionable. Even if such weakly-bound
stakes might occur at $H_{0}\simeq H_{c}$, they would get hybridized with the
itinerant states due to impurity scattering or surface topographic defects
because the energy gap $E_{g}(H_{0})=\Delta_{0}-v_{x}|p_{s}(t)|$ for itinerant
quasiparticles is lower than the bottom of a shallow well in $\Delta(z)$. For
instance, in the clean limit at $\kappa\gg 1$ and $T=0$ the gap
$E_{g}(H_{0})=\Delta_{0}(1-H_{0}/H_{g})$ decreases with $H_{0}$ and vanishes
at $H_{g}=(2/3)^{1/2}H_{c}$, while $\Delta$ remains constant eli . In any
case, the density of itinerant quasiparticles with $E_{g}<E<\Delta_{0}$ and
different orientations of ${\bf v}$ is much greater than the density of bound
states (if any) in which ${\bf v}$ is parallel to the surface, so the
contribution of bound states to rf losses is negligible.

To evaluate the surface resistance $R_{s}(H_{0})$ Deyo et al. used a ”three-
liquid” model in which the dissipative part of conductivity was taken in the
Drude-like form,

$\sigma_{1}=(n_{c}\tau_{c}+n_{b}\tau_{b})e^{2}/m,$ (10)

where $n_{b}$ and $n_{c}$ denote the densities of bound and continuum states,
and $\tau_{b}$ and $\tau_{c}$ are their respective relaxation times. Here
$n_{c}$ and $n_{b}$ were evaluated by summing up the quasiparticle states but
how it was done and what states were taken into account were not specified in
Ref. corn . Moreover, it was stated that $E<\Delta_{0}$ belong to bond states
and $E>\Delta_{0}$ to itinerant states, which contradicts the fact that the
energies $E_{g}(H_{0})<E<\Delta_{0}$ actually belong to continuum states.
Unlike the putative surface bound states, the itinerant quasiparticles with
$E_{g}(H)<E<\Delta_{0}$ do contribute significantly to $R_{s}$. The
phenomenological Eq. (10), which is not obtained from the theory of
electromagnetic response of superconductors mb ; tinkh ; kopnin , gives
incorrect temperature and field dependencies of $\sigma_{1}(T,H_{0})$. For
instance, at low fields Eq. (10) yields $\sigma_{1}\propto n_{c}\propto
T^{1/2}e^{-\Delta/T}$ inconsistent with the Mattis-Bardeen $\sigma_{1}\propto
T^{-1}\ln(T/\hbar\omega)e^{-\Delta/T}$ which has a logarithmic singularity at
$\omega\to 0$.

The problems with the field dependence of $\sigma_{1}(H_{0})$ given by Eq.
(10) are apparent if $\tau_{b}=\tau_{c}$ and the sum $n_{c}+n_{b}$ combines to
the total density of thermally-activated quasiparticles $n_{qp}(H_{0})$. Since
$n_{qp}(H_{0})$ increases with $H_{0}$ agp , Deyo et al. concluded that
itinerant quasiparticles cannot produce the decrease of $R_{s}(H_{0})$ with
$H_{0}$ observed on Nb cavities. However, this just reflects the inadequacy of
Eq. (10) because $\sigma_{1}(H_{0})$ decreasing with $H_{0}$ has been derived
from a microscopic theory ag ; kg . To obtain the decrease of $R_{s}(H_{0})$
with $H_{0}$ in their model Deyo et al. not only ascribed $E>\Delta_{0}$ and
$E<\Delta_{0}$ to continuum and bound states, but also assumed that these
states have different relaxation times $\tau_{c}$ and $\tau_{b}$. As a result,
$\sigma_{1}(H_{0})$ can decrease with $H_{0}$ if $\tau_{b}\ll\tau_{c}$ because
the density of quasiparticles with $E>\Delta_{0}$ and large $\tau_{c}$
decreases while the density of quasiparticles with $E<\Delta_{0}$ and small
$\tau_{b}$ increases as $H_{0}$ increases.

Assigning unequal $\tau_{b}$ and $\tau_{c}$ for itinerant states at
$E_{g}<E<\Delta_{0}$ and $E>\Delta_{0}$, and using various time constants
including inelastic electron-phonon (e-p) collision times $\tau_{ep}$ kaplan
as adjustable parameters in the phenomenological Eq. (10) in Ref. corn  has no
theoretical basis. Generally, $\sigma_{1}$ is determined by the rf power
absorbed by quasiparticles and is independent of $\tau_{ep}$ at low $\omega$
and $H_{0}$, whereas $\sigma_{1}$ at $H_{0}\ll H_{c}$ in the clean limit is
independent of any elastic or inelastic relaxation times mb . The e-p
collisions do control the energy transfer from quasiparticles to phonons,
producing a nonequilibrium correction $\delta f(E,t)$ to the quasiparticle
distribution function which becomes significant at very low $T$ or high
$\omega$ and $H_{0}$. Yet in the range of the parameters, where $\delta
f(E,t)$ is much smaller than the equilibrium
$f_{0}(E)=\left(e^{E/T}+1\right)^{-1}$, the e-p collisions do not affect
$R_{s}$ determined mostly by elastic impurity scattering and the effects of
current on the density of states and the coherence factors. These key effects
are not properly taken into account in Eq. (10).

A microscopic theory of nonlinear surface resistance based on the Keldysh
technique for nonequilibrium superconductors in the dirty limit at
$\lambda\gg\xi$ kopnin  was developed in Refs. ag ; kg . The results of Refs.
ag ; kg  were misinterpreted by Deyo et al. who claimed that they rely
crucially on the ”notion of maintaining $E_{2}=E_{1}+\hbar\omega$”,
”approximation to determine the nonequilibrium quasiparticle distribution
function” or ”The smeared density of states”, implying that the spatial
variation of screening current and temporal oscillations of the
superconducting gap, density of states and coherence factors were disregarded
in Refs. ag ; kg . Neither of these assertions is relevant, as one can see
from the formulas for $R_{s}$ in the case of equilibrium $f_{0}(E)$ ag ; kg :

$\displaystyle\\!\\!R_{s}=R_{1}\\!\\!\int_{0}^{\pi/\omega}\\!\\!\\!\\!dt\\!\int_{0}^{b_{0}}\\!db\\!\int_{E_{g}(t)}^{\infty}\frac{M[E,s(t)]dE}{\cosh(\frac{E}{2T})\cosh[\frac{(E+\hbar\omega)}{2T}]},$
(11) $\displaystyle M=\mbox{Re}G[\hbar\omega+E,s(t)]\mbox{Re}G[E,s(t)]+$
$\displaystyle\mbox{Re}F[\hbar\omega+E,s(t)]\mbox{Re}F[E,s(t)],$ (12)
$\displaystyle s=b\sin^{2}\omega t,\qquad
E_{g}(t)=[\Delta^{2/3}(t)-s^{2/3}(t)]^{3/2},$ (13)

where
$R_{1}=\mu_{0}^{2}\omega^{2}\lambda^{3}\sinh(\hbar\omega/2T)/2\pi\hbar\rho_{n}$,
$E_{g}(t)$ is the quasiparticle gap, $b=b_{0}e^{-2z/\lambda}$,
$b_{0}=(H_{0}/2H_{c})^{2}\Delta_{0}$, $H_{c}=\phi_{0}/2^{3/2}\xi\lambda$,
$\Delta(t)=\Delta_{0}-\pi s(t)/4$, and the Green functions $G^{R}(E,t)$ and
$F^{R}(E,t)$ were obtained by solving the time-dependent Usadel equation at
$\hbar\omega\ll\Delta$ ag . The energy gap $E_{g}(z,t)$ is always lower than
the bottom of the well in $\Delta(z,t)$. Here $E_{g}(t)$, $\Delta(t)$ and the
spectral function $M(E,t)$ which accounts for the density of states and the
coherence factors, oscillate under rf field so no ”smeared density of states”
is involved here. The integration over $b$ accounts for the spatial change of
$J(z,t)$ with the distance from the surface ag . The arguments $E$ and
$E+\hbar\omega$ in $M(E,t)$ have nothing to do with constant quasi-energies
but result from a straightforward calculation of the rf power
$R_{s}H_{0}^{2}/2=\int_{0}^{\infty}\langle J(z,t)\epsilon(z,t)\rangle dz$,
where $\langle...\rangle$ denotes time averaging over the rf period,
$\epsilon(z,t)=(H_{0}\lambda\omega/c)e^{-z/\lambda}\cos\omega t$ is the rf
electric field and the superconducting current $J(z,t)$ kopnin  was calculated
using analytical solutions for $G^{R}(E,t)$ and $F^{R}(E,t)$ given in Refs. ag
; kg . Here $E$ just labels the quasiparticle states in the time-dependent
spectral function $M(E,t)$.

Figure 2:  The surface resistance calculated from Eqs. (11)-(13) for Nb at
$\Delta_{0}=18$ K, $T=2K$ and 1.3 GHz.

Shown in Fig. 2 is an example of $R_{s}(H_{0})$ calculated from Eqs. (11)-(13)
which captures the nonmonotonic field dependence of $R_{s}(H_{0})$ observed on
Nb cavities raise ; ag ; kg . Such generic behavior of $R_{s}(H_{0})$ occurs
even in the case of equilibrium $f_{0}(E)$ of itinerant quasiparticles driven
by a low-frequency rf current at $\omega\tau_{ep}\lesssim 1$, although a
frozen density of nonequilibrium quasiparticles at $\omega\tau_{ep}\gg 1$ may
cause even stronger decrease of $R_{s}(H_{0})$ with $H_{0}$ ag . Developing a
microscopic theory of $R_{s}(H_{0})$ of a nonequlibrium superconductor at low
$T$ and high $\omega$ is a challenging outstanding problem.

In conclusion, the surface bound states shown in Fig. 1 and their contribution
to the nonlinear electromagnetic response are questionable. The
phenomenological Eq. (10) cannot serve as a theoretical basis for the
calculations of $R_{s}(H_{0})$, particularly the field-induced reduction of
$R_{s}(H_{0})$ which has been observed on Nb cavities.

This work was supported by DOE under Grant DE-SC0010081-020.

## References

  * (1) S. Deyo, M. Kelley, N. Sitaraman, T. Oseroff, D.B. Liarte, T. Arias, M. Liepe, and J.P. Sethna, Phys. Rev. B106, 104502 (2022).
  * (2) M. Tinkham Superconductivity (Mc Gray Hill Inc., New York, 2nd ed. 1996).
  * (3) P. Dhakal, Physics Open 5 100034 (2020) and the references therein.
  * (4) K. Maki, In Superconductivity (ed. R. D. Parks, v 2. p. 1035, 1968).
  * (5) J. Bardeen, Rev. Mod. Phys. 34, 667 (1962).
  * (6) N. B. Kopnin, Theory of nonequilibrium superconductivity. (Oxford Univ. Press, 2001)
  * (7) G. Catelani and J. P. Sethna, Phys. Rev. B78, 224509 (2008); P-J. Lin and A. Gurevich, Phys. Rev. B85, 054513 (2012).
  * (8) D. C. Mattis and J. Bardeen, Phys. Rev. 111, 412 (1958).
  * (9) S.B. Kaplan, C.C. Chi, D.N. Langenberg, J.J. Chang, S. Jafarey, and D.J. Scalapino, Phys. Rev. B14, 4854 (1976).
  * (10) A. Gurevich, Physica C 441, 38 (2006).
  * (11) A. Gurevich, Phys. Rev. Lett. 113, 087001 (2014); Supercond. Sci. Technol. 30, 034004 (2017).
  * (12) T. Kubo and A. Gurevich, Phys. Rev. B100, 064522 (2019).


# Feed-forward Disturbance Compensation for Station Keeping in Wave-dominated
Environments
††thanks: This work was supported by the EPSRC under grant No. EP/R513209/1.

Kyle L. Walker1*  —, Adam A. Stokes1*  —, Aristides Kiprakis2*  — and
Francesco Giorgio-Serchi13*  — 1Institute for Integrated Micro and Nano
Systems  2Institute for Energy Systems
School of Engineering
University of Edinburgh
Edinburgh, UK
3Correspondence<EMAIL_ADDRESS>
###### Abstract

When deploying robots in shallow ocean waters, wave disturbances can be
significant, highly dynamic and pose problems when operating near structures;
this is a key limitation of current control strategies, restricting the range
of conditions in which subsea vehicles can be deployed. To improve dynamic
control and offer a higher level of robustness, this work proposes a Cascaded
Proportional-Derivative (C-PD) with Feed-forward (FF) control scheme for
disturbance mitigation, exploring the concept of explicitly using disturbance
estimations to counteract state perturbations. Results demonstrate that the
proposed controller is capable of higher performance in contrast to a standard
C-PD controller, with an average reduction of $\approx 48\%$ witnessed across
various sea states. Additional analysis also investigated performance when
considering coarse estimations featuring inaccuracies; average improvements of
$\approx 17\%$% demonstrate the effectiveness of the proposed strategy to
handle these uncertainties. The proposal in this work shows promise for
improved control without a drastic increase in required computing power; if
coupled with sufficient sensors, state estimation techniques and prediction
algorithms, utilising feed-forward compensating control actions offers a
potential solution to improve vehicle control under wave-induced disturbances.

###### Index Terms:

Feed-forward Control, Disturbance Compensation, State Estimation, Dynamic
Control, Underwater Vehicles.

##  I Introduction

Advanced control of marine vehicles is sharply becoming an industrial
necessity rather than an academic exercise, with the offshore sector seeking
higher levels of autonomy with respect to intervention tasks, inspection tasks
and similar [1]. An industry section where increased autonomy would be highly
beneficial is offshore, in particular the marine renewable sector, as the
transition to cleaner energy sources begins to accelerate and harsher
environments become the subject of exploration for power generation [2].
Remotely Operated Vehicles have become a solidified aspect of many subsea
procedures, but operation near the free-surface or in shallower waters remains
a challenge when the sea state is not calm, owing to the significant influence
of surface waves which often limits the deployment of piloted strategies [3].
Specifically with respect to marine renewable devices, situation within a
turbulent environment is critical to generate sufficient power; thus, a
greater level of robustness, reliability and precision is required with
regards to vehicle control if autonomous operation is to become a reality [4].

Attempts to develop automatic disturbance compensation control algorithms have
varied, with early implementations investigating vision-based solutions [5,
6]. In calmer environments these strategies can be effective, but there is a
lack of robustness due to the dependence on visibility and being feedback
based. To tackle this, alternative solutions have suggested inherently robust
methods in the form of sliding mode and adaptive controllers [7, 8], which
account for hydrodynamic loading by considering this as a set of generalised
system disturbances. For small-magnitude disturbances this can improve
performance sufficiently, but a question remains over stability and
performance guarantees when loading increases, for example low-depth operation
under large wave heights.

Handling time-varying and unsteady disturbances remains an open challenge, but
a proposition which holds promise lies in exploiting forecasting methods [9,
10]; for wave-induced disturbances which are largely predictable in nature,
utilising preview information can assist in disturbance compensation [11].
Given that wave predictions can be deduced through time-history data, this
lends itself to being applicable to a wider range of disturbances than other
proposals. Furthermore, the development of state estimation methods [12]
coupled with a hydrodynamic loading model [14, 13] allows the formulation of a
feed-forward (FF) control action to reduce state error, in conjunction with an
establish feedback controller for set-point regulation. Likewise, the
formulation of a single additional control action ensures computation
requirements remain low, a major benefit in relation to alternative predictive
solutions [15] for these highly dynamic scenarios.

With respect to the above, this work proposes the use of a Feed-Forward (FF)
control element which models the wave disturbances as a product of added
inertia and hydrodynamic drag, coupling this with a Cascaded Position-Velocity
Proportional-Derivative Controller (C-PD) to counteract wave disturbances. The
vehicle state is estimated through an Extended Kalman Filter (EKF), which
determines the corrective control action in the feedback loop. The controller
is simulated under several wave conditions when the disturbance estimation is
both well-known and features significant uncertainty, showing that the
proposed FF+C-PD controller outperforms a standard C-PD controller in all
cases. These findings demonstrate that even coarse predictions of
environmental disturbances can prove useful in mitigating state error of
subsea vehicles, thus improving the accuracy of dynamic control tasks.

##  II Modelling

Throughout this work, the vehicle is considered to possess 3DoF and is
restricted to a planar case; therefore the presented model considers the
surge, heave and pitch motions only. Analogously, the ocean waves are modelled
according to 2nd-order planar theory and are assumed to propagate uni-
directionally.

###  II-A Vehicle Rigid-Body Dynamics

The kinematics of the rigid-body are described by considering two co-ordinate
frames; the earth-fixed and body-fixed. As depicted in Fig. 1, these are
related by a transformation according to [18]:

$\dot{\boldsymbol{\eta}}=\mathbf{J(\boldsymbol{\eta})}\boldsymbol{\nu}$ (1)

where $\boldsymbol{\eta}$ is a state vector describing the position and
orientation of the vehicle, $\boldsymbol{\nu}$ is a state vector of linear and
angular velocities and $\mathbf{J}\in\mathbb{R}^{3\times 3}$ is the
transformation matrix relating the two frames.

Given the above kinematic representation, the vehicle dynamics exhibit
nonlinear behaviour and are defined according to:

$\mathbf{M\dot{\boldsymbol{\nu}}}+\mathbf{C(\boldsymbol{\nu})\boldsymbol{\nu}}+\mathbf{D(\boldsymbol{\nu})\boldsymbol{\nu}}+\mathbf{g(\boldsymbol{\eta})}=\boldsymbol{\tau}+\boldsymbol{\tau_{E}}$
(2)

where $\mathbf{M}=\mathbf{M}_{RB}+\mathbf{M}_{A}\in\mathbb{R}^{3\times 3}$ is
an inertia matrix,
$\mathbf{C(\boldsymbol{\nu})}=\mathbf{C_{RB}(\boldsymbol{\nu})}+\mathbf{C_{A}(\boldsymbol{\nu})}\in\mathbb{R}^{3\times
3}$ is a matrix of Coriolis and centripetal terms,
$\mathbf{D(\boldsymbol{\nu})}=\mathbf{D}_{L}(\boldsymbol{\nu})+\mathbf{D}_{Q}(\boldsymbol{\nu})\odot|\boldsymbol{\nu}|\in\mathbb{R}^{3\times
3}$ is a hydrodynamic damping matrix and
$\mathbf{g(\boldsymbol{\eta})}\in\mathbb{R}^{3}$ is a vector of hydrostatic
restoring forces. In the above, subscripts RB and A relate to contributions
from rigid body and added inertial effects, with:

$\mathbf{M}_{RB}=\begin{bmatrix}m&0&0\\\ 0&m&0\\\
0&0&I_{y}\end{bmatrix}\quad\mathbf{M}_{A}=\begin{bmatrix}X_{\dot{u}}&0&X_{\dot{q}}\\\
0&Z_{\dot{w}}&0\\\ M_{\dot{u}}&0&M_{\dot{q}}\end{bmatrix}$ (3)

where $m$ is the vehicle dry mass, $I_{y}$ is a rotational inertia and
$X_{\dot{q}}$, $Z_{\dot{w}}$, and $X_{\dot{u}}=M_{\dot{u}}$ are added mass
coefficients. The Coriolis and centripetal terms are derived in accordance
with Eq. 3 [18] and the damping matrix is specified as:

$\mathbf{D}_{L}=\begin{bmatrix}X_{u}&0&0\\\ 0&Z_{w}&0\\\
0&0&M_{q}\end{bmatrix}\quad\mathbf{D}_{Q}=\begin{bmatrix}X_{u|u|}&0&0\\\
0&Z_{w|w|}&0\\\ 0&0&M_{q|q|}\end{bmatrix}$ (4)

where $X_{u}$, $Z_{w}$ and $M_{q}$ are linear damping coefficients whilst
$X_{u|u|}$, $Z_{w|w|}$ and $M_{q|q|}$ are quadratic damping coefficients. All
hydrodynamic parameters are defined in Table II.

Finally, $\boldsymbol{\tau}\in\mathbb{R}^{3}$ is a vector of control forces
and moments whilst environmental disturbances are lumped within the vector
$\boldsymbol{\tau_{E}}\in\mathbb{R}^{3}$. As this work is concerned with
mitigating the effects of surface waves on vehicle behaviour, disturbances
arising from ocean currents are assumed negligible and $\boldsymbol{\tau_{E}}$
is formulated to purely describe wave-induced loading, according the model
presented in the following section.

Figure 1: The two frames, earth-fixed and body-fixed, and the relating
transformation for the BlueROV2 Heavy configuration.

###  II-B Wave-Induced Disturbances

Figure 2: Block diagram of the proposed feed-forward disturbance mitigation
technique, with the grey blocks representing the generation of the additional
compensating control actions.

To form a temporal history of the wave elevation, a 2nd-order model was
adopted utilising the principle of superposition. The sea state is therefore
described by a spectrum of $N$ monochromatic components, each with a unique
wave amplitude $A$, wave period, $T$, and phase offset $\epsilon$. It follows
that at a point $x$ and at time $t$ the wave elevation is described by [18]:

$\displaystyle\zeta(x,t)=$
$\displaystyle\sum_{i=0}^{N}A_{i}cos(\kappa_{i}x-\omega_{i}t+\epsilon_{i})$
$\displaystyle+\sum_{i=1}^{N}\frac{1}{2}\kappa_{i}A_{i}^{2}cos2(\kappa_{i}x-\omega_{i}t+\epsilon_{i})$
(5)

where $\kappa$ and $\omega$ are the wave number and circular frequency
respectively. These additional parameters are deduced by solving the
dispersion relation:

$\omega_{i}=\sqrt{g\kappa_{i}\tanh{\kappa_{i}d}}$ (6)

where $g$ and $d$ are the gravitational constant and approximate seabed depth.
In Eq. 5, the wave amplitude is related to the spectral density function
$S(\omega)$ by $A^{2}=2S(\omega)\Delta\omega$ where $\Delta\omega$ is the
difference between successive frequencies. Several models exist to describe
the spectral density function [18]; in this work the JONSWAP spectrum is
adopted owing to the relation with the North Sea, an area of interest for this
application. The spectral density function is therefore:

$S(\omega)=\frac{\alpha
g^{2}}{\omega^{5}}\exp\left[-\frac{5}{4}\left(\frac{\omega_{p}}{\omega}\right)^{4}\right]\gamma^{\Gamma}$
(7)

where $\omega_{p}$ is the spectral peak frequency, $\gamma^{\Gamma}$ is a peak
enhancement factor and

$\displaystyle\alpha=0.0081$
$\displaystyle\sigma=\begin{cases}0.07,&\text{if}\ \omega\leq\omega_{p}\\\
0.09,&\text{if}\ \omega\geq\omega_{p}\end{cases}$
$\displaystyle\Gamma=\exp\left[\frac{(\omega-\omega_{p})^{2}}{2\omega_{p}^{2}\sigma^{2}}\right]$

Spectral information can also be exploited to deduce the fluid particle
motions at a point beneath the surface $z$, which for 2nd-order theory
produces [19]:

$\displaystyle u_{p}(x,z,t)=$
$\displaystyle\sum_{i=1}^{N}\frac{gH_{i}}{2c}\frac{\cosh\kappa_{i}(z+d)}{\cosh\kappa_{i}d}\cos(\kappa_{i}x-\omega_{i}t+\epsilon_{i})+$
$\displaystyle\frac{3}{16}c\kappa_{i}^{2}H_{i}^{2}\frac{\cosh[2\kappa_{i}(z+d)]}{\sinh^{4}\kappa_{i}d}\cos[2(\kappa_{i}x-\omega_{i}t+\epsilon_{i})]$
(8) $\displaystyle w_{p}(x,z,t)=$
$\displaystyle\sum_{i=1}^{N}\frac{gH_{i}}{2c}\frac{\sinh\kappa_{i}(z+d)}{\cosh\kappa_{i}d}\sin(\kappa_{i}x-\omega_{i}t+\epsilon_{i})+$
$\displaystyle\frac{3}{16}c\kappa_{i}^{2}H_{i}^{2}\frac{\sinh[2\kappa_{i}(z+d)]}{\sinh^{4}\kappa_{i}d}\sin[2(\kappa_{i}x-\omega_{i}t+\epsilon_{i})]$
(9)

where $c$ is celerity, defined according to seabed depth to wavelength ratio
[20]. This facilitates deduction of the wave-induced hydrodynamic loads
($X_{E}$, $Z_{E}$ and $M_{E}$ for the surge, heave and pitch respectively)
acting on the body; here, a low-order model is employed which has been
experimentally validated in previous work [13, 14, 22, 23]:

$\boldsymbol{\tau}_{E}=\begin{bmatrix}X_{E}\\\ Z_{E}\\\
M_{E}\end{bmatrix}=\begin{bmatrix}X_{\dot{u}}\dot{\nu}_{p,x}+\\{X_{u}+X_{u|u|}|\nu_{p,x}|\\}\nu_{p,x}\\\
Z_{\dot{w}}\dot{\nu}_{p,z}+\\{Z_{w}+Z_{w|w|}|\nu_{p,z}|\\}\nu_{p,z}\\\
\int_{-L/2}^{L/2}Z_{E}(x^{\prime},z^{\prime},t)x^{\prime}dx\end{bmatrix}$ (10)

where
$\boldsymbol{\nu}_{p}=\left[\nu_{p,x},\nu_{p,z}\right]^{T}=\mathbf{R}_{y}(\theta)\left[u_{p},w_{p}\right]^{T}$
($\mathbf{R}_{y}(\theta)$ is a rotation matrix). Also, $L$ is the vehicle body
length and $(x^{\prime},z^{\prime})$ refers to points along the vehicle axial
length within the local frame.

##  III Control Methodology

Obtaining spectral knowledge of the immediate ocean environment around the
vehicle is key in deducing the magnitude of disturbances when applying the
model described by Eq. 10. Several methods have been proposed for predicting
impending waves and wave loads, including but not limited to _in-situ_ sensor
fusion [16], auto-regressive models [9] and deterministic methods [10], the
latter in particular focusing on utilising spectral information. These have
been applied successfully in the context of wave energy converters, therefore
an analogous deployment with respect to an underwater vehicle holds adjacent
potential. It is therefore postulated here that exploiting knowledge of wave
disturbances to formulate a feed-forward control action can assist in reducing
state perturbations, improving station keeping performance and widening the
range of deployable conditions. The control law is therefore formulated as a
Cascaded Proportional-Derivative (C-PD) with feed-forward (FF) action:

$\displaystyle\boldsymbol{\tau}=$
$\displaystyle\underbrace{\boldsymbol{\tau}_{max}\left[\mathbf{K}_{p,v}\\{\boldsymbol{\nu}-\left(\mathbf{K}_{p}\mathbf{e}+\mathbf{K}_{d}\mathbf{\dot{e}}\right)\\}\right]}_{C-PD}$
(11)
$\displaystyle+\underbrace{\\{\mathbf{M}_{A}\boldsymbol{\dot{\nu}}_{p}+\mathbf{D}(\boldsymbol{\nu_{p}})(-\boldsymbol{\nu}_{p})\\}}_{FF}$

where $\mathbf{K}_{p}$ and $\mathbf{K}_{d}$ are the position PD gains,
$\mathbf{K}_{p,v}$ is the velocity P-gain and $\mathbf{e}$ is the state error.
Also, $\boldsymbol{\tau}_{max}\in\mathbb{R}^{3}$ is a vector describing the
maximum torque available in each DoF. As the heave and pitch states are
controlled via the horizontal thrusters, a thrust allocation algorithm is
embedded within the control architecture to generate the appropriate control
inputs to deliver the required forces and moments. Inclusion of the motor
time-delay as a first-order response yields:

$\boldsymbol{\mu}=\left[\left(1-e^{-\Delta
t/t_{m}}\right)\mathbf{K}_{\boldsymbol{\tau}}^{-1}\right]\boldsymbol{B}_{\boldsymbol{\mu}}^{\dagger}\boldsymbol{\tau}$
(12)

where $\boldsymbol{B}_{\boldsymbol{\mu}}^{\dagger}$ is the Moore-Penrose
pseudo-inverse of the thrust allocation matrix
$\boldsymbol{B}_{\boldsymbol{\mu}}$, $\mathbf{K}_{\boldsymbol{\tau}}$ is a
force co-efficient matrix and $\Delta t_{m}$ is the motor-time constant. Eq.
12 produces a solution $\boldsymbol{\mu}\in\mathbb{R}^{8}$ to be allocated to
each thruster.

(a)

(b)

Figure 3: Map showing the (a) various locations of buoys around Scotland,
circling the buoy selected for this analysis and (b) an enlarged map of the
Moray Firth buoy location.

###  III-A State Estimation

Deploying automatic control in practice requires knowledge of the vehicles
location to be known or at minimum a reasonable estimation to be inferred; to
add an additional layer of realism to the simulations, we employ an EKF as a
state estimator to track the vehicle state. The EKF algorithm has two update
phases; the _predictor_ phase and the _corrector_ phase. The predictor phase
considers the initial estimates and projects the error covariance and state,
$\mathbf{P}$ and $\mathbf{\hat{x}}$, ahead in time such that:

$\mathbf{P}_{k}^{-}=\mathbf{A}_{k}\mathbf{P}_{k-1}\mathbf{A}_{k}^{T}+\mathbf{Q}_{k-1}$
(13)
$\mathbf{\hat{x}}_{k}^{-}=F(\hat{\mathbf{x}}_{k-1},\boldsymbol{\mu_{k-1}},\mathbf{w})$
(14)

where $\mathbf{A}=\frac{dF}{dx}\rvert_{\mathbf{x}=\mathbf{\hat{x}}}$ is the
linearised state transition matrix, $\mathbf{Q}$ is the process error
covariance and $\mathbf{w}$ represents the process noise. Here, $F$ represents
a nonlinear function for the vehicle dynamics. From this, the corrector phase
proceeds to compute the Kalman gain:

$\mathbf{K}_{k}=\mathbf{P}_{k}^{-}\mathbf{H}_{k}^{T}(\mathbf{H}_{k}\mathbf{P}_{k}^{-}\mathbf{H}_{k}^{T}+\mathbf{R}_{k})^{-1}$
(15)

where $\mathbf{H}_{k}$ is a positional measurement and $\mathbf{R}_{k}$ is the
measurement error covariance, before updating the estimate with a measurement
$y_{k}$

$\mathbf{\hat{x}}_{k}=\mathbf{\hat{x}}_{k}^{-}+\mathbf{K}_{k}(y_{k}-\mathbf{H}_{k}\mathbf{\hat{x}}_{k}^{-})$
(16)

Similarly, the error covariance is also corrected in this stage before looping
back to Eq. 13 and repeating for every time-step $k$, where $\mathbf{I}$ is an
identity matrix:

$\mathbf{P}_{k}=(\mathbf{I}-\mathbf{K}_{k}\mathbf{H}_{k})\mathbf{P}_{k}^{-}$
(17)

TABLE I: Case assignments and parameters for the analysed wave spectra.

.  Case Reference Peak Period (s) Significant Wave Height (m) W1 7.1 2.78 W2
9.5 3.47 W3 11.1 3.24

TABLE II: BlueROV2 Heavy dimensions and hydrodynamic parameters; data based on [24, 25, 26]. Parameter | Nomenclature | Value
---|---|---
Weight | $W$ | 112.8 N
Buoyancy | $B$ | 114.8 N
Rotational Inertia, $y$ | $I_{y}$ | 0.253 kgm2
Added Inertia Coeff. |  $X_{\dot{u}}$, $Z_{\dot{w}}$ | 6.36, 18.68 kg
” | $M_{\dot{q}}$ | 0.135 kgm2
” |  $X_{\dot{q}}$, $M_{\dot{u}}$ | 0.67 kgm
Linear Drag Coeff. |  $X_{u}$, $Z_{w}$ | 13.7, 33 kg/s
” | $M_{q}$ | 0.80 kgm2/s
Quadratic Drag Coeff. |  $X_{u|u|}$, $Z_{w|w|}$ | 141, 190 Ns2/m2
” | $M_{q|q|}$ | 0.47 Nms2
Centre of Buoyancy | $r_{B}$ | [0, 0, 0.028]m
Maximum Thrust | $T_{max}$ | 35 N
Thruster Offset | $\alpha$ | 45o

##  IV Scenario Configuration

Given the intended application of the proposed disturbance mitigation method
is for improved performance during inspection and maintenance of
devices/structures in wave-dominated environments, the vehicle was simulated
under three different sea states to analyse performance relative to varying
wave parameters. Spectral data was sourced from the online repository of the
Centre for Environmental Fishes and Aquaculture Science (Cefas) [27],
collected by a wave buoy situated off the coast of Inverness, Scotland in the
Moray Firth. The location of the buoy is shown in Fig. 3 where $d=54m$;
offshore wind farms are typically located within areas of this depth [28, 29],
thus emulating the conditions of a typical inspection or maintenance task. The
selected wave fields and assigned case references are given in Table I with
all vehicle parameters given in Table II.

The controller was tasked with performing station-keeping at a depth of $z=5$m
and for a temporal segment of $600$s with a resolution of $\Delta t=0.05s$,
exposing the vehicle to significant magnitude wave disturbances for a
prolonged period of time. It should be noted here that station-keeping refers
to both positional and attitude regulation, during which the controller
attempts to maintain a reference set-point
$\mathbf{x}_{r}=[x_{r},z_{r},\theta_{r}]^{T}$. Throughout all analysis, sensor
noise is considered which is mitigated by the inclusion of an EKF to estimate
the vehicle state. Comparisons are drawn between the FF controller and a
standard C-PD controller as a baseline reference which also exploits the EKF
to monitor state error.

Figure 4: Temporal segment of the (a) wave in case W3. showing the (b) surge
position, (c) heave position and (d) pitch attitude.

##  V Simulation Results

Performance of the proposed control scheme was analysed by considering the
Root-Mean-Square-Error (RMSE) and absolute maximum error witnessed across the
600s simulation; the power consumed during the station keeping mission was
also recorded and analysed.

###  V-A Station Keeping Performance

A temporal segment for case W3 is displayed in Fig. 4 which evidently shows
the FF compensation having a positive impact on the station keeping error,
relative to the C-PD controller. The preview knowledge of the wave is seen to
reduce RMSE by up to $58.6\%$, with this specific value relating to the pitch
motion for case W1. For the surge and heave, the maximum reductions were also
related to case W1 which returned a $47.5\%$ and $48.2\%$ improvement
respectively. This is an indicator that for lower peak period spectra the
inclusion of FF compensation is more critical, however this could be owed to
the causal effect of larger pitch perturbations in contrast to cases W1 and
W2.

In terms of the maximum displacement recorded during the simulations, the
greatest improvement was also related to the pitch, reducing attitude error by
$50.1\%$ on average. This was much lower for the other DoF, with the surge and
heave only showing $22\%$ and $11.8\%$ reductions; it is suspected that this
is due to a brief section of the simulated waves subjecting the vehicle to
sharp, high magnitude disturbances which the control was unable to correct
for. Behaviour similar to this can be seen in Fig. 4 at $\approx 315$s, where
the traces for the control scheme undergo similar magnitude displacements for
the largest wave height in the segment. Across all cases there was a mean
reduction in RMSE of $\approx 48\%$ and $\approx 28\%$ in maximum error; all
absolute values are shown in Fig. 6.

###  V-B Sensitivity to Noise

Each case was analysed when both the disturbance feed-forward term is deemed
to be accurate (denoted SNR0) and feature imprecisions (denoted SNR15); the
latter attempts to provide insight into controller performance when the
vehicle encounters disturbances that differ from those anticipated by the FF
controller. To achieve this, spectral noise with a SNR of $15$ was injected
directly to the spectral component amplitude and phase offset to alter the
wave (and thus FF compensation calculations) significantly and randomly as
shown in Fig. 5. These results are also compiled within Fig. 6.

When the disturbance is considered to be inaccurate, there is still
improvement in station keeping accuracy across all cases with respect to RMSE
and the majority of cases with respect to maximum error; only case W3 returns
higher maximum error in the surge and heave and in these specific instances
the difference is marginal with increases of $<1\%$ and $5\%$ respectively.
The average reduction in RMSE of $17\%$ supports the claim that these
instances are at isolated periods in time and not a regular occurrence;
similar to above, the pitch experienced the highest improvement in error of
$14.2\%$. Overall these results clearly show that even utilising a spectral
estimation that features inaccuracies can offer a noticeable improvement in
state regulation, implying that well established methods can be applied with
confidence at this end of the control pipeline.

Figure 5: Temporal segment for case W3 during noise analysis; here, SNR0 is
the encountered wave by the vehicle where as SNR15 is the expected wave by the
controller used to generate compensating control actions.  Figure 6: State
error for each control method during the station keeping simulation, showing
the (a)(c)(e) RMSE and (b)(d)(f) absolute maximum error recorded.

###  V-C Station Keeping Efficiency

Given the improvement in performance and additional control action generated
by the FF component, it was anticipated additional power would be required. To
confirm this, the power consumed during the mission was modelled according to
manufacturer data such that [21]:

$\mathbf{P}=0.0011\boldsymbol{\tau}^{3}+0.02078\boldsymbol{\tau}^{2}+0.297\boldsymbol{\tau}$
(18)

where $\mathbf{P}\in\mathbb{R}^{3}$ is the power consumed by the control
actions in each DoF. This was modelled according to a nominal operating volage
of 16V and the results are displayed in Fig. 7.

The power does increase when employing the FF scheme, but interestingly less
power was consumed when considering noisy disturbance estimations. It is
likely that this is attributed to the FF controller anticipating lower
magnitude disturbances for the majority of the simulated mission, which would
explain the increase in error that is shown in Fig. 6. Similarly, the segment
displayed in Fig. 5 demonstrates this behaviour at various instances, for
example at 120s where the disturbance is significantly underestimated. This
supports the connected data in relation to the positional and attitude
RMSE/maximum error. For large magnitude waves it can be argued that this
increased power expenditure is utilised well to reduce station keeping error
significantly, thus in situations where this is paramount it becomes less of a
burdening factor.

Figure 7: Power consumed during the station keeping mission for each case and
control method.

##  VI Conclusions

This paper has proposed active wave-induced disturbance rejection for
underwater vehicles based on the inclusion of a feed-forward control action.
An extensive simulation study demonstrated the capability of the proposed
control to significantly reduce state regulation error in both position and
attitude, highlighting performance improvements of up to $56.8\%$ over three
sea states with varying parameters. Similarly, the maximum error (which is
arguably the critical factor in these scenarios) was also reduced
substantially, with more prominent reductions witnessed for the vehicle pitch.
Given the availability of wave prediction tools, these results provide
evidence of the real potential related to the incorporation of modelled wave-
induced loads directly within the control to improve performance. Considering
this, an interesting avenue for exploration would be to develop a preview of
wave-induced disturbances along a future time-horizon; this would facilitate
the development of different forms of model predictive control, potentially
improving performance even further by evaluating an optimal control sequence,
rather than a one step control action.

## References

  * [1] H. Hastie, K. Lohan, M. Chantler, D. A. Robb, S. Ramamoorthy, R. Petrick, S. Vijayakumar and David Lane, ”The ORCA Hub: Explainable Offshore Robotics through Intelligent Interfaces”, arXiv preprint, 2018.
  * [2] E. Zereik, M. Bibuli, N. Mišković, P. Ridao and A. Pascoal, ”Challenges and future trends in marine robotics”, _Annual Reviews in Control_ , vol. 46, pp. 350-368, 2018.
  * [3] O. Khalid, G. Hao, C. Desmond, H. Macdonald, F. D. McAuliffe, G. Dooly and Weifei Hu, ”Applications of robotics in floating offshore wind farm operations and maintenance: Literature review and trends”, _Wind Energy_ , vol. 25, no. 11, pp. 1880-1899, 2022.
  * [4] S. Sivčev, E. Omerdić, G. Dooly, J. Coleman and Toal, D., ”Towards inspection of marine energy devices using ROVs: Floating wind turbine motion replication”, _Advances in Intelligent Systems and Computing_ , vol 693, 2018.
  * [5] K. N. Leabourne, S. M. Rock, S. D. Fleischer and R. Burton, ”Station keeping of an ROV using vision technology,” Oceans ’97 MTS/IEEE Conference Proceedings, Halifax, NS, Canada, 06-09 Oct. 1997, pp. 634-640.
  * [6] R. L. Marks, H. H. Wang, M. J. Lee and S. M. Rock, ”Automatic visual station keeping of an underwater robot,” Proceedings of OCEANS’94, Brest, France, 13-16 Sept. 1994, pp. 137-142.
  * [7] R. Cui, X. Zhang and D. Cui, ”Adaptive sliding-mode attitude control for autonomous underwater vehicles with input nonlinearities,” Ocean Engineering, vol. 123, 2016, pp. 45-54.
  * [8] T. Elmokadem, M. Zribi and K. Youcef-Toumi, ”Terminal sliding mode control for the trajectory tracking of underactuated Autonomous Underwater Vehicles,”, Ocean Engineering, vol. 129, 2017, pp. 613-625.
  * [9] F. Fusco and J. Ringwood, ”A Study on Short-Term Sea Profile Prediction for Wave Energy Applications”, _Proc. of the 8th European Wave and Tidal Energy Conference_ , Uppsala, Sweden, pp. 756-765, 2009.
  * [10] M. R. Belmont, J. Christmas, J. Dannenberg, T. Hilmer, J. Duncan, J. M. Duncan and B. Ferrier, ”An examination of the feasibility of linear deterministic sea wave prediction in multidirectional seas using wave profiling radar: Theory, simulation, and sea trials”, _Journal of Atmospheric and Oceanic Technology_ , vol. 31, no. 7, pp. 1601-1614, 2014.
  * [11] X. Fang and W.-H. Chen, ”Model Predictive Control With Preview: Recursive Feasibility and Stability”, _IEEE Control Systems Letters_ , vol. 6, pp. 2647-2652, 2022.
  * [12] S. Soylu, A. A. Proctor, R. P. Podhorodeski, C. Bradley and B. J. Buckham, ”Precise trajectory control for an inspection class ROV”, _Ocean Engineering_ , vol. 111, , pp. 508-523, 01 Jan. 2016.
  * [13] K. L. Walker, R. Gabl, S. Aracri, Y. Cao, A. A. Stokes, A. Kiprakis and F. Giorgio-Serchi, ”Experimental Validation of Wave Induced Disturbances for Predictive Station Keeping of a Remotely Operated Vehicle”, _IEEE Robotics and Automation Letters_ , vol. 6, no. 3, pp. 5421-5428, July 2021.
  * [14] K. L. Walker, R. Gabl, S. Aracri, Y. Cao, A. A. Stokes, A. Kiprakis and F. Giorgio-Serchi, ”Experimental Validation of Unsteady Wave Induced Loads on a Stationary Remotely Operated Vehicle”, _Proc. of IEEE International Conference on Robotics and Automation (ICRA)_ , Xi’an, China, pp. 2242-2248, 30 May - 05 June 2021.
  * [15] D. C. Fernández, and Geoffrey A. Hollinger, ”Model Predictive Control for Underwater Robots in Ocean Waves,” IEEE Robotics and Automation Letters, vol. 2, no. 1, pp. 88-95, 2016.
  * [16] J. M. Selvakumar and T. Asokan, ”Station keeping control of underwater robots using disturbance force measurements”, _Journal of Marine Science and Technology_ , vol. 21, pp. 70-85, 2016.
  * [17] C. Long, X. Qin, Y. Bian and M. Hu, ”Trajectory tracking control of ROVs considering external disturbances and measurement noises using ESKF-based MPC”, _Ocean Engineering_ , vol. 241, 2021.
  * [18] T. I. Fossen, ”Guidance and Control of Ocean Vehicles”, _Wiley_ , Chichester, UK, 1994.
  * [19] M. E. McCormick, ”Ocean Engineering Mechanics with Applications”, _Cambridge University Press_ , New York, USA, 2010.
  * [20] D. Reeve, A. chadwick and C. Fleming, ”Coastal Engineering: Processes, Theory and Design Practice”, _Spon Press_ , London, UK, 2004.
  * [21] “T200 Thruster Datasheet,” Tech. Rep. [Online]. Available: http://docs.bluerobotics.com/thrusters/t200/d-model
  * [22] R. Gabl, T. Davey, Y. Cao, Q. Li, B. Li, K. L. Walker, et al., ”Experimental force data of a restrained ROV under waves and current”, _Data_ , vol. 5, no. 3, pp. 57, Jun. 2020.
  * [23] R. Gabl, T. Davey, Y. Cao, Q. Li, B. Li, K. L. Walker, et al., ”Hydrodynamic loads on a restrained ROV under waves and current”, _Ocean Engineering_ , vol. 234, 15 Aug. 2021.
  * [24] M. von Benzon, F. Sørensen, J. Liniger, S. Pedersen, S. Klemmensen and K. Schmidt, ”Integral Sliding Mode Control for a Marine Growth Removing ROV with Water Jet Disturbance,” _Proc. of 2021 European Control Conference (ECC)_ , Delft, Netherlands, pp. 2265-2270, 29 June - 02 July 2021.
  * [25] Blue Robotics, “BlueROV2: The world’s most affordable high-performance ROV”, _BlueROV2 datasheet_ , Jun. 2016 [Revised Jun. 2019].
  * [26] C. J. Wu, ”6-DoF Modelling and Control of a Remotely Operated Vehicle”, M.Eng. Thesis, College of Science and Engineering, Flinders University, 2018.
  * [27] _Wave Spectra Dataset: Feb-May 2021, Moray Firth WaveNet Site_ , Centre For Environment Fisheries and Aquaculture Science (Cefas) Online Repository, downloaded July 2022. [Online]. Available: https://www.cefas.co.uk/data-and-publications/wavenet/.
  * [28] Bailey H, Brookes K and Thompson P (2014) Assessing Environmental Impacts of Offshore Wind Farms: Lessons Learned and Recommendations for the Future. _Aquatic Biosystems_ , 10(8).
  * [29] Oh KY, Nam W, Ryu M S, Kim JY and Epureanu B I (2018) A review of foundations of offshore wind energy converters: Current status and future perspectives. _Renewable and Sustainable Energy Reviews_ , 88: 16-36.