* Marty and Aichelin (2013) R. Marty and J. Aichelin, Phys. Rev. C 87, 034912 (2013), eprint 1210.3476.
  * Lemmon and Mondragon (2010) T. J. Lemmon and A. R. Mondragon (2010), eprint 1012.5438.
  * Acharya et al. (2019) S. Acharya et al. (ALICE), JHEP 10, 084 (2019), eprint 1905.07211.
  * Acharya et al. (2022) S. Acharya et al. (ALICE), Phys. Rev. D 105, L011103 (2022), eprint 2105.06335.
  * Acharya et al. (2021a) S. Acharya et al. (ALICE) (2021a), eprint 2108.02523.
  * Abelev et al. (2012) B. Abelev et al. (ALICE), Phys. Lett. B 718, 295 (2012), [Erratum: Phys.Lett.B 748, 472–473 (2015)], eprint 1203.3641.
  * Acharya et al. (2021b) S. Acharya et al. (ALICE) (2021b), eprint 2109.15240.
  * Aaij et al. (2011) R. Aaij et al. (LHCb), Eur. Phys. J. C 71, 1645 (2011), eprint 1103.0423.
  * Aaij et al. (2021) R. Aaij et al. (LHCb), JHEP 11, 181 (2021), eprint 2109.00220.
  * Acharya et al. (2020b) S. Acharya et al. (ALICE), Phys. Lett. B 805, 135434 (2020b), eprint 1910.14404.
  * Nahrgang et al. (2014) M. Nahrgang, J. Aichelin, P. B. Gossiaux, and K. Werner, Phys. Rev. C 90, 024907 (2014), eprint 1305.3823.
  * Acharya et al. (2017) S. Acharya et al. (ALICE), Phys. Rev. Lett. 119, 242301 (2017), eprint 1709.05260.
  * Durand and O’Raifeartaigh (1976) B. Durand and L. O’Raifeartaigh, PHYSICAL REVIEW D13 (1976).


# Towards Accurate Vehicle Behaviour Classification With Multi-Relational
Graph Convolutional Networks

Sravan Mylavarapu1, Mahtab Sandhu2, Priyesh Vijayan3,
K Madhava Krishna2 , Balaraman Ravindran4, Anoop Namboodiri1 1Authors are with
the Center for Visual Information Technology, KCIS, IIIT Hyderabad2Authors are
with the Robotics Research Center, KCIS, IIIT Hyderabad3 School of Computer
Science, McGill University and Mila. Work done when the author was at Robert
Bosch Center for Data Science and AI.4 Dept. of CSE and Robert Bosch Center
for Data Science and AI, IIT Madras

###### Abstract

Understanding on-road vehicle behaviour from a temporal sequence of sensor
data is gaining in popularity. In this paper, we propose a pipeline for
understanding vehicle behaviour from a monocular image sequence or video. A
monocular sequence along with scene semantics, optical flow and object labels
are used to get spatial information about the object (vehicle) of interest and
other objects (semantically contiguous set of locations) in the scene. This
spatial information is encoded by a Multi-Relational Graph Convolutional
Network (MR-GCN), and a temporal sequence of such encodings is fed to a
recurrent network to label vehicle behaviours. The proposed framework can
classify a variety of vehicle behaviours to high fidelity on datasets that are
diverse and include European, Chinese and Indian on-road scenes. The framework
also provides for seamless transfer of models across datasets without
entailing re-annotation, retraining and even fine-tuning. We show comparative
performance gain over baseline Spatio-temporal classifiers and detail a
variety of ablations to showcase the efficacy of the framework.

##  I INTRODUCTION

Dynamic scene understanding in terms of vehicle behaviours is laden with many
applications from autonomous driving where decision making is an important
aspect to traffic violation. For example, a vehicle cutting into our lane
would require us to decide in moving opposite to it (refer to video for
applications). In this effort, we categorize on-road vehicle behaviours into
one of the following categories: parked, moving away, moving towards us, lane
change from left, lane change from right, overtaking. To realize this, we
decompose a traffic scene video into its spatial and temporal parts. The
spatial part makes use of per-frame semantics, object labels to locate objects
in the 3D world as well as model spatial relations between objects. The
temporal part makes use of flow to track the progress of per frame inter-
object relations over time.

We use Multi-Relational Graph Convolution Networks (MR-GCN)[1] that are now
gaining traction to learn per frame embeddings of an object’s (object of
interest) relation with its surrounding objects in the scene. A temporal
sequence of such embeddings is used by a recurrent network (LSTM)[2] and
attention model to classify the behaviour of the object of interest in the
scene. Apart from its high fidelity and seamless model transfer capabilities,
the attractiveness of the framework stems from its ability to work directly on
monocular stream data bypassing the need for explicit depth inputs such as
from LIDAR or stereo modalities. The end to end nature of the framework limits
explicit dependencies to the entailment of a well-calibrated monocular camera.
Primarily the proposed framework uses inter-object spatial relations and their
evolution over time to classify vehicle behaviours. Intuitively a car’s
temporal evolution of its spatial relations with its adjacent objects (also
called landmarks in this paper) help us, humans, to identify its behaviour. As
seen in Fig. 1 the change of car’s spatial relation with lane marking tells us
that the car is moving away from us. In a similar vein, an object’s relational
change with respect to lane markings also decides its lane change behaviour
while its changing relations with another moving car in the scene governs its
overtaking behaviour. The network’s ability to capture such intuition is a
cornerstone of this effort.

Figure 1: The figure depicts how objects (vehicle) behaviour can be modelled
as a Spatio-temporal scene graph. The evolution pattern of this scene graph
helps us to classify the behaviour. Here the car (red) moves ahead of a lane
marking (yellow) . This spatial relation helps us identify that the car is
moving ahead

.

Figure 2: This illustrates the overall pipeline of our framework. The input
(a) of the system consists of consecutive frames of a monocular video from
cameras mounted on the cars. We use traditional object tracking pipelines, as
shown in (b) to extract features. These extracted features/tracklets are then
transferred into the bird’s-eye view as in (d). Using these pre-processed
features in bird’s-eye view, we generate spatial scene graphs as illustrated
with examples in (e). MR-GCNs are used to encode spatial-relations between
objects from this scene graph. Subsequently, this spatial information from MR-
GCNs over different frames are passed to an LSTM to model the temporal
evolution of such spatial-relations and predict the different vehicle
behaviours.

Specifically the paper contributes as follows

  * •

We showcase Multi-relation GCN (MR-GCN) along with a recurrent network as an
effective framework for modelling Spatio-temporal relations in on-road scenes.
Though Graph Convolutional Networks (GCNs) have gained traction recently in
the context of modelling relational aspects of a scene, to the best of our
knowledge this is perhaps the first such reporting of using an MR-GCN in an
end-to-end framework for on-road temporal-scene modelling, ref III-B1.

  * •

We reproduce the performance benefits of our framework over a variety of
datasets, namely KITTI, Apollo-scape and Indian datasets, to showcase its
versatility. Also, our model can effectively transfer its learning seamlessly
across different datasets preventing any entailment for annotation, retraining
and fine-tuning (as shown in tables III, II)

  * •

Further results on diverse datasets and various on-road participants (not
restricted to cars) such as in Indian roads is another feature of this effort.
We also provide ablation studies that show that the model’s performance
improves with an increase in the number of objects (landmarks) in the scene.
Albeit, the model is yet fairly robust and superior to the baseline when only
a minimal fraction of the landmarks are available.

##  II RELATED WORK

Vehicle Behaviour Understanding : The closest works to our task are [3, 4, 5].
In [3, 5] vehicle’s behaviour is classified using a probabilistic model. The
classification is with respect to the ego-vehicle network, for example, a car
is only to be classified as overtaking if it overtakes the ego car. Whereas in
our proposed method, we classify vehicle behaviours with respect to a global
frame, i.e. we are able to identify and classify overtaking and other
behavioural interactions among any of the objects present in the visible
scene. Methods [6, 7] predict trajectories based on the object’s previous
interactions and trajectory. These methods are concerned with predicting
future motion rather than classifying maneuvers. [8] follows a hard rule-based
approach for the task of maneuver prediction in contrast to our learnable
approach, which is more generic. In [4] traffic scenes are summarized into
textual output. This is in contrast to the proposed framework, which tries to
classify objects of interest in terms of a specific behaviour exhibited over
the last few evolving time instants. To our best knowledge the works [9],
[10], [11] are the closest to our proposed method, but uses various external
sensors such as GPS, face camera and Vehicle dynamics and the classification
is for the ego-vehicle alone. In the Computer Vision community, event
recognition for stationary surveillance video is well researched over standard
datasets [12], [13]. Closer home in the robotics community, motion detection
methods such as [14, 15, 16, 17] some that have included semantics along with
motion [18, 19, 20] bear some resemblance. However, the proposed method goes
beyond the traditional task of motion segmentation. While understanding the
presence of motion in the scene, it presents a scalable architecture that
recognizes much more complex and diverse set of behaviours by just using
monocular videos.

Graph Neural Networks : Scene understanding has been traditionally formulated
as a graph problem [21, 22]. Graphs have been extensively used for scene
understanding, and most of them were solved using belief propagation [23],
[24]. While the success of deep learning approaches into the structured data
prediction problems of scene understanding like semantic segmentation and
motion segmentation have been successful, there has been relatively sparse
literature on scene understanding in the non-euclidean domains like video and
behaviour understanding. Recently there have been multiple successful attempts
at modelling learning algorithms for graph data and specifically, [25, 26, 27,
28, 29] have proposed algorithms to generalize deep learning to unstructured
data. With the success of these methods on the graph classification task,
multiple recent works have extended the methods to address multiple 3D
problems like Shape segmentation [30], 3D correspondence [31], CNN on surfaces
[32]. We model the relations in a dynamic scene using graph networks to
predict the behaviour of the vehicles. The end to end framework wherein the
vehicle maps a temporal sequence of spatial relations to scalable number of
behaviours is a distinguishing novelty of this effort.

##  III PROPOSED METHOD

We propose a pipeline to predict maneuver behaviours of different vehicles in
a video by modeling the evolution of their Spatio-temporal relations with
other objects. The behaviour prediction pipeline is a two-stage process; the
first stage involves modeling a video as a series of spatial scene graphs, and
the second stage involves modeling this Spatio-temporal information to predict
maneuver behaviors of vehicles.

###  III-A Scene Graph Construction

The scene graph construction is itself a multi-step process where we first
identify and track different objects from $T$ consecutive video frames. Next,
we orient these objects in a bird’s eye view to facilitate the identification
of relative positions of various identified objects. Then, for every frame, we
encode the relative positional information of different objects with a spatial
graph. The steps are elaborated below:

####  III-A1 Object Tracking

For a complete dynamic scene understanding, we use three main feature
extraction pipelines as input to identify the different objects and determine
the spatial relationships in a video. The major components required for
dynamic feature extraction are instance segmentation, semantic segmentation
and per-pixel optical flow in a video frame. We follow the method of [33] to
compute the instance segmentation of each object. Specifically, the instance
segmentation helps detect moving objects and for long term tracking of the
object. We follow the method of [34] to compute the semantic understating of
the scene specifically for static objects like poles and lane markings, which
act as landmarks. The optical flow is computed between two consecutive frames
in the video using [35]. Each of the objects in the scene is tracked by
combining optical flow information with features from the instance and
semantic segmentation.

####  III-A2 Monocular to Bird’s-eye view

Once we obtain stable tracks for objects in the image space, we move those
tracks from image space to a bird’s-eye view to better capture the relative
positioning of different objects. For the conversion from image to bird’s-eye
view, we employ Eqn: 1 as described in [36] for object localization in
bird’s-eye view. All the objects are assigned a reference point, for lanes it
is the center point of lane marking, for vehicles and other landmarks it is
the point adjacent to the road that is assigned as a reference point as seen
in Fig. 2 (d). Let $b=(x,y,1)$ be the reference point in homogeneous
coordinates of the image. We transfer these homogeneous co-ordinates $b$ to a
Bird’s-eye view as follows.

$\small B=(B_{x},B_{y},B_{z})=\frac{-hK^{-1}b}{\eta^{T}K^{-1}b}$ (1)

Here $K$ is the camera intrinsic calibration matrix, $\eta$ is the surface
normal, and $h$ is the camera height from the ground plane and $B$ is the 3D
co-ordinate (in the bird’s-eye view) of the point $b$ in the image.

####  III-A3 Spatial Scene Graph

The spatial information in a video frame at time $t$ is captured as a scene
graph $S^{t}$ based on the relative position of different objects in the
bird’s-eye view. The spatial relation, $S^{t}_{i,j}$ between a subject, $i$
and an object, $j$ at time, $t$ is the quadrant in which the object lies with
respect to the subject, i.e $S^{t}_{i,j}\in\mathcal{R}_{s}$, where
$\mathcal{R}_{s}=\\{$top left, top right, bottom left and bottom right$\\}$.
An example of a spatial scene graph is illustrated in Fig. 2 (e).

###  III-B Vehicle Maneuver Behaviour Prediction

We obtain Spatio-temporal representations of vehicles by modeling the temporal
dynamics of their spatial relations with other vehicles and landmarks in the
scenes over time. The Spatio-temporal representation is obtained by stacking
three different neural network layers, each with its own purpose. First, we
use a Multi-Relational Graph Convolutional Network (MR-GCN) layer [1] to
encode the spatial information of each object (node in graph) with respect to
its surroundings at time $t$ into an embedding $E^{t}$. Then, for all nodes,
it’s spatial embedding from all timesteps, $\\{E^{1},E^{2},...,E^{T}\\}$ is
fed into a Long Short-term Memory layer (LSTM) [2] to encode the temporal
evolution of its spatial information. The Spatio-temporal embeddings,
$\\{C^{1},C^{2},...,C^{T}\\}$ obtained at each timestep from the LSTM is then
passed through a (multi-head) self-attention layer to select relevant Spatio-
temporal information necessary to make behavior predictions. The details of
the different stages in the pipeline follow in order below. The overall
pipeline is best illustrated in Fig. 2.

####  III-B1 MR-GCN

The multi-relational interactions ($S^{t}_{i,j}\in\mathcal{R}_{s}$) between
objects are encoded using a Multi-Relational Graph Convolutional Network, MR-
GCN [1], originally proposed for knowledge graphs.

For ease of convenience, we reformulate the representation of each spatial
scene graph as $G$, where $G$ is defined as a set of $|\mathcal{R}_{s}|$
Adjacency matrices, $G=\\{A_{1},..,A_{|\mathcal{R}_{s}|}\\}$ corresponding to
the different relations. Herein, $A_{r}\in\mathbb{R}^{n\times n}$ with
$A_{r}[i,j]=1$ if $S_{i,j}=r$ otherwise $A_{r}[i,j]=0,\forall
r\in\mathcal{R}_{s}$; $n$ denotes number of nodes in the graph.

The output of the $k^{th}$ MR-GCN layer, $h^{k}_{G}$ is defined below in Eqn:
2. The MR-GCN layer convolves over neighbors from different relations
separately and aggregates them by a simple summation followed by the addition
of the node information.

$\displaystyle\small
h^{k}_{G}=ReLU\big{(}\sum_{r\in\mathcal{R}}\hat{A}_{r}h^{k-1}_{G}W_{r}^{k}+h^{k-1}_{G}W_{s}^{k}\big{)}$
(2)

where $\hat{A}_{r}$ is the degree normalized adjacency matrix of relation,
$r$, i.e $\hat{A}_{r}=D_{r}^{-1}A_{r}$ where $D_{r}$ is the degree matrix;
$W_{r}^{k}$ is the weights associated with the $r^{th}$ relation of the
$k^{th}$ MR-GCN layer and $W_{s}^{k}$ is the weight associated with computing
the node information (self-loop) at layer $k$; For each object, the input to
the first MR-GCN layer, are learned object embeddings
$\mathcal{E}_{o}\in\mathbb{R}^{|O|*d}$ corresponding to their object type;
$O=\\{vehicle,landmark\\}$ in our case. Thus, the input layer
$h_{G}^{0}\in\mathbb{R}^{n*d}$. If $d$ is the dimensions of a $k^{th}$ MR-GCN
layer, then $W_{s}^{k},\,W_{r}^{k}\in\mathbb{R}^{d*d}$ and
$h_{G}^{k+1}\in\mathbb{R}^{n*d}$.

The MR-GCN layer defined in Eqn: 2 provides a relational representation for
each node in terms of its immediate neighbors. Multi-hop representations for
nodes can be obtained by stacking multiple MR-GCN layers, i.e., stacking $K$
layers of MR-GCN provides a $K$-hop relational representation, $h^{K}_{G}$. In
essence, MR-GCN processes the spatial information from each timestep $t$
($S^{t}$) and outputs a $K$-hop spatial embedding, $E^{t}$ with
$E^{t}=h^{K}_{G}$.

####  III-B2 LSTM

To capture the temporal dynamics of spatial relations between objects, we use
the popular LSTM to process spatial embeddings of all the objects in the video
over time to obtain Spatio-temporal embeddings, $C$. At each time-step, $t$,
the LSTM takes in the current spatial embedding of the objects, $E^{t}$ and
outputs a Spatio-temporal embedding of the objects, $C^{t}$ conditioned on the
states describing the Spatio-temporal evolution information obtained until the
previous time-step as defined in Eqn: 3.

$\displaystyle\small C^{t}$ $\displaystyle=LSTM(E^{t},C^{t-1})$ (3)

####  III-B3 Attention

We further improve the temporal encoding of each node’s spatial evolution by
learning to focus on those temporal aspects that matter to understand vehicle
maneuvers. Specifically, here we use Multi-Head Self-Attention described in
[37] on the LSTM embeddings. The Self-Attention layer originally proposed in
the NLP literature, obtains a contextual representation of a data at a time,
$t$ in terms of data from other timesteps. The current data is referred to as
the query, $Q$, and the data at other time-steps is called value, $V$, and is
associated with a key, $K$. Attention function outputs a weighted sum of value
based on the context as defined in Eqn: 4.

$\displaystyle\small
Attention_{m}(\,Q,\,K,\,V\,)=softmax(\dfrac{Q\,K^{T}}{\sqrt{d_{k}}})V$ (4)

Herein our task, at every time-step, the query is based on the current LSTM
embedding, and the key and value are representations based on the LSTM
embeddings from other time-steps. For a particular time-step $t$, the
attention function looks out for temporal correlations between the Spatio-
temporal LSTM embeddings at the current time-step and the rest. High scores
are given to those time-steps with relevant information. Herein, we use the
multi-head attention, which aggregates (CONCAT) information from multiple
attention functions, each of which may potentially capture different
contextualized temporal correlations relevant for the end task. A multi-head
attention with $M$ heads is defined below:

$\displaystyle
Z^{t}=\mathop{\parallel}\limits_{m=1}^{m=M}Attention_{m}(C^{t}\,W_{m}^{Q},C^{t}\,W_{m}^{K},C^{t}\,W_{m}^{V})$
(5)

where $\parallel$ represents concatenation, $W_{m}^{Q}\in\mathbb{R}^{d\times
d_{k}}$ , $W_{m}^{K}\in\mathbb{R}^{d\times d_{k}}$ and
$W_{m}^{V}\in\mathbb{R}^{d\times d_{v}}$ are parameter projection matrices for
$m^{th}$ attention head. $Z^{t}\in\mathbb{R}^{n\times Md_{v}}$ is the
concatenated attention vectors for $t^{th}$ time-step from $M$ heads .

Information from all the time-steps, $Z$, are average pooled along the time-
axis and label predictions are made with them. All the three components in the
behavior prediction stage, are trained end-end by minimizing a cross-entropy
loss defined over the true labels, $Y$ and the predicted labels, $\hat{Y}$.
The final prediction components are described below.

$\displaystyle\small\begin{split}&U=Pool_{AVG}(Z)\\\ &\hat{Y}=UW_{l}\\\
&min~{}~{}CrossEntropyLoss\,(Y,\hat{Y})\end{split}$ (6)

where, $W_{l}$ projects $U$ to number of classes.

Figure 3: Figure depicts results on the Indian dataset. We can see Vehicle
class agnosticism of our Pipeline as non standard Vehicles such as bus and
oil-tanker are correctly classified Figure 4: Figures depicts complex Vehicles
maneuvers such as overtaking which can be seen in (a), the car on the right
represented in magenta color and lane change which can be seen in image (b)
and (c) with yellow on the Apollo dataset

##  IV Experiments and Analysis

We evaluate our vehicle behaviour understating framework on multiple publicly
available datasets captured from moving cameras. Code along with hyper-
parameters and additional implementation details are available on the project
web-page [https://ma8sa.github.io/temporal-MR-GCN/]

###  IV-A Training

All the graphs are constructed for a length of 10 time-steps on all datasets.
We empirically found that using 2 layers for MR-GCN provide us with the best
results, with first and second layers being 128, 32 respectively. Other
variations can be found on our project page. The output dimension for
Attention and LSTM are equal to their corresponding input. Training is done on
single Nvidia GeForce GTX 1080 Ti GPU.

###  IV-B Datasets

Our Method takes in monocular videos to understand the vehicle behaviour. We
evaluate our framework on three datasets, two publicly available datasets
namely, KITTI [38] and Apollo-scapes [39] and our proprietary dataset, Indian
dataset which is captured in a more challenging cluttered environment. For all
the three datasets we collected human annotated labels for the vehicle
maneuvers.

####  IV-B1 Apollo-scape Dataset

Apollo-scape harbors a wide range of driving scenarios/maneuvers such as
overtaking, Lane change. We selected this dataset as our primary dataset
because it contains considerable instances of complex behaviours such as
overtaking, Lane change which are negligible in other datasets.

####  IV-B2 KITTI Dataset

This dataset is constrained to a city and regions around the city along with
few highways. We pick sequences 4, 5, 10 for our purpose.

####  IV-B3 Indian Dataset

To account and infer the behaviour of the objects with high variation in
object types, we capture a dataset in the Indian driving conditions. This
dataset contains Non-Standard vehicles which are not previously modeled by the
framework.

###  IV-C Qualitative Analysis

The following color coding is employed across all qualitative results. Blue
denotes vehicles that are parked, red is indicative of vehicles that are
moving away, and green the vehicles that are moving towards us. Yellow
indicates lane changing: left to right, orange indicates lane changing: right
to left behaviour with magenta denoting overtaking. Complex behaviours such as
overtaking and changing lanes can be observed in Fig. 4 which shows
qualitative results obtained on Apollo dataset. Fig. 4 (a) depicts a car
(right most car) in the midst of overtaking (seen with magenta color). Fig. 4
(b) and 4 (c) present us with cases of lane change. In Fig. 4 (b) we see a van
changing lane, along with a moving car, In Fig. 4 (c) we observe a bus
changing lane from left to right along with two cars parked on the right side.
Object class agnosticism can be observed in Fig. 3 (Indian dataset), where
behaviour of several non-standard vehicles are accurately predicted. In Fig. 3
(a) we can see a car on the left side of the road (predicted changing lane)
coming onto the road along with parked cars and one car coming towards us. In
Fig. 3 (b) truck and bus are parked on left side of the road while another bus
is accurately classified as moving forward and not overtaking. Note that
behaviour is classified as overtaking only if a moving vehicle overtakes
(moves ahead of) another moving Vehicle. Fig. 3 (c) provides a case with a
semi-truck and bus moving away from us, along with a truck parked. Fig. 5
shows qualitative results obtained on KITTI sequences 4, 5, 10.

###  IV-D Quantitative Results

In-depth analysis of our model’s performance is available in the form of
confusion matrix, ref table I. To show that our method is sensor invariant we
evaluate our model on KITTI and Indian datasets, results are tabulated in
table II. We provide ablation studies on model as well as land marks and
compare our model with few baselines.

Figure 5: Figure depicts Results obtained on KITTI Tracking dataset.

####  IV-D1 Ablation study on LandMarks

Landmarks are the stationary points (lane marking, poles, etc.) which along
with other vehicles are used by MR-GCN to encode a spatially-aware embedding
for a node. We carry out an ablation study to understand the effects of how
the number of landmarks available in a sequence affects the performance of our
Pipeline. We train three models, one with 50 percent landmarks available,
another with 75 percent available and the last one with all the landmarks
available. The results mentioned in above experiments are tabulated in table
III under their respective rows ours(50%), ours(75%), ours(full). Behaviours
which are directly dependent on landmarks (lane change, moving ahead and back)
show a difference of atleast 10% when comparing the full model and the model
trained with 75% landmarks and a difference of atleast 15% when compared with
a model trained with 50% landmarks. The difference is less significant for
overtaking as its influenced by other moving vehicles. We observed that with
less number of landmarks, the network is biased towards ’parked’ label, hence
the accuracy of it remains high, while others reduce. The observed results
make it evident that the proposed model is highly robust even with fewer
landmarks and also, that performance improvement is achievable with the
utilization of more information in the form of increased landmarks.

####  IV-D2 Ablation study on Model

To study the importance of each component of our system, we perform an
ablation study on the model. We observe the importance of encoding spatial
information along with temporal aspects by comparing variations of our network
with each other, as described in table IV. First, we use MR-GCN followed by a
simple LSTM (G + L) which provides below par results but performs
comparatively better than MR-GCN followed by a single-head attention (G + SA).
Adding a single head Attention to MR-GCN and LSTM (G + L + SA), helps to
weight the spatio-temporal encodings as described in III-B3 providing a huge
improvement in accuracies. The final model with multi head attention (G+L+MA)
attains improvement specially for $overtake$ showing it’s complex behaviour.
Overall, each component in the network provides a meaningful representation
based on the information encoded in it.

Baseline comparison :  We observe the importance of spatial information
encoded as a graph by comparing the performance of our method (G+L+MA) against
a simple LSTM (L) and a LSTM followed by Attention (L + MA) models, which
encode relational information with a 3D location based positional features;
ref table IV. To train these baseline models (L and L + MA), we give 3D
locations of objects obtained from the bird’s eye view (III-A2) as a direct
input. For each object in the video, we create a feature vector consisting of
its distance and angle with all other nodes for $T$ time-steps. The distance
is a simple Euclidean distance. To account for the feature of an object (lane-
markings/Vehicles), we create a 2 dimensional one-hot vector {1,0}
representing Vehicles and {0,1} representing static objects. To create a
feature vector for the $i^{th}$ object, we find distances and angles with all
other nodes for $T$ time-steps in the scene. In both the cases, the above
features are fed to the LSTM at each time-step. We pool the output from LSTM
for L and output from Attention layer for L + MA, along time dimension and
project to number of classes by using a dense layer.
Table IV shows clear improvement in terms of accuracies for models having a
MR-GCN to encode the spatial information when compared to the simple recurrent
models (L and L + MA). Our method outperforms the baselines by a good margin
for all classes.

The complete method and algorithm for creating feature vectors and training
the models are described in detail in the project web page.
[https://ma8sa.github.io/temporal-MR-GCN/].

GT Predicted | MVA | MTU | PRK | LCL | LCR | OVT | total
---|---|---|---|---|---|---|---
MVA | 691 | 1 | 13 | 16 | 4 | 83 | 808
MTU | 1 | 212 | 21 | 0 | 2 | 0 | 236
PRK | 8 | 39 | 1330 | 19 | 7 | 10 | 1413
LCL | 10 | 0 | 6 | 137 | 0 | 8 | 161
LCR | 6 | 0 | 5 | 3 | 112 | 3 | 129
OVT | 18 | 0 | 1 | 0 | 0 | 53 | 72
TABLE I: confusion matrix for model trained on Apollo dataset. Due to space constraints we have used abbreviations for label names, MVA : moving away from us, MTU: moving towards us, PRK: parked, LCL: lane change left to right, LCR: lane change right to left, OVT: overtaking and GT : Ground truth Trained on | Apollo | KITTI | Indian
---|---|---|---
Tested on | KITTI | Indian | KITTI | Indian
Moving away from us | 99 | 99 | 85 | 85
Moving towards us | 98 | 93 | 86 | 74
Parked | 99 | 99 | 89 | 84
TABLE II: This showcases the transfer learning capabilities of our Method. Values provided are accuracies Tested on | Apollo-Scape
---|---
Method |  ours (50%) |  ours (75%) |  ours (full)
Moving away | 72 | 76 | 85.3
Moving towards us | 67 | 75 | 89.5
Parked | 97 | 97 | 94.8
LC left - right | 66 | 74 | 84.1
LC right - left | 71 | 72 | 86.4
Overtaking | 68 | 65 | 72.3
TABLE III: Comparison between models trained with 50% of landmarks and 75% and 100% of landmarks. Values provided are accuracies. Classes | MVA | MTU | PRK | LCL | LCR | OVT
---|---|---|---|---|---|---
Architecture |  |  |  |  |  |
G + L + MA | 85 | 89 | 94 | 84 | 86 | 72
G + L + SA | 84 | 75 | 95 | 51 | 65 | 51
G + L | 78 | 72 | 74 | 37 | 49 | 41
G + SA | 70 | 18 | 54 | 38 | 14 | 40
L | 37 | 35 | 34 | 6 | 24 | 15
L + MA | 56 | 41 | 46 | 8 | 27 | 13
TABLE IV: Baselines and ablation study on the model trained on Apollo.
Abbreviations for model architectures (rows) are, G : MR-GCN, L : LSTM, MA :
Multi-Head Attention, SA : Single-Head Attention. Nomenclature in the
architecture names reflect the components in it. column abbreviations are
similar to table I. All the values provided are accuracies.

####  IV-D3 Transfer Learning

As embeddings obtained from MR-GCN are agnostic to the visual scene, they are
dependent on the scene graph obtained from the visual data (images). To show
this, we trained the model only on Apollo dataset and tested it on Indian and
KITTI. While testing, we removed the lane change and overtaking behaviour as
it is not present in KITTI and Indian dataset. In table II, we observe that
results are better with the model trained on Apollo dataset as compared to
models trained on their respective datasets. This could be attributed to the
presence of more complex behaviours (overtaking, lane change) in Apollo and
shear difference in size of datasets (Apollo has 4K frames as compared to 651
and 722 frames in Indian and KITTI datasets).

##  V Conclusion

This paper is the first to pose the problem of classifying a Vehicle-of-
interest in a scene in terms of its behaviours such as ”Parked”, ”Lane
Changing”, ”Overtaking” and the like by leveraging static landmarks and
changing relations. The main highlight of the model’s architecture is its
ability to learn behaviours in an end-to-end framework as it successfully maps
a temporal sequence of evolving relations to labels with repeatable accuracy.
The framework helps in detecting traffic violations such as overtake on narrow
bridges and in deciding to maintain safe distance from a driver (aggressive)
changing his state frequently. The ability to transfer across datasets with
high fidelity and robustness in the presence of a reduced number of objects
and an improved performance over base-line methods summarizes the rest of the
paper.

##  VI Acknowledgement

The work described in this paper is supported by MathWorks. The opinions and
views expressed in this publication are from the authors, and not necessarily
that of the funding bodies.

## References

  * [1] M. Schlichtkrull, T. N. Kipf, P. Bloem, R. Van Den Berg, I. Titov, and M. Welling, “Modeling relational data with graph convolutional networks,” in _European Semantic Web Conference_. Springer, 2018, pp. 593–607.
  * [2] F. A. Gers, J. Schmidhuber, and F. Cummins, “Learning to forget: Continual prediction with lstm,” 1999.
  * [3] A. Kanazawa, S. Tulsiani, A. A. Efros, and J. Malik, “Learning category-specific mesh reconstruction from image collections,” _CoRR_ , vol. abs/1803.07549, 2018.
  * [4] C.-Y. Chen, W. Choi, and M. Chandraker, “Atomic scenes for scalable traffic scene recognition in monocular videos,” in _2016 IEEE Winter Conference on Applications of Computer Vision (WACV)_. IEEE, 2016, pp. 1–9.
  * [5] S. Sivaraman and M. M. Trivedi, “Looking at vehicles on the road: A survey of vision-based vehicle detection, tracking, and behavior analysis,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 14, no. 4, pp. 1773–1795, 2013.
  * [6] J. Schulz, C. Hubmann, J. Löchner, and D. Burschka, “Interaction-aware probabilistic behavior prediction in urban environments,” in _2018 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 2018, pp. 3999–4006.
  * [7] R. Sabzevari and D. Scaramuzza, “Multi-body motion estimation from monocular vehicle-mounted cameras,” _IEEE Transactions on Robotics_ , vol. 32, no. 3, pp. 638–651, 2016.
  * [8] X. Geng, H. Liang, B. Yu, P. Zhao, L. He, and R. Huang, “A scenario-adaptive driving behavior prediction approach to urban autonomous driving,” _Applied Sciences_ , vol. 7, no. 4, p. 426, 2017.
  * [9] A. Jain, H. S. Koppula, B. Raghavan, S. Soh, and A. Saxena, “Car that knows before you do: Anticipating maneuvers via learning temporal driving models,” in _Proceedings of the IEEE International Conference on Computer Vision_ , 2015, pp. 3182–3190.
  * [10] A. Jain, A. R. Zamir, S. Savarese, and A. Saxena, “Structural-rnn: Deep learning on spatio-temporal graphs,” in _Proceedings of the ieee conference on computer vision and pattern recognition_ , 2016, pp. 5308–5317.
  * [11] A. Narayanan, I. Dwivedi, and B. Dariush, “Dynamic traffic scene classification with space-time coherence,” 2019.
  * [12] X. Wang and Q. Ji, “Video event recognition with deep hierarchical context model,” in _CVPR_ , 2015.
  * [13] A. H. e. a. Sangmin Oh, “A large-scale benchmark dataset for event recognition in surveillance video,” in _CVPR_ , 2011.
  * [14] R. K. Namdev, A. Kundu, K. M. Krishna, and C. Jawahar, “Motion segmentation of multiple objects from a freely moving monocular camera,” in _Robotics and Automation (ICRA), 2012 IEEE International Conference on_. IEEE, 2012, pp. 4092–4099.
  * [15] R. Vidal and S. Sastry, “Optimal segmentation of dynamic scenes from two perspective views,” in _CVPR_ , vol. 2, 2003.
  * [16] A. Kundu, K. Krishna, and J. Sivaswamy, “Moving object detection by multi-view geometric techniques from a single camera mounted robot,” in _IROS_ , 2009\.
  * [17] N. Fanani, M. Ochs, A. Sturck, and R. Mester, “Cnn-based multi-frame imo detection from a monocular camera,” in _2018 IEEE Intelligent Vehicles Symposium (IV)_. IEEE, 2018, pp. 957–964.
  * [18] N. D. Reddy, P. Singhal, V. Chari, and K. M. Krishna, “Dynamic body vslam with semantic constraints,” in _Intelligent Robots and Systems (IROS), 2015 IEEE/RSJ International Conference on_. IEEE, 2015, pp. 1897–1904.
  * [19] J. Vertens, A. Valada, and W. Burgard, “Smsnet: Semantic motion segmentation using deep convolutional neural networks,” in _IROS_ , 2017.
  * [20] T. Chen and S. Lu, “Object-level motion detection from moving cameras,” _IEEE Transactions on Circuits and Systems for Video Technology_ , 2016.
  * [21] M. A. Fischler and R. A. Elschlager, “The representation and matching of pictorial structures,” _IEEE Transactions on computers_ , vol. 100, no. 1, pp. 67–92, 1973.
  * [22] D. Marr and H. K. Nishihara, “Representation and recognition of the spatial organization of three-dimensional shapes,” _Proc. R. Soc. Lond. B_ , vol. 200, no. 1140, pp. 269–294, 1978.
  * [23] P. F. Felzenszwalb and D. P. Huttenlocher, “Pictorial structures for object recognition,” _International journal of computer vision_ , vol. 61, no. 1, pp. 55–79, 2005.
  * [24] L. Sigal, M. Isard, H. Haussecker, and M. J. Black, “Loose-limbed people: Estimating 3d human pose and motion using non-parametric belief propagation,” _International journal of computer vision_ , vol. 98, no. 1, pp. 15–48, 2012.
  * [25] D. K. Duvenaud, D. Maclaurin, J. Iparraguirre, R. Bombarell, T. Hirzel, A. Aspuru-Guzik, and R. P. Adams, “Convolutional networks on graphs for learning molecular fingerprints,” in _Advances in Neural Information Processing Systems 28_ , C. Cortes, N. D. Lawrence, D. D. Lee, M. Sugiyama, and R. Garnett, Eds. Curran Associates, Inc., 2015, pp. 2224–2232.
  * [26] T. N. Kipf and M. Welling, “Semi-supervised classification with graph convolutional networks,” _arXiv preprint arXiv:1609.02907_ , 2016.
  * [27] J. Bruna, W. Zaremba, A. Szlam, and Y. LeCun, “Spectral networks and locally connected networks on graphs,” _CoRR_ , vol. abs/1312.6203, 2013.
  * [28] M. Henaff, J. Bruna, and Y. LeCun, “Deep convolutional networks on graph-structured data,” _CoRR_ , vol. abs/1506.05163, 2015.
  * [29] M. Defferrard, X. Bresson, and P. Vandergheynst, “Convolutional neural networks on graphs with fast localized spectral filtering,” _CoRR_ , vol. abs/1606.09375, 2016.
  * [30] L. Yi, H. Su, X. Guo, and L. J. Guibas, “Syncspeccnn: Synchronized spectral cnn for 3d shape segmentation.” in _CVPR_ , 2017, pp. 6584–6592.
  * [31] O. Litany, T. Remez, E. Rodola, A. Bronstein, and M. Bronstein, “Deep functional maps: Structured prediction for dense shape correspondence,” in _2017 IEEE International Conference on Computer Vision (ICCV)_. IEEE, 2017, pp. 5660–5668.
  * [32] H. Maron, M. Galun, N. Aigerman, M. Trope, N. Dym, E. Yumer, V. G. Kim, and Y. Lipman, “Convolutional neural networks on surfaces via seamless toric covers,” _ACM Transactions on Graphics (TOG)_ , vol. 36, no. 4, p. 71, 2017\.
  * [33] K. He, G. Gkioxari, P. Dollár, and R. Girshick, “Mask r-cnn,” in _ICCV_ , 2017.
  * [34] S. Rota Bulò, L. Porzi, and P. Kontschieder, “In-place activated batchnorm for memory-optimized training of dnns,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , 2018, pp. 5639–5647.
  * [35] D. Pathak, R. Girshick, P. Dollár, T. Darrell, and B. Hariharan, “Learning features by watching objects move,” in _CVPR_ , 2017.
  * [36] S. Song and M. Chandraker, “Robust scale estimation in real-time monocular SFM for autonomous driving,” in _2014 IEEE Conference on Computer Vision and Pattern Recognition, CVPR 2014, Columbus, OH, USA, June 23-28, 2014_. IEEE Computer Society, 2014, pp. 1566–1573. [Online]. Available: https://doi.org/10.1109/CVPR.2014.203
  * [37] A. Vaswani, N. Shazeer, N. Parmar, J. Uszkoreit, L. Jones, A. N. Gomez, Ł. Kaiser, and I. Polosukhin, “Attention is all you need,” in _Advances in neural information processing systems_ , 2017, pp. 5998–6008.
  * [38] A. Geiger, P. Lenz, C. Stiller, and R. Urtasun, “Vision meets robotics: The kitti dataset,” _International Journal of Robotics Research (IJRR)_ , 2013\.
  * [39] X. Huang, X. Cheng, Q. Geng, B. Cao, D. Zhou, P. Wang, Y. Lin, and R. Yang, “The apolloscape dataset for autonomous driving,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops_ , 2018, pp. 954–960.


# A non-intrusive neural-network based gradient algorithm for parameter
estimation in non-stationary elasticity

Stefan Frei  Jan Reichle and Stefan Volkwein Department of Mathematics and
Statistics, University of Konstanz, 78457 Konstanz, Germany, {stefan.frei,
jan.reichle<EMAIL_ADDRESS>
###### Abstract

We present a non-intrusive gradient algorithm for parameter estimation
problems in non-stationary elasticity. To avoid multiple (and potentially
expensive) solutions of the underlying partial differential equation (PDE), we
approximate the PDE solver by a neural network within the gradient algorithm.
The network is trained offline for a given set of parameters. The algorithm is
applied to an unsteady linear-elastic contact problem; its convergence and
approximation properties are investigated numerically.

##  1 Introduction

We are interested in the solution of parameter estimation problems subject to
computationally expensive partial differential equations (PDE). Typically, the
available algorithms to solve such problems require multiple solutions of the
underlying PDE to compute state and sensitivities. To fix ideas, we study a
model from elasticity describing the deformation of a metal plate under
applied forces, where the material parameters $p\in{\cal P}^{\text{ad}}$ are
unknown and have to be estimated. In industrial practice commercial software
is frequently used to predict, for instance, the resulting deformations and
stresses. In many cases, the software is not accessible and does not allow for
a computation of sensitivities with respect to the parameters.

In this work we apply a neural network to approximate the solution of a non-
stationary linear-elastic contact problem. For an overview of data-based
approaches to solve inverse problems, we refer, e.g., to [3], but see also [2,
10]. The approximation quality of neural networks is studied in [7], for
instance, and the many references cited therein.

The paper is organized as follows: In Section 2 the parameter identification
problem is introduced as a PDE-constrained minimization problem. The finite
element (FE) discretization and the neural network are briefly explained in
Sections 3 and 4, respectively. Section 5 is devoted to the gradient based
optimization method, where the derivatives are computed by backpropagation. In
Section 6 a numerical example is presented. Finally, a conclusion is made in
Section 7.

##  2 Problem formulation

In the present work we consider a non-stationary linear-elastic contact
problem, see, e.g., [5], in a time interval $[0,T]$ and in a two-dimensional
bounded spatial domain $\Omega$. However, our approach can also be applied to
other evolution equations. We assume that a plate lies in between a fixed
rigid geometry in the form of a ridge, see Figure 1.

Figure 1: Spatial domain $\Omega$ and its finite element triangulation
$\mathscr{T}_{h}$.

The plate is pulled from the right boundary (${\Gamma_{\mathrm{r}}}$), while
it is fixed on the left one (${\Gamma_{\mathrm{l}}}$). On top and bottom
(${\Gamma_{\mathrm{tb}}}$) we apply frictionless contact conditions which
allow the plate to deform only towards its interior. Then, we have
$\partial\Omega={\Gamma_{\mathrm{tb}}}\cup{\Gamma_{\mathrm{lr}}}$, where
${\Gamma_{\mathrm{lr}}}={\Gamma_{\mathrm{l}}}\cup{\Gamma_{\mathrm{r}}}$.
Moreover, we set

$\displaystyle\bm{H}=L^{2}(\Omega;\mathbb{R}^{2})\quad\text{and}\quad{\bm{V}_{\mathrm{d}}}=\big{\\{}{\bm{\varphi}}\in
H^{1}(\Omega;\mathbb{R}^{2})\,\big{|}\,{\bm{\varphi}}=0\text{ on
}{\Gamma_{\mathrm{lr}}}\big{\\}}.$

For $T>0$ let $Q=(0,T)\times\Omega$,
${\Sigma_{\mathrm{tb}}}=(0,T)\times{\Gamma_{\mathrm{tb}}}$ and
${\Sigma_{\mathrm{lr}}}=(0,T)\times{\Gamma_{\mathrm{lr}}}$. The resulting
deformation vector ${\bm{u}}=(u_{1},u_{2})$ is governed by the initial-
boundary value problem

$\displaystyle\rho{\bm{u}}_{tt}-\nabla\cdot\sigma({\bm{u}})$
$\displaystyle=\bm{f}$ $\displaystyle\text{in }Q,$ (1a)
$\displaystyle{\bm{u}}\cdot{\bm{n}}\leq 0,\hskip
2.84526pt{\bm{n}}\cdot\big{(}\sigma({\bm{u}}){\bm{n}}\big{)}$
$\displaystyle\leq 0,\hskip
2.84526pt({\bm{u}}\cdot{\bm{n}})\big{(}({\bm{n}}\cdot(\sigma({\bm{u}}){\bm{n}})\big{)}=0$
$\displaystyle\text{on }{\Sigma_{\mathrm{tb}}},$ (1b)
$\displaystyle{\bm{\tau}}\cdot\big{(}\sigma({\bm{u}}){\bm{n}}\big{)}$
$\displaystyle=0$ $\displaystyle\text{on }{\Sigma_{\mathrm{tb}}},$ (1c)
$\displaystyle{\bm{u}}$ $\displaystyle={{\bm{u}}_{\mathrm{d}}}$
$\displaystyle\text{on }{\Sigma_{\mathrm{lr}}},$ (1d)
$\displaystyle{\bm{u}}(0,\cdot)={\bm{u}}_{t}(0,\cdot)$ $\displaystyle=0$
$\displaystyle\text{on }\Omega$ (1e)

with an inhomogeneity $\bm{f}\in L^{2}(0,T;\bm{H})$, the matrices
$\sigma({\bm{u}})=2\mu\epsilon(u)+\lambda\mathrm{tr}(\epsilon({\bm{u}}))I$ and
$\epsilon({\bm{u}})=\nicefrac{{1}}{{2}}\,(\nabla{\bm{u}}+(\nabla{\bm{u}})^{\top})$,
where $I$ stands for the $2\times 2$ identity matrix. The (positive) Lamé
parameters $\lambda$ and $\mu$ correspond to the elasticity modulus $E$ and to
the Poisson ratio $\nu$ (which will be the parameters we wish to determine by
means of experiments) by the relations

$\displaystyle\lambda=\frac{E\nu}{(1+\nu)(1-2\nu)},\qquad\mu=\frac{E}{2(1+\nu)}.$

The convex and compact set of admissible parameters is given as

$\displaystyle{\mathscr{P}_{\mathsf{ad}}}\colonequals\big{\\{}p=(E,\nu)\in{\mathbb{R}}^{2}\,\big{|}\,{E_{\mathsf{a}}}\leq
E\leq{E_{\mathsf{b}}},~{}{\rm
NewA}\leq\nu\leq{\nu_{\mathsf{b}}}\big{\\}}\equalscolon\big{[}{E_{\mathsf{a}}},{E_{\mathsf{b}}}\big{]}\times\big{[}{\rm
NewA},{\nu_{\mathsf{b}}}\big{]},$

where $0<{E_{\mathsf{a}}}\leq{E_{\mathsf{b}}}$ and $0<{\rm
NewA}\leq{\nu_{\mathsf{b}}}<\nicefrac{{1}}{{2}}$. The density is fixed to
$\rho=2700\,\nicefrac{{\mathrm{kg}}}{{\mathrm{m}^{3}}}$ and the (spatially
independent) Dirichlet conditions are given by

$\displaystyle{{\bm{u}}_{\mathrm{d}}}=0\text{ on
}(0,T]\times{\Gamma_{\mathrm{l}}},\quad{{\bm{u}}_{\mathrm{d}}}(t)=(18t^{2}-12t^{3},0)\text{
on }(0,T]\times{\Gamma_{\mathrm{r}}}.$

In (1b) and (1d) the vectors ${\bm{n}}$ and ${\bm{\tau}}$ denote the outer
unit normal vector and the corresponding tangential vector, respectively.

To deal with the contact conditions (1b) we use Nitsche’s method [1]. It can
be shown [1] that the three conditions (1b) are equivalent to the following
equality for arbitrary $\gamma>0$:

${\bm{n}}\cdot\big{(}\sigma({\bm{u}}){\bm{n}}\big{)}=-\frac{1}{\gamma}\left[{\bm{u}}\cdot{\bm{n}}-\gamma{\bm{n}}\cdot\big{(}\sigma({\bm{u}}){\bm{n}}\big{)}\right]_{+},$
(2)

where $[g]_{+}\colonequals\max\\{g,0\\}$ denotes the positive part of any
function $g$. Condition (2) can be easily incorporated in a variational
formulation; cf. [5], for instance. For that purpose we introduce the velocity
vector ${\bm{v}}={\bm{u}}_{t}$. Then the variational formulation of (1) is
given as

$\displaystyle\begin{aligned}
\frac{\mathrm{d}}{\mathrm{d}t}\,{\langle{\bm{v}},{\bm{\varphi}}\rangle}_{\bm{H}}+\int_{\Omega}\sigma({\bm{u}}):\nabla{\bm{\varphi}}\,\mathrm{d}{\bm{x}}\hskip
56.9055pt&\\\
+\frac{1}{\gamma}\,\Big{\langle}\left[{\bm{u}}\cdot{\bm{n}}-\gamma{\bm{n}}\cdot\big{(}\sigma({\bm{u}}){\bm{n}}\big{)}\right]_{+},{\bm{\varphi}}\cdot{\bm{n}}\Big{\rangle}_{L^{2}({\Gamma_{\mathrm{tb}}})}&=0\quad\forall{\bm{\varphi}}\in{\bm{V}_{\mathrm{d}}}\text{
in }(0,T],\\\
\frac{\mathrm{d}}{\mathrm{d}t}\,{\langle{\bm{u}},{\bm{\psi}}\rangle}_{\bm{H}}-{\langle{\bm{v}},{\bm{\psi}}\rangle}_{\bm{H}}&=0\quad\forall{\bm{\psi}}\in\bm{H}\text{
in }(0,T],\end{aligned}$ (3)

where the symbol “$\,:\,$” stands for the Frobenius product of matrices.

For known results on the well-posedness of dynamic contact problems, we refer
to [6, Chapter 4]. Throughout, we assume that the solution ${\bm{u}}$ to (3)
belongs to $C(\overline{Q})$.

Next we formulate the parameter identification as an optimization problem.
Suppose we are given five spatial points ${\bm{x}}_{i}=(x_{i},y_{i})\in\Omega$
for $i=1,\ldots,5$. We set $t_{1}=\nicefrac{{T}}{{2}}$, $t_{2}=T$ and
${y_{\mathsf{a}}^{i}}\colonequals\inf\\{y\in{\mathbb{R}}\,|\,(x_{1}^{i},y)\in\Omega\\}$,
${y_{\mathsf{b}}^{i}}\colonequals\sup\\{y\in{\mathbb{R}}\,|\,(x_{1}^{i},y)\in\Omega\\}$
for $i=1,\ldots,5$. Assume that we have the following available measurements:

  * 1)

point values

$\displaystyle
C^{1}({\bm{u}})\colonequals\big{(}{\bm{u}}({\bm{x}}_{1},t_{1}),{\bm{u}}({\bm{x}}_{2},t_{1}),\ldots,{\bm{u}}({\bm{x}}_{5},t_{1}),{\bm{u}}({\bm{x}}_{1},t_{2}),\ldots,{\bm{u}}({\bm{x}}_{5},t_{2}\big{)}\in{\mathbb{R}}^{2\times
10}$

of the deformation vector;

  * 2)

mean values

$\displaystyle C^{2}({\bm{u}})$
$\displaystyle\colonequals\textstyle\Big{(}\kappa_{1}\int_{y_{\mathsf{a}}^{1}}^{y_{\mathsf{b}}^{1}}{\bm{u}}(x_{1},y,t_{1})\,\mathrm{d}y,\ldots,\kappa_{5}\int_{y_{\mathsf{a}}^{5}}^{y_{\mathsf{b}}^{5}}{\bm{u}}(x_{5},y,t_{1})\,\mathrm{d}y,$
$\displaystyle\hskip
19.91692pt\textstyle\kappa_{1}\int_{y_{\mathsf{a}}^{1}}^{y_{\mathsf{b}}^{1}}{\bm{u}}(x_{1},y,t_{2})\,\mathrm{d}y,\ldots,\kappa_{5}\int_{y_{\mathsf{a}}^{5}}^{y_{\mathsf{b}}^{5}}{\bm{u}}(x_{5},y,t_{2})\,\mathrm{d}y\Big{)}\in{\mathbb{R}}^{2\times
10}$

with the factor
$\kappa_{i}=\nicefrac{{1}}{{({y_{\mathsf{b}}^{i}}-{y_{\mathsf{a}}^{i}})}}$;

  * 3)

mean values

$\displaystyle C^{3}({\bm{u}})$
$\displaystyle\colonequals\textstyle\Big{(}\kappa_{1}\int_{y_{\mathsf{a}}^{1}}^{y_{\mathsf{b}}^{1}}\sigma_{\text{vM}}(x_{1},y,t_{1})\,\mathrm{d}y,\ldots,\kappa_{5}\int_{y_{\mathsf{a}}^{5}}^{y_{\mathsf{b}}^{5}}\sigma_{\text{vM}}(x_{5},y,t_{1})\,\mathrm{d}y,$
$\displaystyle\hskip
19.91692pt\textstyle\kappa_{1}\int_{y_{\mathsf{a}}^{1}}^{y_{\mathsf{b}}^{1}}\sigma_{\text{vM}}(x_{1},y,t_{2})\,\mathrm{d}y,\ldots,\kappa_{5}\int_{y_{\mathsf{a}}^{5}}^{y_{\mathsf{b}}^{5}}\sigma_{\text{vM}}(x_{5},y,t_{2})\,\mathrm{d}y\Big{)}\in{\mathbb{R}}^{1\times
10}$

of the van-Mises stress

$\displaystyle\sigma_{\text{vM}}\colonequals\big{(}\nicefrac{{3}}{{2}}\,\sigma^{\text{dev}}({\bm{u}}):\sigma^{\text{dev}}({\bm{u}})\big{)}^{1/2},\text{
where
}\sigma^{\text{dev}}({\bm{u}})\colonequals\sigma({\bm{u}})-\nicefrac{{1}}{{3}}\,\mathrm{tr}(\sigma({\bm{u}}))I.$

These $m\colonequals 50$ measurements define the nonlinear observation
operator

$\displaystyle{\bm{u}}\mapsto{\mathcal{C}}({\bm{u}})={\big{(}C^{1}_{1,\cdot}({\bm{u}}),C^{1}_{2,\cdot}({\bm{u}}),C^{2}_{1,\cdot}({\bm{u}}),C^{2}_{2,\cdot}({\bm{u}}),C^{3}({\bm{u}})\big{)}^{\top}}\in{\mathbb{R}}^{m}.$

For given measurement vector ${{\bm{u}}^{\mathrm{obs}}}\in{\mathbb{R}}^{m}$ we
consider the minimization problem

$\displaystyle\begin{aligned} &\min
J({\bm{u}},(E,\nu))\colonequals\frac{1}{2}\,{\|{\mathcal{C}}({\bm{u}})-{{\bm{u}}^{\mathrm{obs}}}\|}^{2}_{{\mathbb{R}}^{m}}+\frac{\alpha}{2}\,\bigg{(}\Big{(}\frac{E}{E_{0}})^{2}+\big{(}\frac{\nu}{\nu_{0}}\Big{)}^{2}\bigg{)}\\\
&\hskip 1.42262pt\text{subject to }({\bm{u}},(E,\nu))\text{
satisfies~{}\eqref{Varform} and
}(E,\nu)\in{\mathscr{P}_{\mathsf{ad}}},\end{aligned}$ (4)

where the parameters $E,\nu$ are normalized by characteristic values
$E_{0},\nu_{0}$ of elasticity modulus and Poisson ratio, respectively.

ptLayer 0(Input Layer)ptLayer 1ptLayer 2ptLayer 3(Output Layer)

Figure 2: Example of a dense neural network with two hidden layers.

##  3 Finite element (FE) discretization

For discretization, we use the space of $P_{1}$ finite elements on a
triangular mesh $\mathscr{T}_{h}$:

$\displaystyle V^{h}\colonequals\big{\\{}v_{h}\in
C(\overline{\Omega})\,\big{|}\,v_{h}|_{T}\in P_{1}(T)\,\forall
T\in\mathscr{T}_{h}\big{\\}},~{}{\bm{V}^{h}}=V^{h}\times
V^{h},~{}{\bm{V}^{h}_{\mathrm{d}}}\colonequals{\bm{V}^{h}}\cap{\bm{V}_{\mathrm{d}}}.$

The FE mesh is illustrated in Figure 1, and the discrete formulation of (3)
reads: Find ${\bm{u}}^{h}\in{{\bm{u}}_{\mathrm{d}}}+{\bm{V}^{h}_{\mathrm{d}}}$
and ${\bm{v}}^{h}\in{\bm{V}^{h}}$, such that
${\bm{u}}^{h}(0)={\bm{v}}^{h}(0)=0$ in $\bm{H}$ and

$\displaystyle\begin{aligned}
\frac{\mathrm{d}}{\mathrm{d}t}\,{\langle{\bm{v}}^{h},{\bm{\varphi}}^{h}\rangle}_{\bm{H}}+\int_{\Omega}\sigma({\bm{u}}^{h}):\nabla{\bm{\varphi}}^{h}\,\mathrm{d}{\bm{x}}\hskip
42.67912pt&\\\
+\frac{1}{\gamma}\,\Big{\langle}\left[{\bm{u}}^{h}\cdot{\bm{n}}-\gamma{\bm{n}}\cdot\big{(}\sigma({\bm{u}}^{h}){\bm{n}}\big{)}\right]_{+},{\bm{\varphi}}^{h}\cdot{\bm{n}}\Big{\rangle}_{L^{2}({\Gamma_{\mathrm{tb}}})}&=0\quad\forall{\bm{\varphi}}^{h}\in{\bm{V}^{h}_{\mathrm{d}}}\text{
in }(0,T],\\\
\frac{\mathrm{d}}{\mathrm{d}t}\,{\langle{\bm{u}}^{h},{\bm{\psi}}^{h}\rangle}_{\bm{H}}-{\langle{\bm{v}}^{h},{\bm{\psi}}^{h}\rangle}_{\bm{H}}&=0\quad\forall{\bm{\psi}}\in{\bm{V}^{h}}\text{
in }(0,T].\end{aligned}$ (5)

The Nitsche parameter is chosen as
$\nicefrac{{1}}{{\gamma}}=\nicefrac{{\gamma^{0}}}{{h}}$ for a parameter
$\gamma^{0}>0$, see e.g., [5]. System (5) is discretized in time using the
backward Euler method.

##  4 Neural network approximation

To solve the optimization problem (4) numerically, we have to compute solution
pairs $({\bm{u}}^{h},{\bm{v}}^{h})$ to (5) for many parameters
$p=(E,\nu)\in{\mathscr{P}_{\mathsf{ad}}}$. For that reason we introduce a
neural network with the goal to approximate the FE solution map
$({\bm{u}}^{h},{\bm{v}}^{h})=\mathcal{S}^{h}(p)$. We assume that (5) can be
solved offline for a certain number of parameter pairs
$p^{(i)}=(E^{(i)},\nu^{(i)})$, $i=1,\ldots,n$, which is used to train the
network.

We will use a fully connected dense network, which can be defined by the
following map

$\displaystyle
N(W,B;p)=W^{L}(\sigma(W^{L-1}+\dots\sigma(W^{1}p+b^{1})\dots+b^{L-1})+b^{L},$

where $B=\\{b^{l}\in{\mathbb{R}}^{n_{l}}\,|\,l\in\\{1,\dots,L\\}\\}$ is a set
of biases and $W=\\{W^{l}\in{\mathbb{R}}^{n_{l}\times
n_{l-1}}\,|\,l\in\\{1,\dots,L\\}\\}$ a set of matrices, which are the
parameters of the neural network. An illustration of a (small) network is
shown in Figure 2. The network receives a set
$\\{p^{(i)}\\}_{i=1}^{n}\subset{\mathscr{P}_{\mathsf{ad}}}$ of parameters as
input and shall deliver the observation
$C^{(i)}=\mathcal{C}({\bm{u}}^{h}_{i})\in{\mathbb{R}}^{m}$,
${\bm{u}}^{h}_{i}=\mathcal{S}^{h}(p^{(i)})$, as an output. In the numerical
examples, we will use the softplus activation function
$\sigma(s)=\ln(e^{s}+1)$.

To train the network, we use a least-squares loss function, defined by

$\displaystyle\mathcal{L}(W,B)=\frac{1}{2n}\sum_{i=1}^{n}{\|C^{(i)}-N(W,B;p^{(i)})\|}_{{\mathbb{R}}^{m}}^{2}.$

Both input and output data are normalized around the respective mean values
and with respect to the variances before they are used to train the network,
e.g., for the $j$-parameter $p_{j}$ by
$\tilde{p}_{j}^{(i)}\colonequals\nicefrac{{(p_{j}^{(i)}-\mathrm{mean}_{j})}}{{\delta_{j}}}$
with $\mathrm{mean}_{j}=\sum_{i=1}^{n}p^{(i)}/n$ and variance $\delta_{j}$.

Training. We apply the neural network to the model problem (5). To train the
neural network, (5) has been solved for $100\,000$ equally distributed
parameter pairs ($E,\nu$) in $[5\cdot 10^{10},10^{11}]\times[0.3,0.4]$ using
$1000$ equally distributed values for $E$ and $100$ equally distributed values
for $\nu$. For validation, we create another set of $19\,483$ data pairs
corresponding to randomly created parameters
$p\in{\mathscr{P}_{\mathsf{ad}}}$. We use the deep-learning API keras, which
is based on tensorflow [4]. A dense network with two hidden layers of size
$200$ and $100$ is applied. For the training, we use a combination of the
stochastic gradient algorithm and the Adam algorithm, see [9], for in total
500 epochs. More precisely, we execute both algorithms for 50 epochs each and
use the result of the one, which performed better in terms of a reduction in
the loss function evaluated in the validation data. We start with a learning
rate of $10^{-3}$, which is increased gradually up to $5\cdot 10^{-2}$ in the
last epochs. As a result the loss function decreased by a factor of roughly
$10^{4}$ within 500 epochs. The total time required for training was
approximately $107,6$ seconds.

Approximation properties. In order to test the neural network approximation we
compare the two functionals

$\displaystyle\begin{aligned}
\mathcal{F}^{h}(E,\nu)&\colonequals\frac{1}{2}\,{\|{\mathcal{C}}({\bm{u}}^{h})-{{\bm{u}}^{\mathrm{obs}}}\|}^{2}_{{\mathbb{R}}^{m}}&&\text{for
}{\bm{u}}^{h}=\mathcal{S}^{h}(E,\nu),\\\
\mathcal{F}^{N}(E,\nu)&\colonequals\frac{1}{2}\,{\|{\mathcal{C}}({\bm{u}}^{N})-{{\bm{u}}^{\mathrm{obs}}}\|}^{2}_{{\mathbb{R}}^{m}}&&\text{for
}{\bm{u}}^{N}=N(W,B;(E,\nu)).\end{aligned}$ (6)

We assume that the parameters we would like to find in our parameter
estimation problem are given by

$\displaystyle p_{1}\colonequals(E_{1},\nu_{1})\quad\text{with
}E_{1}=7.513\cdot 10^{9},~{}\nu_{1}=0.3547$ (7)

and define the observed data by means of the FE solution
${\bm{u}}^{h}_{1}=\mathcal{S}^{h}(p_{1})$ of (5) as
${{\bm{u}}^{\mathrm{obs}}}=\mathcal{C}({\bm{u}}_{1}^{h})$. In Figure 3, we
compare $\mathcal{F}^{h}$ and its neural network approximation
$\mathcal{F}^{N}$ on ${\mathscr{P}_{\mathsf{ad}}}$. We observe that both
graphs nearly coincide. Moreover, $\mathcal{F}^{N}$ seems to be convex despite
the non-linear activation function $\sigma$ of the network, which is
beneficial for numerical optimization. The minimum is attained very close to
$p_{1}$, see also the numerical results for the parameter estimation problem
given below.

Figure 3: Comparison of the FE objective $\mathcal{F}^{h}(p_{1})$ and its
neural network approximation $\mathcal{F}^{N}(p_{1})$ for
${{\bm{u}}^{\mathrm{obs}}}=\mathcal{C}({\bm{u}}_{1}^{h})$ and
${\bm{u}}^{h}_{1}=\mathcal{S}^{h}(p_{1})$ with $p_{1}$ defined in (7).

##  5 Neural-network based gradient algorithm

Gradient-type algorithms typically require frequent solutions of the PDE (5)
and its sensitivities. In order to avoid these, we approximate the operator
$p\mapsto\mathcal{C}(\mathcal{S}^{h}(p))$) by means of a neural network, i.e.,
we seek the minimizer of the function $\mathcal{F}^{N}(p)$ defined in (6).
Therefore, we need to compute the derivative of the network output $N(W,B;p)$
with respect to $p\in{\mathscr{P}_{\mathsf{ad}}}$. This can be computed by
means of a backpropagation:

###### Lemma 5.1

Let $N(W,B;p)$ be the dense network defined in Section 3. Moreover, let
$z^{l}=W^{l}a^{l-1}+b^{l}$ and $a^{l}=\sigma_{l}(z^{l})$ be the output of the
$l$-the layer of the network. It holds for the derivatives
$\rho^{l}=(\rho^{l}_{j})_{j=1}^{n_{l}}$,
$\rho^{l}_{j}=\partial_{a^{l}_{j}}\mathcal{F}^{N}$:

$\displaystyle\rho^{L}=a^{L}-y\quad\text{ and
}\quad\rho^{l}=(W^{l+1})^{\top}\big{(}\rho^{l+1}\circ\sigma_{l+1}^{\prime}(z^{l+1})\big{)}\,\text{
for }l=0,\dots,L-1,$

where $(u\circ v)_{j}=u_{j}v_{j}$. In particular, we obtain the derivative
$\rho_{j}^{0}\colonequals\partial_{p_{j}}\mathcal{F}^{N}$ of $\mathcal{F}^{N}$
with respect to the network input $p=(p_{1},p_{2})$.

Proof. The proof for $l\geq 1$ is standard and can be found in any textbook on
neural networks, see, e.g., [8]. The extension to $l=0$ is completely
straight-forward.

We use Lemma 5.1 to apply a gradient algorithm to find the minimum of
$\mathcal{F}^{N}(p)$. The resulting algorithm is given as Algorithm 1. For the
line search, we use an Armijo-Goldstein line search in the numerical results
given below.

Algorithm 1 (Neural-network based gradient algorithm)

1:procedure Gradient algorithm(network $N$, $p^{0}$)

2: $p\leftarrow p^{0}$

3: while Termination criterion is not satisfied do

4: $\Delta p\leftarrow$ calculate $\nabla_{p}\mathcal{F}^{N}(p)$ % using Lemma
5.1

5: $\eta\leftarrow$ Line_Search$(\mathcal{F}^{N}(\cdot),p,\Delta p)$

6: $p\leftarrow p-\eta\Delta p$ return $p$

##  6 Numerical examples

We use again the example introduced in Section 4 with
${{\bm{u}}^{\mathrm{obs}}}=\mathcal{C}({\bm{u}}^{h}_{1})$ and
${\bm{u}}_{1}^{h}=\mathcal{S}^{h}(p_{1})$ for $p_{1}$ specified in (7). The
neural network $N(W,B;\cdot)$ is trained offline as described in Section 4. An
evaluation of the neural network $N(W,B;p)$ for a new parameter pair $p$
requires then only $29.8$ $\mu s$ in average, while our FE solver needs $243$
seconds to solve (5).

Due to the convexity of the functional $\mathcal{F}^{N}$ (see Fig. 3) a
regularization is not required, i.e., we set $\alpha=0$ in (4). We stop the
gradient algorithm, when the initial value of the norm of the gradient
$\|\nabla_{p}{\cal F}^{N}(p)\|_{{\mathbb{R}}^{m}}$ is reduced by a factor of
$10^{-8}$. The behavior of the gradient algorithm is visualized in Figure 4,
left sketch. We observe that the gradient algorithm converges quickly towards
the minimizer. This is also confirmed in Figure 5, where the objective
function $\mathcal{F}^{N}$ and its gradient are visualized within the
iteration.

Figure 4: Convergence behavior of the neural-network based gradient algorithm
applied to problem (6) to the correct solutions $p_{1}$ (7) (left sketch) and
$p_{2}$ (right sketch). The correct solutions are marked with blue dots, the
iterates in red.

Figure 5: Left: Objective function ${\cal F}^{N}$, right: Norm of the gradient
$\|\nabla_{p}{\cal F}^{N}\|_{{\mathbb{R}}^{2}}$ for the first example
corresponding to parameters (7).

The resulting parameters $p$ are (after re-normalisation) $\bar{E}_{1}\approx
7.51385\cdot 10^{10}$, $\bar{\nu}_{1}\approx 0.354768$, which yields a
relative error of approximately $1.1\cdot 10^{-3}$. This seems to be the error
resulting from the network approximation.

Finally, we apply the algorithm in the same way to a second point
$p_{2}=(E_{2},\nu_{2})$, where $E_{2}=7.45\cdot 10^{10}$ and $\nu_{2}=0.3481$.
The convergence behavior is shown in Figure 4, right sketch. The resulting
parameters are $\bar{E}_{2}\approx 7.45071\cdot 10^{10}$,
$\bar{\nu}_{2}\approx 0.348164$, which gives a relative error of around
$3.8\cdot 10^{-3}$.

##  7 Conclusion

We have presented a non-intrusive gradient algorithm based on a neural-network
approximation of the underlying PDE. First results for a dynamic contact
problem look very promising. Despite the non-linear activation functions
involved in the network approximation, the approximated cost function ${\cal
F}^{N}$ seems to conserve the convexity of the function ${\cal F}^{h}$. We
think that the algorithm can be particularly interesting in cases, when the
PDE solvers are computationally very expensive or not accessible at all for
the computation of sensitivities.

## References

  * [1] P. Alart and A. Curnier.  A mixed formulation for frictional contact problems prone to Newton like solution methods.  Computer Methods in Applied Mechanics and Engineering, 92(3):353–375, 1991.
  * [2] H. Antil, H.C. Elman, A. Onwunta, and D. Verma.  A deep neural network approach for parameterized pdes and bayesian inverse problems.  Machine Learning: Science and Technology, 4(3):035015, 2023.
  * [3] S. Arridge, P. Maaß, O. Öktem, and C.B. Schönlieb.  Solving inverse problems using data-driven models.  Acta Numerica, 28:1–174, 2019.
  * [4] F. Chollet et al.  Keras.  https://keras.io, 2015.
  * [5] F. Chouly and P. Hild.  A Nitsche-based method for unilateral contact problems: numerical analysis.  SIAM Journal on Numerical Analysis, 51(2):1295–1307, 2013.
  * [6] C. Eck, J. Jarusek, and M. Krbec.  Unilateral Contact Problems: Variational Methods and Existence Theorems, volume 270.  CRC Press, 2005.
  * [7] P. Grohs, F. Hornung, A. Jentzen, and P. Zimmermann.  Space-time error estimates for deep neural network approximations for differential equations.  Advances in Computational Mathematics, 49(4), 2023.
  * [8] C.F. Higham and D.J. Higham.  Deep learning: An introduction for applied mathematicians.  SIAM Review, 61(4):860–891, 2019.
  * [9] D.P. Kingma and J. Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2014.
  * [10] M. Raissi, P. Perdikaris, and G.E. Karniadakis.  Physics-informed neural networks: a deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations.  Journal of Computational Physics, 378:686–707, 2019.


# Corrected Weight Functions for Stellar Oscillation Eigenfrequencies

R. H. D. Townsend Department of Astronomy, University of Wisconsin-Madison,
475 N Charter St, Madison, WI 53706, USA<EMAIL_ADDRESS>S. D.
Kawaler Department of Physics and Astronomy, Iowa State University, 2323
Osborn Drive, Ames, IA 50011, USA<EMAIL_ADDRESS>
###### Abstract

Kawaler et al. (1985) present a variational expression for the
eigenfrequencies associated with stellar oscillations. We highlight and
correct a typographical error in the weight functions appearing in these
expressions, and validate the correction numerically.

Asteroseismology (73) — Pulsating variable stars (1307) — Stellar oscillations
(1617) — Astronomical Methods (1043)

††facilities: We have made extensive use of NASA’s Astrophysics Data System
Bibliographic Services.††software: Astropy (Astropy Collaboration et al.,
2013, 2018, 2022), GYRE (Townsend & Teitler, 2013; Townsend et al., 2018;
Goldstein & Townsend, 2020), Matplotlib (Hunter, 2007)

tablenum

##  1 Introduction

Under the assumption of a vanishing surface density, the equations governing
linear, adiabatic, non-radial stellar oscillations are self-adjoint and
therefore obey a variational principle (see, e.g. Unno et al., 1989).
Variational expressions for mode eigenfrequencies serve as a valuable
computational check on numerical solutions to the oscillation equations. Via
so-called weight functions, they also quantify how eigenfrequencies are
sensitive to different parts of a star’s internal structure; they are
therefore key tools in determining the asteroseismic probing potential of each
mode.

Kawaler et al. (1985) present a variational expression that has been adopted
in many studies of stellar oscillations (e.g., Kawaler et al., 1986; Brassard
et al., 1992a, b; Córsico & Benvenuto, 2002; Timmes et al., 2018). In the
present paper we highlight and correct a typographical error in one of the
weight functions appearing in their expression, and use numerical calculations
to demonstrate the validity of the correction.

##  2 Weight Functions

The Kawaler et al. (1985) variational expression for the eigenfrequency
$\sigma(\mathbf{y})$ of a mode is

$\sigma^{2}(\mathbf{y})=\frac{\int_{0}^{R}\left[\mathcal{C}(\mathbf{y},r)+\mathcal{N}(\mathbf{y},r)+\mathcal{G}(\mathbf{y},r)\right]\rho
r^{2}\operatorname{d}\\!{r}}{\int_{0}^{R}\mathcal{T}(\mathbf{y},r)\rho
r^{2}\operatorname{d}\\!{r}}.$ (1)

Here, $\mathbf{y}(r)=\\{y_{1}(r),\ldots,y_{4}(r)\\}$ is the vector of
eigenfunctions introduced by Dziembowski (1971), $r$ is the radial coordinate
(with $R$ the stellar radius), and $\rho$ the density. The three weight
functions are

$\displaystyle\mathcal{C}(\mathbf{y},r)$
$\displaystyle=g^{2}\ell(\ell+1)S_{\ell}^{-2}(y_{2}-y_{3})^{2},$ (2)
$\displaystyle\mathcal{N}(\mathbf{y},r)$ $\displaystyle=r^{2}N^{2}y_{1}^{2},$
(3) $\displaystyle\mathcal{G}(\mathbf{y},r)$
$\displaystyle=-\frac{gr}{U}\left[y_{4}+\ell(\ell+1)y_{3}\right]^{2}.$ (4)

where $\ell$ is the model harmonic degree, $g$ is the gravity, $N$ the Brunt-
Väisälä frequency, $S_{\ell}$ the Lamb frequency, and $U$ the usual homology
invariant. The function

$\mathcal{T}(\mathbf{y},r)=r^{2}\left[y_{1}^{2}+\ell(\ell+1)\left(\frac{g}{r\sigma^{2}}\right)^{2}y_{2}^{2}\right]$
(5)

appearing in the denominator is proportional to the differential inertia of
the mode.

By re-deriving these expressions from equations (13.13) and (13.15) of Unno et
al. (1979)111The equivalent equations in Unno et al. (1989) are (14.16) and
(14.19)., we find that the correct expression for $\mathcal{G}$ is

$\mathcal{G}(\mathbf{y},r)=-\frac{gr}{U}\left[y_{4}+(\ell+1)y_{3}\right]^{2};$
(6)

comparing this against equation (4), there is an extra factor of $\ell$
multiplying the $y_{3}$ term in the latter. We have confirmed that this error
is typographical in nature; the spurious factor was not included in the
calculations presented in Kawaler et al. (1985) and follow-up papers by the
same authors.

###  2.1 Numerical Validation

Figure 1: The relative frequency error $\Delta$ between numerical and
variational frequencies of the $n=3$ polytrope, plotted as a function of mode
radial order $n_{\rm pg}$ using the original (equation (4); left) and
corrected (equation (6); right) expressions for $\mathcal{G}$.

To validate our analysis we explore how well equation (1) can reproduce the
numerical eigenfrequencies of an $n=3$ polytrope model, using both original
and corrected expressions for $\mathcal{G}$. We apply the GYRE oscillation
code (Townsend & Teitler, 2013; Townsend et al., 2018; Goldstein & Townsend,
2020) to calculate adiabatic eigenfrequencies and eigenfunctions for low-
radial order, $\ell=0,\ldots,4$ modes of this model. GYRE’s VACUUM boundary
condition (where the Lagrangian pressure perturbation vanishes at the surface)
is adopted, and the oscillation equations are discretized using the COLLOC_GL2
(second-order Gauss-Legendre collocation) scheme on a radial grid of
$\sim$1000$$ points.

For each mode found we evaluate the relative frequency error

$\Delta=\left|\frac{\sigma(\mathbf{y})-\sigma}{\sigma}\right|,$ (7)

where $\sigma$ is the numerical eigenfrequency reported by GYRE. Figure 1
plots $\Delta$ as a function of mode radial order $n_{\rm pg}$ (defined in the
Eckart-Osaki-Scuflaire classification scheme; see Unno et al., 1989), for
calculations based on the original (equation 4) and corrected (equation 6)
expressions for $\mathcal{G}$. In the former case (left panel) the errors are
significant for $\ell\geq 2$ modes, reaching $\Delta\gtrsim 0.1$ at small
$|n_{\rm pg}|$. In the latter (right panel) the frequency errors are small for
all modes, $\Delta\lesssim${10}^{-4}$$, and arise primarily from the
truncation error of GYRE’s numerical scheme. This confirms the validity of the
corrected expression.

The original expression for $\mathcal{G}$ happens to yield small $\Delta$ for
both radial ($\ell=0$) and dipole ($\ell=1$) modes. When $\ell=1$ the spurious
factor of $\ell$ in equation (4) makes no difference to the value of
$\mathcal{G}$. When $\ell=0$, this factor introduces an additive offset

$\delta\mathcal{G}(\mathbf{y})=\frac{gr}{U}y_{3}\left(y_{3}+2y_{4}\right)$ (8)

to $\mathcal{G}$; however, integration by parts, together with the fact that
$y_{3}=0$ at the outer boundary for radial modes, can be used to demonstrate
that this offset makes no net contribution to the numerator of the variational
expression (1).

## Acknowledgments

This work has been supported by NSF grants ACI-1663696, AST-1716436 and
PHY-1748958, and NASA grant 80NSSC20K0515.

## References

  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123
  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167
  * Brassard et al. (1992a) Brassard, P., Fontaine, G., Wesemael, F., & Hansen, C. J. 1992a, ApJS, 80, 369
  * Brassard et al. (1992b) Brassard, P., Fontaine, G., Wesemael, F., & Tassoul, M. 1992b, ApJS, 81, 747
  * Córsico & Benvenuto (2002) Córsico, A. H., & Benvenuto, O. G. 2002, Ap&SS, 279, 281
  * Dziembowski (1971) Dziembowski, W. A. 1971, Acta Astron., 21, 289
  * Goldstein & Townsend (2020) Goldstein, J., & Townsend, R. H. D. 2020, ApJ, 899, 116
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science & Engineering, 9, 90
  * Kawaler et al. (1985) Kawaler, S. D., Winget, D. E., & Hansen, C. J. 1985, ApJ, 295, 547
  * Kawaler et al. (1986) Kawaler, S. D., Winget, D. E., Hansen, C. J., & Iben, I., J. 1986, ApJ, 306, L41
  * Timmes et al. (2018) Timmes, F. X., Townsend, R. H. D., Bauer, E. B., et al. 2018, ApJ, 867, L30
  * Townsend et al. (2018) Townsend, R. H. D., Goldstein, J., & Zweibel, E. G. 2018, MNRAS, 475, 879
  * Townsend & Teitler (2013) Townsend, R. H. D., & Teitler, S. A. 2013, MNRAS, 435, 3406
  * Unno et al. (1989) Unno, W., Osaki, Y., Ando, H., Saio, H., & Shibahashi, H. 1989, Nonradial oscillations of stars, 2nd edn. (University of Tokyo Press)
  * Unno et al. (1979) Unno, W., Osaki, Y., Ando, H., & Shibahashi, H. 1979, Nonradial oscillations of stars, 1st edn. (University of Tokyo Press)