# Exploring descriptors for titanium microstructure via digital fingerprints
from variational autoencoders

Michael D. White1,111Corresponding author
Email address<EMAIL_ADDRESS>(M.D. White), Gowtham
Nimmal Haribabu1,3, Jeyapriya Thimukonda
Jegadeesan3, Bikramjit Basu3, Philip J. Withers1,2 and Chris P. Race2,4

1Department of Materials, University of Manchester, Manchester, UK, M13 9PL
2Henry Royce Institute, University of Manchester, Manchester, UK, M13 9PL
3Materials Research Centre, Indian Institute of Science, Bangalore, India,
560012
4Department of Materials Science and Engineering, University of Sheffield,
Sheffield, UK, S1 3JD

###### Abstract

Microstructure is key to controlling and understanding the properties of
metallic materials, but traditional approaches to describing microstructure
capture only a small number of features. To enable data-centric approaches to
materials discovery, to allow efficient storage of microstructural data and to
assist in quality control in metals processing, we require more complete
descriptors of microstructure. The concept of _microstructural fingerprinting_
, using machine learning (ML) to develop quantitative, low-dimensional
descriptors of microstructures, has recently attracted significant attention.
However, it is difficult to interpret conclusions drawn by ML algorithms,
which are often referred to as “black boxes”. For example, convolutional
neural networks (CNNs) can be trained to make predictions about a material
from a set of microstructural image data, but the feature space that is
learned is often used uncritically and adopted without any validation.

Here we explore the use of variational autoencoders (VAEs), comprising a pair
of CNNs, which can be trained to produce microstructural fingerprints in a
continuous latent space. The VAE architecture also permits the reconstruction
of images from fingerprints, allowing us to explore how key features of
microstructure are encoded in the latent space of fingerprints. We develop a
VAE architecture based on ResNet18 and train it on two classes of Ti-6Al-4V
optical micrographs (bimodal and lamellar) as an example of an industrially
important alloy where microstructural control is critical to performance. The
latent/feature space of fingerprints learned by the VAE is explored in several
ways, including by supplying interpolated and randomly perturbed fingerprints
to the trained decoder and via dimensionality reduction to explore the
distribution and correlation of microstructural features within the latent
space of fingerprints.

We show that the fingerprints generated via the trained VAE exhibit smooth,
interpolable behaviour with stability to local perturbations, supporting their
suitability as general purpose descriptors for microstructure. We also show
that key properties of the microstructures (volume fraction and grain size)
are strongly correlated with position in the encoded feature space, supporting
the use of VAE fingerprints for quantitative exploration of process-structure-
property relationships.

##  1 Introduction

In the field of metallurgy and materials science, process-structure-property
(PSP) linkages are instrumental in guiding material design for targeted
applications [1, 2]. Despite the importance of PSP linkages, a rigorous
mathematical framework is not currently available for systematic analysis in
this context [3]. The key problem is that, whilst compositional processing
parameters and measured properties are inherently described by numbers,
microstructure lacks a comprehensive numerical descriptor. A central
impediment is that the characteristic microstructural features exhibit
heterogeneity over a wide range of size, spatial, and temporal scales. From
the application standpoint, it is important to identify a subset of salient
measures of internal structure that can be tracked through a material’s
processing history, that capture the dominant influences on the targeted
properties. In conventional microstructural analysis, some quantitative
methods are used, but these generally rely on metrics applied to image data,
such as grain size, phase fraction and correlation coefficients. This can
provide some crucial information, but the metrics which are suitable in each
case will depend on the morphology of the microstructure. For example,
interlamellar spacing is useful to describe lamellar microstructures, but is
redundant when considering an equiaxed microstructure. In any case, such
measures embody only a tiny fraction of the information contained in a
microstructural image. A quantitative description of microstructure that is
independent of morphology (a microstructural fingerprint) and embodies a full
range of features would enable direct comparisons between microstructures with
different morphologies, and provide new methods for constructing PSP
relationships [3, 4], or for quantifying the deviation of a given
microstructure from some ideal standard in a quality control procedure.

In the past two decades, reasonable progress has been achieved in the
quantification and low-dimensional representation of microstructures, over
larger processing and material compositional windows [5, 6, 7]. These advances
have been possible with the use of concepts and toolsets from data science and
informatics [8, 9]. The first potential benefit is realised with the use of
automated feature engineering of the hierarchical material structures, e.g.,
establishing low-dimensional representations that provide maximum value in
capturing high-fidelity PSP linkages. A systematic and comprehensive
quantification of the material structure is now possible, for example by
performing statistical analysis of image features [10] or by combining the
physics-inspired framework of $n$-point spatial correlations ($n$-point
statistics) with machine learning approaches, such as principal component
analysis (PCA) [11].

A machine learning approach that has not yet been fully explored, in the
context of microstructural fingerprinting, is variational autoencoders (VAEs),
which comprise a pair of convolutional neural networks (CNNs), referred to
individually as the encoder and the decoder. VAEs were first introduced by
Kingma and Welling [12]. Initial applications focussed on generating images of
individual objects, particularly human faces, utilising datasets such as the
CelebA dataset [13]. More recently, attention has shifted towards machine
learning applications in materials science, such as quantification of
microstructure and the development of new PSP relationships. The key
difference is that microstructural image data are such that the entire image
field contains potentially meaningful information, rather than an image of a
foreground object of interest and some arbitrary background. Attempts have
been made to construct VAEs for texture embedding. One example is TextureVAE,
which consists of a VGG19 encoder and a decoder comprising 4 convolution
blocks [14]. The network was tested on various microstructures and latent
dimensions were artificially varied to visualise their effect on the resulting
reconstructions. Further exploration of the encoded space has also been shown
to carry the potential for material property prediction [15, 16]. Conversely,
the ability to generate synthetic microstructures from specified properties
was demonstrated [17]. Dimensionality reduction of the encoded space has also
been shown to provide meaningful visualisation of the space, enabling property
prediction from microstructural image data [18] and identification of new
processing routes for target orientation distributions [19].

Here, we employ a deep residual block VAE architecture, based on ResNet18
[20], to encode optical micrograph data from two classes of Ti-6Al-4V; a
lamellar microstructure and a bimodal equiaxed microstructure. The encoded
space is explored via paths and Gaussian permutations, as a tool for
generating artificial microstructures. The latent space in further explored
through dimensionality reduction via $t$-SNE, with analysis of morphological
metric distributions across the space. Support vector regression (SVR) is then
applied to correlate fingerprints contained in the latent space with the same
morphological metrics to quantify the trends visualised by the $t$-SNE.

Our titanium alloy dataset is representative of the microstructures of an
important class of industrial materials in which control of microstructure is
key to delivering the required in-service performance and in which quality
assurance of material (and hence microstructure) is a critical part of
manufacturing processes. Furthermore, the nature of the dataset allows us to
evaluate the performance of the VAE with respect to several aspects common to
a broad range of metallurgical challenges:

  1. 1.

the descriptors should transparently handle a range of microstructures (here
we have widely varying grain morphology);

  2. 2.

the space of descriptors (the latent space of the VAE) should be well-behaved,
in the sense that the fingerprints should vary smoothly with changes to the
microstructure and vice versa;

  3. 3.

the representation in the latent space should be interpretable in terms of key
features of the microstructure (or properties of the material), which is to
say that the fingerprints should show statistical correlation with features
and properties.

We explore these aspects of the VAE performance for our titanium dataset.

##  2 Materials and Methods

###  2.1 Dataset

The LightForm Ti-6Al-4V alloy bimodal/lamellar dataset (LFTi64BL) is an open
access dataset, curated within the LightForm group at the University of
Manchester, and can be accessed via Zenodo [21]. The dataset comprises 40
optical 8 bit micrographs of resolution $1292\times 968$ pixels, containing
equal numbers of bimodal equiaxed and lamellar microstructures. An example
from each classification is shown in Figure 1.

(a) Bimodal

(b) Lamellar

Figure 1: Representative optical micrographs of Ti-6Al-4V from LFTi64BL
dataset.

This dataset contains a relatively small number of images for machine learning
purposes. To expand the dataset, patches were extracted with random rotations
and reflections applied. The expanded LFTi64BL dataset contains 1,000 patches
from each image, resulting in 40,000 patches of resolution $256\times 256$
pixels.

###  2.2 Greyscale Normalisation

Each image was white balanced by clipping greyscale intensities outside the
90${}^{\text{th}}$ percentile and remapping to the range $[0,1]$. Figure 2
shows distributions of mean greyscale intensities across the bimodal and
lamellar datasets within LFTi64BL separately, before and after normalisation.

(a) Bimodal

(b) Lamellar

Figure 2: Distribution of mean greyscale intensities across the LFTi64BL
dataset before and after normalisation, for bimdoal and lamellar
microstructures separately.

Prior to white-balancing, the variation in mean greyscale intensity is heavily
influenced by fluctuations in lighting conditions during image capture.
However, after white-balancing, greyscale intensity is more normally
distributed across the dataset. This is now indicative of the distribution of
phase fractions, with lower mean greyscale intensity corresponding to a higher
volume fraction of the $\upbeta$ phase.

###  2.3 Variational Autoencoders (VAEs)

Variational autoencoders (VAEs) are a tool for encoding image data into a
compressed format (or fingerprint) that preserves morphological features. They
comprise a pair of convolutional neural networks (CNNs), referred to
individually as the encoder and decoder. Fingerprints output by the trained
encoder can be fed into downstream tasks such as image classification and
property prediction. Figure 3 provides a visual interpretation of the general
VAE architecture.

Figure 3: Schematic of general VAE architecture that takes an image $x$ as
input, encodes the image into a fingerprint $e(x)$ and is trained to
reconstruct the input image with the mapping $d(e(x))$.

Suppose we have an image, $x\in\mathbb{R}^{m_{1}\times m_{2}}$. The encoder
takes the image $x$ as an input and computes a compressed representation
$z\in\mathbb{R}^{D}$, such that $D\ll m_{1}m_{2}$, where $z=e(x)$ and $D$ is a
tunable parameter that denoted the dimensionality of the encoded space. The
decoder then takes the encoded vector $z$ as an input and aims to reconstruct
the input image, $x$. We denote the reconstruction as $\hat{x}=d(e(x))$, where
$\hat{x}\in\mathbb{R}^{m_{1}\times m_{2}}$. The encoded space learned is
continuous and normally distributed.

To train a VAE, we must have a set of images,
$X=\\{x_{1},x_{2},\dots,x_{N}\\}$, where $N$ is the total number of images in
the dataset. The dataset is split into two subsets,
$X_{\text{train}}\in\mathbb{R}^{N_{\text{train}}\times m_{1}\times m_{2}}$ and
$X_{\text{eval}}\in\mathbb{R}^{N_{\text{eval}}\times m_{1}\times m_{2}}$,
where $X_{\text{train}}$ is used to train the VAE, $X_{\text{eval}}$ is used
for evaluating the VAE on unseen data, $N_{\text{train}}+N_{\text{eval}}=N$,
$X_{\text{train}}\cup X_{\text{eval}}\equiv X$ and $X_{\text{train}}\cap
X_{\text{eval}}\equiv\varnothing$. The training set is split into batches and
fed into the VAE, one batch at a time. Ultimately, $d(e(x_{i}))$ for $x_{i}\in
X_{\text{train}}$ is computed and some loss function is used as a metric for
assessing the encoded representation and reconstruction quality. The loss
function is then used as an input to an optimiser. The Adam optimiser [22] is
currently the state of the art and is utilised throughout all models discussed
herein. The optimiser updates weights and biases in both the encoder and
decoder from a single loss function after each batch operation. The model can
then be evaluated on $X_{\text{eval}}$ to measure the potential for transfer
learning, but metrics calculated on these images in the evaluation set are not
utilised for updating any weights or biases. Encoded representations,
generated from the trained VAE on microstructural image data, can be
considered as a signature, or microstructural fingerprint (as described in
[3]), and will be referred to as such throughout this paper.

###  2.4 Loss Functions

During VAE training, a loss function is periodically computed on the model, to
determine performance and provide the inputs for the optimiser to update the
weights and biases in the networks, with the aim of minimising the loss
function. This is typically a combination between a measure in the spacial
domain of the reconstruction accuracy and imposing a restriction on the
encoded space to be normally distributed [12]. To measure how normally
distributed the encoded space is, the Kullback-Leibler (KL) divergence,
$D_{\text{KL}}$, can be determined between each encoded representation and the
unit normal distribution $\mathcal{N}(0,1)$ [23].

The reconstruction accuracy can be measured in several ways. A popular method
is the mean squared error (MSE) [24], which is a distance metric in the
spacial domain, given by

$\text{MSE}=\frac{1}{m}\sum_{j=1}^{m}(x_{j}-\hat{x}_{j})^{2},$

where $m=m_{1}m_{2}$ is the total number of pixels in the image and the
$x_{j},\hat{x}_{j}$ denote pixels in the input and reconstruction,
respectively.

Another metric that operates in the spacial domain is the binary cross-entropy
(BCE), denoted here as $\mathcal{L}_{b}$, which is a measure of negative log
likelihood and is given by

$\mathcal{L}_{b}=-\sum_{i=1}^{m}x_{i}\log\hat{x}_{i}-\sum_{i=1}^{m}(1-x_{i})\log(1-\hat{x}_{i}).$

The issue with minimising these loss functions is that they generally result
in blurry reconstructions. With the aim of minimising blur, a loss function on
the frequency domain was proposed in [25]. This requires calculating the 2D
Fourier transform of both the input and reconstruction. The spectral loss,
$\mathcal{L}_{f}(x,\hat{x})$, can then be given by the MSE between the 2D
Fourier transforms, i.e.,

$\mathcal{L}_{f}(x,\hat{x})=\frac{1}{m}\sum_{j=1}^{m}\left((\text{Im}\\{\mathcal{F}(x)_{j}\\}-\text{Im}\\{\mathcal{F}(\hat{x})_{j}\\})^{2}+(\text{Re}\\{\mathcal{F}(x)_{j}\\}-\text{Re}\\{\mathcal{F}(\hat{x})_{j}\\})^{2}\right),$

where $\mathcal{F}$ denotes the 2D FFT, $\text{Im}\\{\mathcal{F}\\}$ denotes
the imaginary part of $\mathcal{F}$ and $\text{Re}\\{\mathcal{F}\\}$ the real
part. The total loss function to be minimised by the optimiser is then given
by

$\mathcal{L}=\alpha\mathcal{L}_{b}(x,\hat{x})+(1-\alpha)\mathcal{L}_{f}(x,\hat{x}),$
(1)

where $\alpha\in[0,1]$ is a tunable hyperparameter.

###  2.5 ResNet18-VAE

ResNet [20] is a deep CNN composed of residual blocks and was initially
proposed for classification of the ImageNet dataset [26], which is a dataset
containing over 1 million natural images with 1000 classifications. As such,
the final layer of a standard ResNet is a 1000-dimensional fully connected
layer, where the output from each node corresponds to a probability for each
ImageNet classification. The depth of the network can be controlled by varying
the number of layers within each block and the total number of blocks. Here,
we consider ResNet18, which contains 8 residual blocks, each comprised of 2
convolution layers with subsequent ReLU activation and batch normalisation.
Each block has the potential to be effectively skipped by summing the input
received at each block with the output from that block, after convolution.

To convert ResNet18 into an encoder, it was modified by replacing the average
pool and 1000-dimensional fully connected layers at the end of the network
with a flattening of the final convolution output. This was followed by a
512-dimensional fully connected layer with tanh activation, which is then
simultaneously fed into two separate fully connected layers. Each fully
connected layer consists of $|z|$ neurons, where $|z|$ is the specified
dimension of the encoded space. This is somewhat arbitrary, but 256 was used
to generate the results presented in this paper. One of these fully connected
layers is utilised as a set of means, $\mu$, whilst the other is treated as a
set of standard deviations, $\sigma$. These are then combined into the output
$z$ as

$z=\mu+\exp(\sigma\mathcal{N}(0,1)),$

where $\mathcal{N}(0,1)$ denotes the standard normal distribution, with mean 0
and standard deviation 1. This branching into $\mu$ and $\sigma$, with
subsequent combination of the two, is what sets variational autoencoders apart
from standard autoencoders. The decoder then essentially mirrors the encoder
with transpose convolution layers replacing the standard convolution layers.
The architecture for this ResNet18 VAE is provided in Figure 4.

(a) Encoder

(b) Decoder

Figure 4: ResNet18 VAE architecture for the $256\times 256$ inputs used in the
present study.

###  2.6 Morphological Analysis

Morphological measurements were automated for the LFTi64BL dataset as metrics
for correlation with fingerprints produced by the VAE. These metrics can then
also be plotted as colour maps over across a dimensionality reduction of the
fingerprint space to visualise how such features are distributed. Each metric
requires the images to be binarised prior to measurement.

####  2.6.1 Image Binarisation

Binarised images are required to compute certain metrics on the images, such
as phase fraction and grain size. A high-pass Gaussian filter was applied to
each image, in the Fourier domain, to normalise illumination across the
images. The images were then binarised with Otsu’s thresholding method [27],
before applying area closing to remove noise from all images and
$\upalpha_{\text{s}}$ laths from the bimodal images. Pixels corresponding to
$\upalpha_{\text{p}}$ grains are assigned a label, whilst all other pixels are
labelled as 0. Figure 5 provides representative examples of the resulting
binary images. The binarised lamellar microstructures are highly accurate, but
there is some retention of $\upalpha_{\text{s}}$ laths in the bimodal
microstructures due to connectivity with the $\upalpha_{\text{p}}$ grains.

(a) Original lamellar

(b) Binary lamellar

(c) Original bimodal

(d) Binary bimodal

Figure 5: Representative optical micrograph images patches of Ti-6Al-4V before
and after binarisation.

####  2.6.2 Phase Fraction

Due to the way in which the images were binarised (described in Section
2.6.1), the phase fraction is simply determined as the sum of all the pixel
values in the binary image divided by the total number of pixels in the image.

####  2.6.3 Lamellae Direction

For the lamellar images in the LFTi64BL dataset, another descriptor we can
consider is the dominant direction, or orientation, of the lamellae, relative
to the bottom edge of the image. This can be quantified as the mean angle
subtended between the elongation direction of each lamella and the bottom edge
of the image. An erosion is applied on the binary image to isolate overlapping
grains before each grain is labelled with a unique integer. Watershed
segmentation is then applied with markers taken from the labelled image and
the initial binary image as a mask. Each grain in the watershed image is
isolated and eigenvectors are computed for each grain. The primary eigenvector
describes the direction of each lamella. The angle between the primary
eigenvector and the bottom edge of the image is then calculated for each grain
and averaged to provide the metric for lamellae direction of a given image.
Figure 6 shows the erosion and watershed method applied to a representative
micrograph and Figure 7 shows the corresponding eigenvectors for an individual
grain.

(a) Original micrograph

(b) Binary with area closing

(c) Eroded and labelled

(d) Watershed segmentation

Figure 6: Erosion and watershed methods applied to a representative micrograph
for grain isolation. Figure 7: Eigenvectors plotted on an individual lamella
for direction measurement, with the primary eigenvector shown in red and the
secondary eigenvector plotted in blue.

####  2.6.4 Bimodal Grain Size

As a metric for quantifying the bimodal microstructures in LFTi64BL, average
grain size measurements were automated following ASTM E1382-97 [28]. Random
line scans are applied on the binarised bimodal micrographs and peaks in the
derivative of the profile along the line scans are used to detect grain edges.
The distance in pixels between peaks in the derivative of the profile line is
then converted into a distance in $\mu$m.

###  2.7 Support Vector Regression (SVR)

Once fingerprints have been extracted from the ResNet18 VAE, regression
algorithms can be trained to predict morphological metrics. Here, we use
support vector regression (SVR), which is an extension to support vector
machines for continuous variables [29]. Gaussian process regression is a
popular alternative to SVR, but relies on relatively low-dimensional inputs
compared to the dimensionality of the ResNet18 VAE fingerprints. Fingerprints
are randomly split into a training set containing 90% of the fingerprints and
a test set containing 10% for input into the SVR. This is repeated 10 times
with random splits to perform 10-fold cross-validation.

###  2.8 $t$-Stochastic Neighbour Embedding ($t$-SNE)

The encoded space learned is high-dimensional (256 dimensions for ResNet18-VAE
architecture, outlined in Section 2.5). This makes it difficult to visualise
the encoded space in its entirety. To aid in visualising this high-dimensional
space, we can perform dimensionality reduction down to 2 or 3 dimensions,
which allows us to plot the encoded space and to assess how the
microstructures in our training set are distributed. Here, we consider
$t$-stochastic neighbour embedding ($t$-SNE) [30].

The fingerprints are represented as a similarity matrix, $S$, where entries,
$S_{i,j}$, denote the probability that $z_{j}$ is a nearest neighbour of
$z_{i}$. The aim is then to minimise the distance between $S_{i,j}$ and
$S_{j,i}$ [31]. Standard SNE utilises a Gaussian distribution to determine
similarity between fingerprints, whereas $t$-SNE relies on the Student’s $t$
distribution, hence its name. Due to the longer tail of the $t$ distribution,
relative to a Gaussian distribution, the use of the $t$ distribution results
in more separated embeddings and helps alleviate the crowding issue often
encountered with SNE [32]. Principal component analysis (PCA) is used to
initialise the $t$-SNE.

###  2.9 Traversing the Encoded Space

One way in which the encoded space can be traversed is through construction of
a specific path. Here, we consider a linear path. Two images, $x_{1}$ and
$x_{2}$, are randomly selected and $e(x_{1}),e(x_{2})$ are computed. A linear
path is then constructed from $e(x_{1})$ to $e(x_{2})$ in the encoded space,
according to the following equation.

$e_{n}(x_{1},x_{2})=e(x_{1})+\frac{n+1}{N}(e(x_{2})),\hskip
10.00002ptn=1,\dots,N,$

where $e_{n}(x_{1},x_{2})$ denotes each fingerprint along the path and N is
the number of steps along the path. Fingerprints along the path are then
supplied to the decoder to output reconstructions.

Another method for generating potentially valid encoded representations of
microstructure is to form a random cloud, centered at a known valid encoded
representation, $e(x)$, for some image $x$. This can be achieved by sampling
random Gaussian noise from $\mathcal{N}(0,1)$ and summing with $e(x)$, i.e.,

$e_{n}(x)=e(x)+\gamma\mathcal{N}(\mu,\sigma),\hskip 10.00002ptn=1,\dots,N,$

where $e_{n}(x)$ is a random neighbour of $e(x)$, $\gamma$ is a tunable scale
factor that controls the noise level, $\mu,\sigma\in\mathbb{R}^{D}$ denote the
element-wise mean and standard deviation for each dimension in the encoded
space and $N$ is the number of neighbouring fingerprints to be output.

##  3 Results

###  3.1 Reconstruction Accuracy

The ResNet18 VAE was trained under two separate regimes. The first included
exclusively either lamellar or bimodal microstructures and the second included
the full LFTi64BL dataset, containing both lamellar and bimodal
microstructures. Figure 8 shows some example reconstructions from the lamellar
dataset after 1000 epochs. When training is restricted to the lamellar
dataset, grain boundaries are accurately identified and reconstructed, but
there is some smoothing evident in the reconstructions and the interlamellar
$\upbeta$ appears coarsened. This becomes clearer upon inspection of the
morphological metric distributions shown in Figure 9. The
$\upalpha_{\text{p}}$ volume fraction is slightly reduced in the
reconstructions, with some anomalies between 30% and 50%. The mean lamellae
direction, relative to the bottom edge of each image patch, is consistent.
This confirms that the $\upalpha_{\text{p}}$ grains are oriented correctly in
the reconstructions, but there is an increase in the mean lamellae width,
which is likely due to the smoothing effect removing smaller grains. There is
also an increase in the mean lamellae aspect ratio, owed to the coarsening of
the interlamellar $\upbeta$ in the reconstructions which results in
$\upalpha_{\text{p}}$ grains appearing more elongated.

Figure 8: Representative examples of (a) $256\times 256$ patches from the
lamellar microstructures and (b) their corresponding reconstructions from the
ResNet18 VAE architecture, after 1000 epochs.

Figure 9: Morphological metric distributions across the original image patches
and their corresponding reconstructions for the lamellar dataset.

Training the VAE on exclusively bimodal microstructures results in a similar
smoothing effect as with the lamellar images, although this appears more
pronounced in the bimodal case due to the nature of the fine scale features
present in this case. Figure 10 shows some example reconstructions for the
bimodal micrographs. All $\upalpha_{\text{s}}$ laths retained in the prior
$\upbeta$ grains are completely absent from the reconstructions and the
outputs are effectively a mask for the $\upalpha_{\text{p}}$ grains. This is
confirmed to be a feature of the entire set of reconstructions in Figure 11,
which shows a drastic increase in both the $\upalpha_{\text{p}}$ volume
fraction and mean $\upalpha_{\text{p}}$ grain size.

Figure 10: Representative examples of (a) $256\times 256$ patches from the
bimodal microstructures and (b) their corresponding reconstructions from the
ResNet18 VAE architecture, after 1000 epochs.

Figure 11: Morphological metric distributions across the original image
patches and their corresponding reconstructions for the bimodal dataset.

###  3.2 Traversing the Encoded Space

The encoded space was explored with the methods discussed in Section 2.9.
Figure 12 shows reconstructions along a linear path between a pair of
fingerprints from the training set. A linear path between the fingerprints was
constructed and 8 equispaced fingerprints were determined along the path.
These fingerprints were supplied to the trained decoder to generate the images
in Figure 12.

Figure 12: Reconstructions along linear path through the encoded space learned
during training on exclusively lamellar image data from LFTi64BL.

The fingerprints along this path show a smooth transition along the linear
path, which is a direct result of the continuity of the latent space learned
by the VAE. These microstructures are synthetic and are not included in the
training set, although they do have the same characteristics as the lamellar
microstructures in the training set.

Figure 13 shows the output from the same method applied to the full LFTi64BL
dataset, with the linear path defined between a lamellar and a bimodal
microstructure. The same smooth transition between the input microstructures
is observed, however, intermediate microstructures along the path stray
considerably away from the training set, particularly towards the centre of
the path. This is to be expected, but confirms that, for VAEs to generate
convincing artificial micrographs, it is crucial that the training set be
cohesive and not contain drastic variations in microstructure. Otherwise, the
latent space constructed is likely to contain microstructures that are not
representative of the training set.

Figure 13: Reconstructions along linear path through the encoded space learned
during training on the full LFTi64BL dataset.

To visualise neighbouring microstructures localised around an individual
fingerprint, a sample microstructure was randomly selected and supplied to the
trained encoder. The fingerprint obtained was perturbed with a small amount of
Gaussian noise, as described in Section 2.9 with $\gamma=0.2$, and provided to
the trained decoder to generate a synthetic microstructure. This process was
repeated and Figure 14 shows 10 realisations of microstructures produced.

Figure 14: Synthetic microstructures generated from the addition of unit
Gaussian noise, with a scale factor of $\gamma=0.2$, to a known valid encoded
representation, from which morphology is inherited.

Each artificial microstructure possesses similar features to the input
micrographs in terms of grain morphology and direction of lamellae with
respect to the bottom edge of the image. This shows that local fingerprints
within the latent space are likely to possess similar microstructural
features. Figure 15 shows example micrographs constructed with increasing
values of $\gamma$. As $\gamma$ increases, images generated are no longer
representative of the training dataset. PCA was trained on the fingerprints
constructed by the VAE and then used to transform the noise-perturbed
fingerprints alongside the original fingerprints (see Figure 16). This shows
the noise-perturbed fingerprints emanating from a single point, which
corresponds to the fingerprint to which the noise was applied, and we can see
that noise added with a scale factor of $\gamma>0.5$ are all completely
outside the encoded space learned by the VAE.

(a) Original

(b) $\gamma=0.05$

(c) $\gamma=0.1$

(d) $\gamma=0.2$

(e) $\gamma=0.5$

(f) $\gamma=1$

Figure 15: Synthetic micrograph examples generated via the Gaussian cloud
method with various $\gamma$ values, compared with original micrograph. Figure
16: PCA applied to fingerprints leanred by the VAE alongside fingerprints
constructed with the Gaussian cloud method. Original fingerprints are shown in
grey and the colour map denotes the $\gamma$ values applied to generate the
artifical fingerprints.

###  3.3 Metric Distributions Across the Encoded Space

Dimensionality reduction, in the form of $t$-SNE (described in Section 2.8),
was also applied to the original 256-dimensional fingerprints to reduce them
down to 2-dimensional vectors. This enables fingerprints to be plotted in a
2-dimensional scatter plot to visualise the entire latent space and
distribution of metrics across the space. Figure 17 shows the full LFTi64BL
dataset reduced to 2-dimensions. The colour map in this figure that
illustrates the microstructure classification, with 0 denoting a bimodal
microstructure and 1 denoting lamellar. There is a strong clustering of each
class, with only a small overlap between them, allowing fingerprints from each
class to be easily separated. Training an SVM with 90% of the fingerprints
allocated for training and 10% reserved for testing yields a mean
classification accuracy of 99.9% $\pm$ 0.001 after 10-fold cross-validation.

Figure 17: $t$-SNE with 2 components applied to encoded representations, where
the colour map denotes classifications, with bimodal microstructures in red
and lamellar microstructures in blue.

The same methodology for dimensionality reduction was applied to the lamellar
and bimodal fingerprints separately to map distributions of various
microstructural features. Figures 18(a) and 18(b) shows greyscale intensity
from the lamellar micrographs, before and after normalisation.

Prior to normalisation, there is a strong clustering between two groups of
images, heavily influenced by illumination during image capture. After
normalisation, greyscale intensity is more closely linked to volume fraction
and is normally distributed across the dataset. Dimensionality reduction then
shows a smooth gradient of greyscale intensity across the latent space. The
plots with volume fraction yield similar results, although less pronounced.
Figures 18(c) and 18(d) shows plots of the encoded spaces for bimodal and
lamellar microstructures, trained separately.

Finally, we look at directionality (discussed in Section 2.6.3) for lamellar
microstructures and grain size (Section 2.6.4) for bimodal microstructures as
morphological features of interest. Figure 18(e) shows the $t$-SNE plot with a
colour map corresponding to lamellae direction.

This appears to show random scatter across the encoded space, in contrast to
the Gaussian cloud for image generation that seems to reconstruct images with
similar direction, when perturbing fingerprints with a small amount of
Gaussian noise. Nearest neighbours are shown to have similar directionality,
but this is contained within small regions of the latent space and not
universal. This random distribution of directionality may be useful in
practice, though, as this implies that when encoding the microstructural
information in such a manner, directionality can be effectively ignored and
sample orientation when imaging would not be of any concern.

Figure 18(f) shows the $t$-SNE plot for the bimodal fingerprints with a colour
map denoting grain size. Here, we see a gradient of grain size across the
latent space, suggesting that nearest neighbours will share similar
morphologies.

(a) Lamellar greyscale prior to normalisation

(b) Lamellar greyscale after normalisation

(c) Lamellar $\upalpha_{\text{p}}$ volume fraction

(d) Bimodal $\upalpha_{\text{p}}$ volume fraction

(e) Lamellae directionality.

(f) Bimodal mean $\upalpha_{\text{p}}$ grain size.

Figure 18: $t$-SNE with 2 components applied to encoded representations of the
LFTi64BL dataset, with colour maps denoting various morphological metrics.

In each $t$-SNE reduction, there is a small group of fingerprints that are
separated away from the main cluster. There are no discernible differences
between the microstructures that result in these fingerprints and those that
form the main cluster. Despite these microstructures being included in the
training set, their reconstructions after training are effectively just noise.
This is a caveat of the VAE training process and not the dimensionality
reduction. Removing these outliers from the set of fingerprints and rerunning
the $t$-SNE with the remaining fingerprints results in retention of the main
cluster and a complete absence of the secondary cluster. Removing the images
that reside in the secondary cluster from the training set and retraining the
VAE from scratch results in a new cluster forming, with a different
morphology, that contains a different set of images. Nonetheless, we are able
to explore the encoded space of valid reconstructions, using the $t$-SNE as a
guide to identify suitably encoded microstructures.

###  3.4 Metric Predictions from Encoded Representations

If an approach to microstructural fingerprinting is successfully encoding the
essence of the microstructure, then we might reasonably expect to be able to
recover key features of the microstructure, such as grain size, from the
encoded fingerprint. A strong correlation between the fingerprint and a given
property also opens the possibility of reversing the inferential process and
asking what fingerprints (or which regions of latent space) would correspond
to a microstructure exhibiting a given property of interest (with the
possibility with a VAE of then reconstructing an image of the corresponding
microstructure).

The predictability of morphological metrics across the encoded space was
validated with Support Vector Regression (SVR), described in Section 2.7. The
fingerprints output from the RestNet18 VAE are randomly split into 90%
training data and 10% test data. This is repeated 10 times to perform 10-fold
cross-validation. The SVR is then trained on the training set and outputs for
the test set are compared with the true values to measure the prediction
accuracy. Figure 19 shows scatter plots for the SVR predictions across the 10
train-test splits combined, against the true values, with the line $y=x$
plotted in red to illustrate deviation from exact predictions. Table 1 then
shows percentage error and standard deviation of the predictions relative to
the true values.

(a) Lamellar greyscale intensity

(b) Bimodal greyscale intensity

(c) Lamellar volume fraction $\upalpha_{\text{p}}$

(d) Bimodal volume fraction $\upalpha_{\text{p}}$

(e) Lamellae direction

(f) Bimodal grain size

Figure 19: Scatter plots to illustrate combined 10-fold cross validation SVR predictions of various morphological metrics with the line $y=x$ added to highlight where true predictions should lie. Class | Metric | Mean Percentage Error
---|---|---
Lamellar | Greyscale Intensity | 0.60 % $\pm$ 0.01
Volume Fraction $\upalpha_{\text{p}}$ (%) | 1.83 % $\pm$ 0.03
Directionality (∘) | 148 % $\pm$ 29
Bimodal | Greyscale Intensity | 0.75 % $\pm$ 0.01
Volume Fraction $\upalpha_{\text{p}}$ (%) | 2.50 % $\pm$ 0.03
Grain Size ($\mu$m) | 7.91 % $\pm$ 0.12
Table 1: Quantitative analysis of error for SVR predictions.

##  4 Discussion

Reconstructions from the proposed ResNet18 VAE were shown to accurately
identify grain boundaries, albeit with some smoothing. This information loss
is acceptable for correlating with the morphological features discussed, as
these features are still captured. However, prediction of properties that rely
upon small scale features, such as fine $\upalpha_{\text{s}}$ laths, may
become hindered. This could be alleviated with higher resolution images and
patches with a higher magnification.

The encoded space learned by the VAE is continuous and contains valid
microstructural fingerprints that are not included in the training set. Two
methods were used to explore the encoded space as a tool for generating
artificial microstructures. The linear path method illustrates the continuity
of the space and shows a smooth transition between microstructures, when
linearly interpolating between two known fingerprints. With a training set
containing only lamellar microstructures, reconstructions along this path are
completely artificial but still resemble lamellar microstructures. Once the
training set is expanded to include both bimodal and lamellar microstructures,
the linear interpolation between a bimodal fingerprint and a lamellar
fingerprint yields a blend between the two morphologies. These can still be
perceived as valid microstructures, but it is important to note that these are
not representative of any data in the training set and VAEs trained in such a
way should be interpolated carefully.

Gaussian noise was also used to explore the encoded space more locally around
an individual fingerprint. Iteratively perturbing a learned fingerprint and
supplying the output to the trained decoder yields another set of artificially
generated microstructures. In this case, the decoder generates microstructures
with similar features to the input image, including the lamellae thickness and
direction relative to the bottom of the image, provided that a suitably small
amount of noise is applied. This behaviour is only observed locally for some
features, such as lamellae direction, which is randomly distributed across the
encoded space. This may be due to the fact that the VAE was trained on patches
with random rotations applied, making the VAE rotationally invariant in this
case. This removes any bias towards sample orientation during image capture.
Adding large amounts of noise would result in variations of such features.

Dimensionality reduction via $t$-SNE was performed to reduce the fingerprints
to 2 dimensions. This enables them to be plotted in a 2-dimensional scatter
plot to visualise the distribution of morphological features across the
encoded space. The fingerprint position in the encoded space appears to
dependent significantly on the volume fraction of $\upalpha_{\text{p}}$,
relative to other morphological metrics. This was confirmed by SVR
predictions, which was trained to predict volume fraction with an average
percentage error of 1.83 % $\pm$ 0.03 for the lamellar dataset and 2.50 %
$\pm$ 0.03 for the bimodal dataset. Grain size was also highly correlated with
the encoded representations and SVR was able to predict grain size with
reasonable accuracy, giving an average percentage error of 7.91 % $\pm$ 0.12.
Directionality of lamellae appeared to be distributed randomly across the
encoded space from the $t$-SNE plots and this was also confirmed with SVR,
which was unable to learn any trends in the data and repeatedly predicted the
mean direction for each fingerprint, resulting in an average percentage error
of 148 % $\pm$ 29\. In practice, the random nature of the distribution of
lamellae direction could be a useful feature, as this implies that sample
orientation under the microscope can be effectively ignored when capturing a
dataset for VAE without any bias being introduced.

##  5 Conclusions

The latent space of microstructural fingerprints generated from variational
autoencoders (VAEs) has been explored to further our understanding of how VAEs
encode feature information from microstructural image data.

  * •

We show that a variational autoencoder (VAE) architecture based on ResNet18 is
able to produce accurate reconstructions of microstructures (up to some
smoothing).

  * •

Interpolation of fingerprints along a linear path in latent space and random
perturbations about fingerprints of input microstructures resulted in the
generation of plausible synthetic microstructures, demonstrating the
suitability of the VAE for smooth representation of microstructure.

  * •

Fingerprints constructed by the trained VAE are shown to correlate with
morphological features, including $\upalpha{\text{p}}$ volume fraction and
grain size, using support vector regression (SVR) with 10-fold cross-
validation.

  * •

Principal component analysis (PCA) is shown to provide useful insight into the
amount of noise that can be added to known fingerprints before the perturbed
fingerprints become disconnected from the learned latent space.

Our study based on a set of micrographs of titanium alloy microstructure
demonstrates that a VAE can encode material microstructures to produce
fingerprints that exhibit several of the key features required in a general
purpose fingerprint. Important next steps will be to test this approach on a
broader class of microstructures and to apply it to datasets exhibiting
variation in both microstructure and measured properties, to allow exploration
of the suitability of VAEs for predicting microstructure-process-property
relationships in an explainable framework. We also note that generative
adversarial networks (GANs) are found to generate high-fidelity synthetic
microstructures, with strong statistical similarity to the training data, even
in a limited data regime [33]. It would be interesting to consider synthetic
images output from a GAN as inputs to the VAE to determine how the generated
images are distributed through the encoded space and assess their potential
for bolstering morphology prediction.

## Acknowledgements

This research was supported through funding from the Scheme for Promotion of
Academic and Research Collaboration (SPARC) grant MHRD-18-0015. MDW was
supported by the University of Manchester and Rolls-Royce plc under the
Engineering and Physical Sciences Research Council (EPSRC) grant EP/S022635/1
and Science Foundation Ireland (SFI) grant 18/EPSRC-CDT/3584. CPR was
supported by a University Research Fellowship of the Royal Society. PJW was
supported by the Henry Royce Institute for Advanced Materials, funded through
EPSRC grants EP/R00661X/1, EP/S019367/1, EP/P025021/1 and EP/P025498/1.

## References

  * [1] B. Basu, N. H. Gowtham, Y. Xiao, S. R. Kalidindi, K. W. Leong, Biomaterialomics: Data science-driven pathways to develop fourth-generation biomaterials, Acta Biomaterialia (2022).  doi:10.1016/j.actbio.2022.02.027.
  * [2] E. A. Holm, R. Cohn, N. Gao, A. R. Kitahara, T. P. Matson, B. Lei, S. R. Yarasi, Overview: Computer vision and machine learning for microstructural characterization and analysis, Metallurgical and Materials Transactions A 51 (2020) 5985–5999.  doi:10.1007/s11661-020-06008-4.
  * [3] T. L. Burnett, P. J. Withers, Completing the picture through correlative characterization, Nature Materials 18 (2019) 1041–1049.  doi:10.1038/s41563-019-0402-8.
  * [4] B. L. DeCost, E. A. Holm, A computer vision approach for automated analysis and classification of microstructural image data, Computational Materials Science 110 (2015) 126–133.  doi:10.1016/j.commatsci.2015.08.011.
  * [5] M. I. Latypov, S. R. Kalidindi, Data-driven reduced order models for effective yield strength and partitioning of strain in multiphase materials, Journal of Computational Physics 346 (2017) 242–261.  doi:10.1016/j.jcp.2017.06.013.
  * [6] Y. C. Yabansu, P. Steinmetz, J. Hötzer, S. R. Kalidindi, B. Nestler, Extraction of reduced-order process-structure linkages from phase-field simulations, Acta Materialia 124 (2017) 182–194.  doi:10.1016/j.actamat.2016.10.071.
  * [7] E. Popova, T. M. Rodgers, X. Gong, A. Cecen, J. D. Madison, S. R. Kalidindi, Process-structure linkages using a data science approach: Application to simulated additive manufacturing data, Integrating Materials and Manufacturing Innovation 6 (1) (3 2017).  doi:10.1007/s40192-017-0088-1.
  * [8] K. Rajan, Materials informatics: The materials “gene” and big data, Annual Review of Materials Research 45 (1) (2015) 153–169.  doi:10.1146/annurev-matsci-070214-021132.
  * [9] D. L. McDowell, R. A. LeSar, The need for microstructure informatics in process–structure–property relations, MRS Bulletin 41 (8) (2016) 587–593.  doi:10.1557/mrs.2016.163.
  * [10] M. D. White, A. Tarakanov, P. J. Withers, C. P. Race, K. J. H. Law, Digital fingerprinting of microstructures, Computational Materials Science 218 (2023) 111985\.  doi:10.1016/j.commatsci.2022.111985.
  * [11] D. T. Fullwood, S. R. Niezgoda, B. L. Adams, S. R. Kalidindi, Microstructure sensitive design for performance optimization, Progress in Materials Science 55 (6) (2010) 477–562.  doi:10.1016/j.pmatsci.2009.08.002.
  * [12] D. P. Kingma, M. Welling, Auto-Encoding Variational Bayes, arXiv (May 2014).  doi:10.48550/arXiv.1312.6114.
  * [13] Z. Liu, P. Luo, X. Wang, X. Tang, Deep learning face attributes in the wild, in: Proceedings of International Conference on Computer Vision (ICCV), 2015.
  * [14] A. Sardeshmukh, S. Reddy, B. P. Gautham, P. Bhattacharyya, TextureVAE: Learning Interpretable Representations of Material Microstructures Using Variational Autoencoders, in: Proceedings of the AAAI 2021 Spring Symposium on Combining Artificial Intelligence and Machine Learning with Physical Sciences, Stanford, CA, USA, March 22nd - to - 24th, 2021, Vol. 2964 of CEUR Workshop Proceedings, CEUR-WS.org, 2021.
  * [15] R. Cang, H. Li, H. Yao, Y. Jiao, Y. Ren, Improving direct physical properties prediction of heterogeneous materials from imaging data via convolutional neural network and a morphology-aware generative model, Computational Materials Science 150 (2018) 212–221.  doi:10.1016/j.commatsci.2018.03.074.
  * [16] Y. Kim, H. K. Park, J. Jung, P. Asghari-Rad, S. Lee, J. Y. Kim, H. G. Jung, H. S. Kim, Exploration of optimal microstructure and mechanical properties in continuous microstructure space using a variational autoencoder, Materials & Design 202 (2021) 109544.  doi:10.1016/j.matdes.2021.109544.
  * [17] H. S. Stein, D. Guevarra, P. F. Newhouse, E. Soedarmadji, J. M. Gregoire, Machine learning of optical properties of materials – predicting spectra from images and images from spectra, Chemical Science 10 (2019) 47–55.  doi:10.1039/C8SC03077D.
  * [18] Y. Pathak, K. S. Juneja, G. Varma, M. Ehara, U. D. Priyakumar, Deep learning enabled inorganic material generator, Physical Chemistry Chemical Physics 22 (2020) 26935–26943.  doi:10.1039/D0CP03508D.
  * [19] S. Sundar, V. Sundararaghavan, Database development and exploration of microstructure versus process relationships using variational autoencoders, arXiv (2020).  doi:10.48550/ARXIV.2001.09171.
  * [20] K. He, X. Zhang, S. Ren, J. Sun, Deep residual learning for image recognition, arXiv (2015).  doi:10.48550/ARXIV.1512.03385.
  * [21] M. White, C. Daniel, J. Quinta da Fonseca, X. Zeng, N. Byers, B. Karnasiewicz, Lamellar and bi-modal ti-64 microstructure images (Nov. 2021).  doi:10.5281/zenodo.5714384.
URL https://zenodo.org/records/5714384

  * [22] D. P. Kingma, J. Ba, Adam: A Method for Stochastic Optimization, arXiv (2014).  doi:10.48550/ARXIV.1412.6980.
  * [23] D. P. Kingma, M. Welling, An introduction to variational autoencoders, Foundations and Trends® in Machine Learning 12 (4) (2019) 307–392.  doi:10.1561/2200000056.
  * [24] R. Yu, A tutorial on VAEs: From Bayes’ rule to lossless compression, arXiv (2020).  doi:10.48550/arXiv.2006.10273.
  * [25] S. Björk, J. N. Myhre, T. H. Johansen, Simpler is better: spectral regularization and up-sampling techniques for variational autoencoders, arXiv (2022).  doi:10.48550/ARXIV.2201.07544.
  * [26] J. Deng, W. Dong, R. Socher, L.-J. Li, K. Li, L. Fei-Fei, Imagenet: A large-scale hierarchical image database, in: 2009 IEEE conference on computer vision and pattern recognition, Ieee, 2009, pp. 248–255.  doi:10.1109/CVPR.2009.5206848.
  * [27] S. L. Bangare, A. Dubal, P. S. Bangare, S. T. Patil, Reviewing Otsu’s method for image thresholding, International Journal of Applied Engineering Research 10 (9) (2015) 21777–21783.  doi:10.37622/IJAER/10.9.2015.21777-21783.
  * [28] E04 Committee, Test Methods for Determining Average Grain Size Using Semiautomatic and Automatic Image Analysis, Tech. rep., ASTM International.  doi:10.1520/E1382-97R15.
  * [29] M. Awad, R. Khanna, Support Vector Regression, In: Efficient Learning Machines: Theories, Concepts, and Applications for Engineers and System Designers, Apress, Berkeley, CA, 2015, pp. 67–80.  doi:10.1007/978-1-4302-5990-9_4.
  * [30] R. Silva, P. Melo-Pinto, t-SNE: A study on reducing the dimensionality of hyperspectral data for the regression problem of estimating oenological parameters, Artificial Intelligence in Agriculture 7 (2023) 58–68.  doi:10.1016/j.aiia.2023.02.003.
  * [31] G. E. Hinton, S. Roweis, Stochastic neighbor embedding, in: S. Becker, S. Thrun, K. Obermayer (Eds.), Advances in Neural Information Processing Systems, Vol. 15, MIT Press, 2002.
  * [32] L. van der Maaten, G. Hinton, Visualizing data using t-sne, Journal of Machine Learning Research 9 (86) (2008) 2579–2605.
  * [33] N. H. Gowtham, T. J. Jeyapriya, C. Bhattacharya, B. Basu, A deep adversarial approach for the generation of synthetic titanium alloy microstructures with limited training data, Computational Materials Science 230 (2023) 112512.  doi:10.1016/j.commatsci.2023.112512.


# Occurrence of a Turning Point in the Dynamic Solution for a
Reissner–Nordström Black Hole

Alexander Shatskiy<EMAIL_ADDRESS>
Abstract:

A self-consistent exact solution for a Reissner–Nordström black-and-white hole
formed as a result of accretion has been considered. Prior to the formation of
a black-and-white hole, there is a bulk charged sphere at the center of the
system. The occurrence of a turning (bounce) point in the newly formed black-
and-white hole is investigated. The occurrence of a turning point was
investigated by an example of solution for the accretion of a neutral
spherical dust shell. The model equations are written with allowance for the
cosmological $\Lambda$-term. Within the model under consideration, both the
black-and-white hole and its turning point are formed in the already existing
Universe; therefore, the black-and-white hole in this model is not
<<eternal>>.

##  I NTRODUCTION

The existence of black holes (BHs), as well as the mechanisms of their
formation, have long been considered as investigated. The mathematical
solutions for real (rotating) Kerr BHs and charged Reissner–Nordström BHs show
the existence of a turning point. The turning (or bounce) point is considered
to be a point at which the free-fall trajectory of a particle changes from the
direction along which the particle radial coordinate decreases to a direction
along which this coordinate increases. Here, of interest are specifically the
turning points located under both black hole horizons; i.e., under the Cauchy
horizon: in the internal R-region of space.

The existence of internal turning point (turning point below) in the steady-
state solution for BH is related to the existence of black-and-white hole
rather than a conventional BH. However, the mechanism of turning-point
occurrence in the self-consistent solution of collapsing matter has not yet
been studied (see, e.g., Ch. Misner, K. Thorne, and J. Wheeler (1977), §34.6).
Therefore, the models describing the occurrence of black-and-white holes still
contain more questions than answers. It is agreed that black-and-white holes
can be only <<eternal>> BHs, i.e., the ones formed jointly with the creation
of the Universe. Therefore, one of the purposes of this study was to justify
the model in which a black-and-white hole itself and its turning point can be
formed in the already existing Universe. Note that the white hole is formed in
the new expanding universe rather than in the same Universe in which the black
hole is formed.

We will consider a self-consistent exact solution for a spherically symmetric
Reissner–Nordström black-and-white hole during its formation. Prior to the
formation of a black-and-white hole, there is a bulk charged sphere at the
center of the system. The occurrence of a turning point was investigated by an
example of exact solution for the accretion of a thin neutral spherical dust
shell. The model equations will be written with allowance for the cosmological
$\Lambda$-term.

Previously many researchers have studied the dynamics of spherical dust shells
in a centrally symmetric electric field (see, e.g., D.A. Tretyakova, A.A.
Shatskiy, I.D. Novikov, S. Alexeyev (2012); A.A. Shatskiy and Novikov (2010);
V.A.Berezin and V.I.Dokuchaev (2014)). In particular, since the topology of
Reissner–Nordström BH was investigated by plotting Carter–Penrose diagrams in
V.A.Berezin and V.I.Dokuchaev (2014), we omit these diagrams here and refer
the reader to V.A.Berezin and V.I.Dokuchaev (2014) for the information on this
issue. The topology of real rotating BHs was investigated by the author in
Shatskiy (2020).

##  II DESCRIPTION OF THE MODEL

Let there be no BH and horizons at the initial instant. The matter
distribution at the initial instant is described by the following model. There
is a sphere at the center of the system, which has an electrical charge $q$,
radius $r_{q}$ , and mass $m_{q}$. Let the mass $m_{q}$ be insufficient for
the formation of a BH and gravitational collapse of the sphere.

Let also a spherical dust shell with a mass ${m_{dust}}$ be at rest on some
radius $R_{1}$ at the initial instant. During further evolution this dust
shell begins to fall freely towards the center of the system. The mass
${m_{dust}}$ turns out to be sufficient for the further formation of a
Reissner–Nordström BH before reaching the shell radius $r_{q}$. Let the total
mass of the system be M. The BH horizon radius $r_{h}^{+}$ and the Cauchy
horizon radius ${r_{h}^{-}}$, corresponding to the charge $q$ and mass ${M}$,
can be written in the Reissner–Nordström coordinates as111In this study, if
not specified separately, we use the theoretical system of units, in which the
speed of light and gravitational constant are, respectively, ${c=1}$ and
${G=1}$.:

$r_{h}^{\pm}=M\pm\sqrt{M^{2}-q^{2}}$ (1)

For convenience we introduce the charge-to-mass ratio so as to satisfy the
relations

$r_{h}^{-}:=\beta^{2}r_{h}^{+}\,,\quad\beta=const\leq 1\,.$ (2)

According to (1), the necessary conditions for this situation are

$q=\kappa M\,,\quad\kappa:=\frac{2\beta}{1+\beta^{2}}\leq 1\,.$ (3)

We introduce also a factor $\gamma$ to link the charge and mass with the
radius $r_{q}$:

$r_{q}=\gamma r_{h}^{-}=\gamma\beta^{2}r_{h}^{+}\,,\quad\gamma=const<1\,.$ (4)

Obviously, the factor $\gamma$ should to smaller than unity, because both BH
horizons must arise before the possible collision of the falling dust shell
with the charged sphere; i.e., the condition ${r_{q}<r_{h}^{-}\leq r_{h}^{+}}$
must be satisfied.

It follows from the above expressions that

$r_{q}=\beta\gamma q=\beta\gamma\kappa M\,,\quad r_{h}^{-}=\beta q\,,\quad
r_{h}^{+}=q/\beta\,.$ (5)

Below we will prove the possibility of existence of turning point ${r_{t}}$,
which arises in solution in the gap between the radii $r_{q}$ and
${r_{h}^{-}}$.

To proceed further, we should solve the problem with limitation on the sphere
size $r_{q}$. Here, the main limiting factor is the electric field strength
$E$. The maximum strength, ${E=q/r^{2}}$, should be below the vacuum
breakdown: ${\sim 10^{4}\,\makebox{V/sm}}$. Since the maximum strength is
achieved near the sphere surface (i.e., at $r_{q}$), the value
${E_{q}=q/r_{q}^{2}}$ should be smaller than ${10^{4}\,\makebox{V/sm}}$.

Let it be

$E_{q}=10^{3}\,\makebox{V/sm}=3.(3)\makebox{ <<CGSE>> units}$ (6)

It follows from (5) that

$E_{q}:=\frac{q}{r_{q}^{2}}=\frac{\beta^{-2}\gamma^{-2}}{q}$ (7)

Based on (7), we can express the electrical charge $q$ as

$q=\frac{\beta^{-2}\gamma^{-2}}{E_{q}}$ (8)

We assume that ${\beta^{2}=0.5}$ and ${\gamma=0.5}$; then
${\kappa=\sqrt{8}/3}$. Having rewritten expression (8) in CGSE units, we will
have (to the end of this section) the following expression for the charge in
these units:

$q=\frac{\beta^{-2}\gamma^{-2}c^{4}}{GE_{q}}\approx 2.915\cdot
10^{49}\makebox{<<CGSE>> units}\approx 6.073\cdot 10^{58}\makebox{<<e>>}$ (9)

Correspondingly, for the total BH mass $M$ we have

$M=\frac{q}{\kappa\sqrt{G}}\approx 1.197\cdot 10^{53}\makebox{g}\approx
6.02\cdot 10^{20}\,M_{\odot}$ (10)

Based on (5) and (9), we find the charged sphere radius ${r_{q}}$ to be

$r_{q}=\frac{q\beta\gamma\sqrt{G}}{c^{2}}=\frac{c^{2}}{\beta\gamma\sqrt{G}E_{q}}\approx\,2.957\cdot
10^{24}\,\makebox{sm}\approx 958\,000\makebox{pc}\sim 1\makebox{Mpc}$ (11)

Let us estimate the density of the excess (or deficit) elementary charges
comprising the charge $q$ of a sphere of radius $r_{q}$. According to (9) and
(11), this surface charge density will be as low as ${\sigma_{q}\sim 5.46\cdot
10^{8}\makebox{e/sm}^{2}\approx 8.74\cdot 10^{-7}\makebox{C/m}^{2}}$.

On the assumption that the charged-sphere matter has a density on the order of
${1\makebox{g/sm}^{3}}$, at ${q^{2}=Gm_{q}^{2}}$, the thickness $d_{q}$ of
this charged spherical shell, according to (9) and (11), turns out to be about
${70}$ m. As compared with the radius $r_{q}$ from (11), this is an ultrathin
film; however, in comparison with the thickness of the layer in which excess
(or deficit) electrons with a charge $q$ are concentrated (i.e., as compared
with atomic sizes), this is an immense value.

One can also calculate the surface tension $N$ formed by the electrical force
of repulsion for a surface charge $q$ on a sphere of radius $r_{q}$:

$N=\frac{q^{2}}{4\pi r_{q}^{3}}=\frac{\beta^{-2}\gamma^{-2}c^{4}}{4\pi
Gr_{q}}\approx 2.616\cdot 10^{24}\makebox{din/sm}$ (12)

With the gravitational force of attraction neglected, we have to divide the
tension N by the specific strength of the film material in order to calculate
the minimum thickness ${d_{min}}$ for a film retaining charge by its surface
tension. Let this material be steel with a specific strength
${P_{s}=10^{10}\makebox{d/sm}^{2}}$. Then we find from (12) that

$d_{min}=\frac{N}{P_{s}}=\frac{q^{2}}{4\pi P_{s}r_{q}^{3}}\approx 2.616\cdot
10^{14}\makebox{sm}$ (13)

However, the mass of this steel will then be about 13 orders of magnitude
larger than the total mass $M$, which is impossible. Thus, equilibrium cannot
be reached by compensating the electrical force of repulsion with the elastic
force of material.

Note that the surface tension force is almost tangential to the small element
of the charged-sphere surface. Therefore, to compensate for the radial
electrical force, the surface tension force should greatly exceed the
electrical force; the larger the radius ${r_{q}}$, the larger this excess must
be.

A calculation of the pressure $P$ exerted by the electrical force of repulsion
on the charged-sphere surface yields

$P=\frac{q^{2}}{4\pi r_{q}^{4}}\approx 0.88\,\makebox{d/sm}^{2}$ (14)

At such large radii this is a very low pressure; practically all materials can
withstand it. Such a large difference between the surface tension $N$ and
pressure $P$ is due to the fact that ${N\propto 1/r_{q}^{3}}$, while
${P\propto 1/r_{q}^{4}}$; hence, ${P/N\propto 1/r_{q}}$. Therefore, if there
is another radial force acting in the opposite direction (gravity force), this
force can balance theoretically the electrical force (see below).

##  III POSSIBILITY OF EQUILIBRIUM ON A CHARGED SPHERE AT THE INITIAL MOMENT

Let us consider the possibility of establishing equilibrium for a charged
sphere of mass $m_{q}$, radius $r_{q}$, and charge $q$. This equilibrium will
be maintained due to the compensation of the electrostatic repulsion of
charges on the sphere by the sphere self-gravity. This equilibrium is possible
only beyond the gravitational radius of the sphere. As will be seen from
further calculations, the BH horizon may not form at all at these masses and
charges.

It is known from the general-physics course that the electric field potential
$A_{t}$ must be constant inside a charged sphere and decrease inversely
proportionally to radius beyond the sphere. Therefore, the electric field
should be zero within the sphere and equal to ${q/r^{2}}$ beyond it. However,
we will be interested in the electric field on the sphere itself. Obviously,
it is equal to the arithmetic mean of the field limits within the sphere and
beyond it when the radii tend to $r_{q}$; i.e., the electric field on the
sphere is ${q/(2r_{q}^{2})}$.

The same holds true for the gravitational charge (mass); however, we will
assume that the mass is uniformly distributed over a layer of finite thickness
$d$ (the physical thickness of the sphere film). It will be shown below that
the thickness $d$ should be much smaller than the sphere radius $r_{q}$;
nevertheless, $d$ is a macroscopic size. At the same time, the charge $q$ is
concentrated only on the sphere surface (in a layer whose thickness has atomic
sizes). We assume that the charge is determined by the deficit of electrons on
the sphere; i.e., the positive charge of the nuclei of sphere lattice is not
completely compensated by the elementary charge. The sizes of electrons and
nuclei are much smaller than the sphere thickness $d$. Therefore, the
gravitational field behaves identically at different sides of the sphere
external surface, and the electric field under the external sphere surface
becomes zero, in correspondence with the aforesaid.

We consider a charge element $e$ on the sphere and the corresponding mass
element $\mu$, which obey the relations

$e:=\xi\cdot q\,,\quad\mu:=\xi\cdot m_{q}\,,\quad\xi=const<<1\,.$ (15)

Let us write now the equations of motion of a charged particle in
gravitational and electromagnetic fields (see Landau1 , §90):

$\mu\frac{du^{i}}{ds}=eF^{ij}u_{j}-\mu\Gamma^{i}_{km}u^{k}u^{m}\,,\quad\Gamma^{i}_{km}=\frac{1}{2}g^{in}\left(\frac{\partial
g_{nk}}{\partial x^{m}}+\frac{\partial g_{nm}}{\partial x^{k}}-\frac{\partial
g_{km}}{\partial x^{n}}\right)\,.$ (16)

Here, ${F_{ik}}$ is the electromagnetic field tensor, ${u^{k}}$ is the
4-vector of velocity of a particle with a mass $\mu$ and a charge $e$, and
$g_{in}$ are the components of the metric tensor.

In this section, we choose a static frame of reference and the
Reissner–Nordström metric for it:

$ds^{2}=fdt^{2}-f^{-1}\,dr^{2}\,-\,r^{2}\,(d\theta^{2}+\sin^{2}\theta\,d\varphi^{2})\,,\quad
f(r):=\left(1-\frac{2m}{r}+\frac{q^{2}}{r^{2}}\right)\,.$ (17)

For this metric we have

$u^{k}=\delta^{k}_{t}u^{t}\,,\quad u^{t}=f^{-1/2}\,,\quad
F_{ik}:=\partial_{i}A_{k}-\partial_{k}A_{i}\,,\quad
A_{i}=\delta_{i}^{t}A_{t}\,.$ (18)

Here, ${A_{i}}$ is the 4-vector potential of electromagnetic field;
${A_{t}=q/r}$ and ${A_{t}=q/r_{q}}$ at ${r\geq r_{q}}$ and ${r\leq r_{q}}$,
respectively.

All quantities depend on only the radius, and Eq. (16) can be rewritten in the
form

$\mu\frac{du^{r}}{ds}=-e\,\sqrt{f}\cdot
F_{rt}-\frac{\mu}{2}\cdot\frac{df}{dr}$ (19)

In correspondence with the aforesaid, we have the following expression for the
single tensor component ${F_{ik}}$ on the sphere:
${\left.F_{rt}\right|_{r_{q}}=-0.5q/r_{q}^{2}}$.

In addition, a factor of ${0.5}$ should arise in the term with $q^{2}$ for the
${\partial_{r}f}$ in formula (19), which can be rewritten for the equilibrium
point in the form

$\left.\mu\frac{du^{r}}{ds}\right|_{r=r_{q}}=0.5\cdot\frac{eq\sqrt{f(r_{q})}}{r_{q}^{2}}-\frac{\mu}{2}\,\left(\frac{2m_{q}}{r_{q}^{2}}-0.5\cdot\frac{2q^{2}}{r_{q}^{3}}\right)=0$
(20)

It follows from (20) that this equilibrium point is stable, because the
positive force of repulsion decreases as ${\propto 1/r_{q}^{3}}$, whereas the
negative force of attraction decreases as${\propto 1/r_{q}^{2}}$.

Figure 1:  Dependence of the charged-sphere radius $r_{q}$ on ${q^{2}}$ in
units of mass $m_{q}$ (see (22)).

Using expressions (15), we replace the quantities ${(e,\,\mu)}$ with
${(q,\,m_{q})}$ in Eq. (20) and arrive at

$\displaystyle(4m_{q}^{4}-q^{4})r_{q}^{2}-2q^{2}m_{q}(2m_{q}^{2}-q^{2})r_{q}+q^{4}(m_{q}^{2}-q^{2})=0$
(21) $\displaystyle
r_{q}(q,m_{q})=\frac{m_{q}q^{2}(2m_{q}^{2}-q^{2})+q^{4}\sqrt{2m_{q}^{2}-q^{2}}}{4m_{q}^{4}-q^{4}}$
(22)

As follows from (22), the function ${r_{q}(q)}$ is increasing (see Fig. 1).
According to (5), we are interested in only ${r_{q}<q}$; hence, the limitation
on the charge $q$ takes the form

$m_{q}^{2}<q^{2}<2m_{q}^{2}$ (23)

Therefore, when determining the root of $r_{q}$ in Eq. (21), we chose the
<<$+$>> sign. It follows from formula (23) that a sphere with a mass $m_{q}$
and charge $q$ cannot form a BH; otherwise, expressions (1) for horizons would
be complex values.

##  IV MODEL EQUATIONS

The frame of reference comoving synchronously with matte222 The frame of
reference synchronously comoving with matter is applicable to dust matter,
because pressure is absent in dust (the metric time component is 1).,
described by coordinates $\tau$ and $R$,

$ds^{2}=d\tau^{2}-e^{\lambda(\tau,R)}\,dR^{2}\,-\,r^{2}(\tau,R)\,(d\theta^{2}+\sin^{2}\theta\,d\varphi^{2}),$
(24)

turns out to be more convenient for further analysis. Let the accreting matter
be a gravitating dust with an energy density ${\varepsilon}$, whose
distribution depends on only the radius $r$. We will also take into account
the influence of the cosmological $\Lambda$-term on the matter dynamics in our
model equations. Then the total energy–momentum tensor beyond a sphere of
radius $r_{q}$ can be presented in the form

$T^{n}_{m}=\left(\begin{tabular}[]{c c c c}$+\frac{q^{2}}{8\pi
r^{4}}$&0&0&0\\\ 0&$+\frac{q^{2}}{8\pi r^{4}}$&0&0\\\ 0&0&$-\frac{q^{2}}{8\pi
r^{4}}$&0\\\ 0&0&0&$-\frac{q^{2}}{8\pi r^{4}}$\\\
\end{tabular}\right)+\left(\begin{tabular}[]{c c c c}$\varepsilon$&0&0&0\\\
0&0&0&0\\\ 0&0&0&0\\\ 0&0&0&0\\\ \end{tabular}\right)+\delta^{n}_{m}\Lambda$
(25)

The first term on the right-hand side of (25) corresponds to energy–momentum
tensor of the electromagnetic field of a point charge $q$, the second term
corresponds to the energy–momentum tensor of dust, and the term
${\delta^{n}_{m}\Lambda}$ is the cosmological ${\Lambda}$-term. Inside a
sphere of radius $r_{q}$ the energy–momentum tensor of electromagnetic field
should be zero, because electric field is absent in a charged sphere.
Concerning the energy–momentum tensor of electromagnetic field on the sphere
of radius $r_{q}$, in correspondence with the considerations reported in
Section III, all tensor components on this sphere must be divided by 2,
because the electric field on the sphere is an average of the fields inside
the sphere and beyond it.

Due to the hydrodynamic independence of dust layers relative to each other,
one can integrate the equations of motion for dust similarly to the solution
of Tolman’s problem (see Tolman (1934); Oppenheimer (1939)). In essence, this
is the same Tolman’s problem in a centrally symmetric electric field for
neutral dust.

The Einstein equations corresponding to metric (24) can be written as333
Derivation of Eqs. (26)-(29) can be found, e.g., in Landau1  (§100, problem
5).:

$\displaystyle 8\pi
T_{\tau}^{\tau}=8\pi\varepsilon+q^{2}/r^{4}+8\pi\Lambda=\frac{1}{r^{2}}\left[1+rr_{,_{\tau}}\lambda_{,_{\tau}}+r^{2}_{,_{\tau}}-e^{-\lambda}\left(2rr_{,_{RR}}+r^{2}_{,_{R}}-rr_{,_{R}}\lambda_{,_{R}}\right)\right]$
(26) $\displaystyle 8\pi
T_{R}^{R}=q^{2}/r^{4}+8\pi\Lambda=\frac{1}{r^{2}}\left(1+2rr_{,_{\tau\tau}}+r^{2}_{,_{\tau}}-e^{-\lambda}r^{2}_{,_{R}}\right)$
(27) $\displaystyle 8\pi
T^{R}_{\tau}=0=\left(2r_{,_{R\tau}}-r_{,_{R}}\lambda_{,_{\tau}}\right)e^{-\lambda}/r$
(28) $\displaystyle 8\pi T^{\theta}_{\theta}=8\pi
T^{\varphi}_{\varphi}=-q^{2}/r^{4}+8\pi\Lambda=\frac{r_{,_{\tau\tau}}}{r}+\frac{\lambda_{,_{\tau\tau}}}{2}+\frac{\lambda_{,_{\tau}}^{2}}{4}+\frac{r_{,_{\tau}}\lambda_{,_{\tau}}}{2r}-\left(r_{,_{RR}}-\frac{1}{2}r_{,_{R}}\lambda_{,_{R}}\right)e^{-\lambda}/r$
(29)

Integrating Eq. (28) over time, we obtain

$e^{-\lambda}r^{2}_{,_{R}}=F_{1}(R)$ (30)

The function ${F_{1}(R)}$, as will be seen below, determines the initial
conditions for the dust velocity distribution.

Substituting (30) into Eq. (27), we find that

$\frac{q^{2}}{r^{2}}+8\pi\Lambda
r^{2}=\left(rr^{2}_{,_{\tau}}\right)_{,\tau}/r_{,_{\tau}}-F_{1}+1$ (31)

Multiplying this equation by $r_{,\tau}$ and integrating the result over time,
we arrive at

$\frac{q^{2}}{r}-\frac{8\pi}{3}\Lambda
r^{3}+r(1-F_{1})+rr^{2}_{,_{\tau}}=F_{2}(R)$ (32)

The function ${F_{2}(R)}$, as will be seen below, determines the initial
conditions for the dust density distribution.

##  V INITIAL CONDITIONS

Furthermore the subscript <<i>> will be assigned to all quantities at
${\tau=0}$.

We choose the scale of coordinate $R$ in (24) at the initial instant such that

$r_{i}:=R$ (33)

This can be done if the function r changes monotonically along the coordinate
$R$ at the initial instant.

We assume also that the entire matter is at rest at the initial instant and
that all derivatives with respect to $\tau$ are zero:

$\left.r_{,_{\tau}}\right|_{i}=0$ (34)

Then, proceeding from (32), we obtain

$F_{2}(R)=\frac{q^{2}}{R}-\frac{8\pi}{3}\Lambda R^{3}+R(1-F_{1})$ (35)

or

$F_{1}(R)=1-\frac{F_{2}}{R}+\frac{q^{2}}{R^{2}}-\frac{8\pi}{3}\Lambda R^{2}$
(36)

Multiplying Eq. (26) by ${r^{2}r_{,_{R}}}$ and expressing
${\lambda_{,_{\tau}}}$ from Eq. (28), we arrive at

$8\pi\varepsilon r^{2}r_{,_{R}}=\left[q^{2}/r-\frac{8\pi}{3}\Lambda
r^{3}+r(1-F_{1})+rr_{,_{\tau}}^{2}\right]_{,_{R}}$ (37)

A comparison of this expression with (32) shows that

$8\pi\varepsilon r^{2}r_{,_{R}}=\frac{dF_{2}}{dR}$ (38)

Based on Eq. (38), we can conclude that the quantity ${8\pi\varepsilon
r^{2}r_{,_{R}}}$ is independent of time:

$8\pi\varepsilon r^{2}r_{,_{R}}=8\pi\varepsilon_{i}r_{i}^{2}$ (39)

Let us introduce the mass function ${m(R)}$:

$m(R):=\int\limits_{0}^{R}4\pi\varepsilon_{i}R^{2}\,dR=\int\limits_{0}^{r_{i}}4\pi\varepsilon
r^{2}\,dr$ (40)

Within our model, there is neither matter nor electric field at radii smaller
than $r_{q}$ (except for the $\Lambda$-term), and there is a charged thin film
with a mass of $m_{q}$ at the radius $r_{q}$. Therefore444The superscripts
¡¡${-}$¿¿ and ¡¡${+}$¿¿ for the radius $r_{q}$ indicate the lower and upper
limits, i.e., respectively, under and above the charged film with a mass
$m_{q}$., ${m(r_{q}^{-})\equiv 0}$ и ${m(r_{q}^{+})\equiv m_{q}}$ and
expression (40) can be integrated not from zero but from ${r_{q}^{-}}$ or from
${r_{q}^{+}}$:

$\displaystyle
m(R>r_{q}^{-})=\int\limits_{r_{q}^{-}}^{R}4\pi\varepsilon_{i}R^{2}\,dR$ (41)
$\displaystyle
m(R>r_{q}^{+})=m_{q}+\int\limits_{r_{q}^{+}}^{R}4\pi\varepsilon_{i}R^{2}\,dR$
(42)

Integrating expression (38) in the limits from ${r_{q}^{-}}$ to $R$ and from
${r_{q}^{+}}$ to $R$, we find from (41-42) that

$\displaystyle 2m(R>r_{q}^{-})=F_{2}(R>r_{q}^{-})-F_{2}(r_{q}^{-})$ (43)
$\displaystyle 2m(R>r_{q}^{+})-2m_{q}=F_{2}(R>r_{q}^{+})-F_{2}(r_{q}^{+})$
(44)

To obtain the $F_{1}$ and $F_{2}$ values at the point $r_{q}$, we will take
into account that the charged sphere is in equilibrium; therefore, the
following equalities must be fulfilled: ${r_{,_{\tau}}=0}$ и
${r_{,_{\tau\tau}}=0}$. In addition, the electric field is absent under the
charged sphere (at $r_{q}^{-}$) (i.e., the terms with $q^{2}$ must be excluded
in that domain), whereas directly above the charged sphere the electric field
is ${q/r_{q}^{2}}$; therefore, the terms with $q^{2}$ on the charged sphere
must be multiplied by ${\frac{1}{2}}$. Then, proceeding from Eqs. (27), (30)
and (35), we find the ${F_{1}}$ and ${F_{2}}$ values on the charged sphere:

$\displaystyle F_{1}(r_{q})=1-\frac{q^{2}}{2r_{q}^{2}}-8\pi\Lambda r_{q}^{2}$
(45) $\displaystyle
F_{2}(r_{q})=\frac{q^{2}}{r_{q}}+\frac{16\pi\Lambda}{3}r_{q}^{3}$ (46)

As was noted in Section II, the matter mass ${m_{q}}$ has accumulated
completely at the charge radius; therefore, expressions for ${r_{q}^{+}}$ must
be used in the calculations. Then, proceeding from (42), (44) and (46), we
obtain the following expression for ${F_{2}(R>r_{q}^{+})}$:

$\displaystyle
F_{2}(R>r_{q}^{+})=\frac{q^{2}}{r_{q}}+\int\limits_{r_{q}^{+}}^{R}8\pi\varepsilon_{i}R^{2}\,dR+\frac{16\pi\Lambda}{3}r_{q}^{3}=\frac{q^{2}}{r_{q}}+2m_{dust}(R)+\frac{16\pi\Lambda}{3}r_{q}^{3}$
(47)

Let us introduce also the mass ${m^{\Lambda}(R)}$, related to the
$\Lambda$-term:

$m^{\Lambda}(R):=\int\limits_{0}^{R}4\pi\Lambda
R^{2}\,dR=\frac{4\pi}{3}\Lambda R^{3}$ (48)

It follows from (40) that the consideration of the cosmological
${\Lambda}$-term can be reduced to subtraction of ${4\pi\Lambda R^{3}/3}$ from
the mass${m(R)}$. According to the cosmological observations, ${\Lambda\sim
10^{-29}\makebox{g/sm}^{3}}$; therefore, on the scales on the order of
$r_{q}$, which are determined by expression (11), the mass related to the
$\Lambda$-term is about 8 orders of magnitude smaller than the total mass $M$
(see (10)). This means that the cosmological $\Lambda$-term can be disregarded
at ${R\sim r_{q}}$.

##  VI INTEGRATION OF THE EQUATIONS

The horizon radii ${r_{h}^{\pm}}$ in comoving coordinates are determined by
the equality of the modulus of invariant velocity $V$ to unity:

$\left.V^{2}\right|_{r_{h}^{\pm}}=1\,,\quad
V^{2}(\tau,R):=\frac{r_{,\tau}^{2}\,e^{\lambda}}{r_{,R}^{2}}\,.$ (49)

Here, the function ${V^{2}(\tau,R)}$ is the squared invariant velocity (see
A.A. Shatskiy, A.Yu. Andreev (1999); A.A. Shatskiy and Novikov (2010)). Using
expression (30), one can rewrite Eq. (49) in the form

$\left.r_{,\tau}^{2}\right|_{r_{h}^{\pm}}=F_{1}(R)$ (50)

Having expressed (using (36)) the function $F_{1}$ in Eq. (32) in terms of
$F_{2}$, we rewrite (32) in the form convenient for integrating and obtaining
the function ${r(\tau,R)}$:

$r_{,_{\tau}}=\pm\sqrt{\left(\frac{F_{2}}{r}-\frac{q^{2}}{r^{2}}+\frac{8\pi\Lambda
r^{2}}{3}\right)-\left(\frac{F_{2}}{R}-\frac{q^{2}}{R^{2}}+\frac{8\pi\Lambda
R^{2}}{3}\right)}$ (51)

This expression demonstrates, in particular, that the particle dynamics
depends on only the distribution of matter under the radius determining the
particle location and is independent of matter above this radius.

Proceeding from (51) and taking into account (36) and (50), we obtain an
equation for calculating the horizon radii ${r_{h}^{\pm}}$:

$1-\frac{F_{2}}{r_{h}^{\pm}}+\frac{q^{2}}{{r_{h}^{\pm}}^{2}}-\frac{8\pi}{3}\Lambda{r_{h}^{\pm}}^{2}=0$
(52)

This expression, with the $\Lambda$-term disregarded, is practically
equivalent to the equality ${f(r_{h}^{\pm})=0}$ from definition (17); i.e.,
the meaning of function $F_{2}$ is the doubled mass $m$ in the matter-comoving
coordinates. Therefore, we neglect below the $\Lambda$-term and derive an
analogue of expression (52) from (1) for the horizon radii in matter-comoving
coordinates:

$r_{h}^{\pm}=\frac{1}{2}F_{2}\pm\sqrt{\frac{1}{4}F_{2}^{2}-q^{2}}$ (53)

Hence, in correspondence with (3), we have

$q=\frac{\kappa F_{2}}{2}\,,\quad\kappa:=\frac{2\beta}{1+\beta^{2}}\,.$ (54)

Based on Eq. (51), one can obtain two stopping points for the dust shell;
these are turning points, at which ${r_{,_{\tau}}=0}$. Expression (51) yields
directly the trivial stopping point: the initial point ${r_{i}:=R}$, at which
the dust shell begins to move. Concerning the desired second stopping point,
we find from (51) that

$\displaystyle r_{t}=\frac{q^{2}}{F_{2}-q^{2}/R}=\frac{\kappa q}{2-\kappa
q/R}$ (55)

Hence, one can see that ${r_{t}\geq 0.5\kappa q}$; the equality is implemented
only at ${q/R\to 0}$.

Let us first consider the motion of probe particles under an external dust
shell. Probe particles cannot form a horizon; according to (47) and (55), we
have the following relations for them:

$\displaystyle F_{2}^{0}=\frac{q^{2}}{r_{q}}\,,\quad
r_{t}^{0}=\frac{r_{q}R}{R-r_{q}}\,.$ (56)

Therefore, in the case ${r_{q}<R<2r_{q}}$ for probe particles, the stopping
radius $r_{t}^{0}$ even exceeds the initial radius $r_{i}$, which indicates
the outside direction of their motion. Correspondingly, in the range
${2r_{q}<R<\infty}$ the probe particles will begin to move inside. For the
stopping point of a probe particle with ${r_{t}^{0}\to r_{q}}$, we find from
(56) that ${R\to\infty}$. Thus, the spherical layers of probe particles will
constantly cross each other during their motion.

Let us now consider the non-probe particles of a dust sphere of mass
${m_{dust}}$, whose radius is $R_{1}$ at the initial instant. In accordance
with (5), the condition for location of turning point radius $r_{t}$ for dust-
sphere particles between the charge-sphere radius $r_{q}$ and Cauchy horizon
radius $r_{h}^{-}$ takes the form

$\displaystyle r_{q}=\beta\gamma q<r_{t}<r_{h}^{-}=\beta q\,,\quad m_{q}<q\,.$
(57)

Introducing the designation

$\displaystyle h:=\frac{\kappa q}{2R_{1}}$ (58)

and taking into account (22) and (55), we rewrite (57) as

$\displaystyle\frac{r_{q}}{q}=\frac{m_{q}q(2m_{q}^{2}-q^{2})+q^{3}\sqrt{2m_{q}^{2}-q^{2}}}{4m_{q}^{4}-q^{4}}<\frac{r_{t}}{q}=\frac{\beta}{\left(1+\beta^{2}\right)\left(1-h\right)}<\beta\,,\quad\frac{q}{m_{q}}>1\,.$
(59)

The range of parameters at which these inequalities are satisfied can be
determined from the plots presented in Fig. 2.

Figure 2:  (Left) Dependence of ${r_{q}/q}$ on ${q^{2}/m_{q}^{2}}$.
(Right) Dependences of ${r_{t}/q}$ on $\beta$ at different $h$ for inequality
(59).
.

Here, we should also consider a seeming contradiction related to possible
intersection of dust layers. As was shown for probe particles, they do not
reach the radius $r_{q}$ at finite $R$ values. The same holds true for the
dust-shell inner layers, because the gravitational potential for these layers
is still insufficient even to form BH horizons. Therefore, the dust-shell
inner layers will begin to move from the center at some instant. At the same
time, the dust-shell outer layers have already a sufficiently high
gravitational potential to form BH horizons and reach the turning point
${r_{t}}$. Thus, the inner layers will cross the outer layers and their roles
will interchange at some instant. Theoretically, our model becomes incorrect
after the intersection of non-probe dust layers; however, we deal in practice
with a fairly thin dust shell, whose thickness was chosen by us
incommensurately smaller than the radius $r_{q}$ (see (11)). Therefore, both
under and above this thin dust shell, we can neglect the inner intersection of
its layers, whereas inside the dust shell there will be constant mixing of its
inner and outer layers. On the whole, this pattern should not affect the model
dynamics. In addition, we can slightly modify our model and assume that bulk
dust shells consist of only one layer; the size of each <<dust particle>> in
this layer is on the order of the entire dust shell thickness ${d_{dust}}$
(${d_{dust}<<r_{q}}$). In this case, there will be no intersection of dust
layers, and each bulk <<dust particle>> will move in a self-consistent
gravitational field. According to the estimates reported in Section II,
${m_{dust}\sim m_{q}\sim M}$ in our model; therefore, at a density of <<dust
particle>> material of about ${\sim 1\makebox{g/sm}^{3}}$, one can obtain a
limitation on the <<dust particle>> size: ddust ${d_{dust}\lesssim
d_{q}\approx 70}$. Therefore, the <<dust particles>> in the dust layer under
consideration can be small asteroids or boulders.

##  VII CHOICE OF APPROPRIATE PARAMETERS

Taking into account (23), we introduce the following designation:

$\displaystyle\alpha:=\frac{q}{m_{q}}>1$ (60)

To implement the electrostatic and gravitational equilibrium of charge $q$ on
a sphere with a mass $m_{q}$ and radius $r_{q}$, one must choose the
parameters $\alpha$, $\beta$ and $h$ such as to satisfy inequalities (59).
This can be done using, e.g., the parameters

$\displaystyle\alpha^{2}=1.1\,,\quad\beta=0.9\,,\quad h=0.4$
$\displaystyle\Rightarrow\,\frac{r_{q}}{q}\approx 0.731<\frac{r_{t}}{q}\approx
0.829\,,\quad\kappa\approx 0.997\,,\quad\frac{R_{1}}{r_{q}}\approx 1.706\,.$
(61)

or

$\displaystyle\alpha^{2}=1.01\,,\quad\beta=0.8\,,\quad h=0.3$
$\displaystyle\Rightarrow\,\frac{r_{q}}{q}\approx 0.673<\frac{r_{t}}{q}\approx
0.697\,,\quad\kappa\approx 0.988\,,\quad\frac{R_{1}}{r_{q}}\approx 2.447\,.$
(62)

Thus, the problem of calculating the equilibrium state of charge $q$ on a
sphere can be considered as solved.

##  VIII ACCRETION OF MATTER ON ALREADY EXISTING REISSNER–NORDSTRÖM BH

We considered above the accretion of a dust shell on a charged sphere in
equilibrium, in the absence of BH. Let us now consider the accretion of a dust
shell on a Reissner–Nordström BH.

In accordance with (1), any charged BH with a mass $M$ and charge $q$ obeys
the inequality

$r_{h}^{-}\leq q\leq M\leq r_{h}^{+}$ (63)

Hence, if a BH already exists at the initial instant, its charge remains
constant, and only its mass increases, the following inequality will always be
valid for any two instants <<$a$>> and <<$b$>>:

$r_{ha}^{-}\leq r_{hb}^{+}$ (64)

Here, the subscripts <<$a$>> and <<$b$>> correspond to the radii of horizons
the BH will have at the corresponding (arbitrary) instants. In other words,
the case where the inner Cauchy horizon at some instant <<$a$>> exceeds the BH
outer horizon at some other instant <<$b$>> is impossible. Hence, independent
of the specific features of neutral matter accretion on the BH, the latter
will have only one turning point. Therefore, the situation where another
turning point would arise above the initial BH horizon (as a result of the
formation of another outer $T$-region) is impossible.

##  IX DISCUSSION

Thus, the formation of another turning point for a charged BH is impossible
after the formation instant of this BH. Therefore, a necessary initial
condition for the turning point formation is the absence of BH at the initial
instant. Hence, the occurrence of a new BH and a new universe (along with the
occurrence of a white hole) is possible only at the instant when horizons of a
new Reissner–Nordström BH (or, possibly, Kerr BH) arise.

Therefore, the above-considered mechanism of turning-point formation for a
Reissner–Nordström BH is the only possible one for this BH. In other words,
within our model, there is a charged sphere at the initial instant; the charge
and entire matter of this sphere are concentrated in a thin charged spherical
film. Beyond this sphere, there is no matter, up to some radius
${r_{i}=R_{1}}$. At this radius, $R_{1}$, a spherical dust shell with a mass
${m_{dust}}$ is at rest; the thickness ${m_{dust}}$ of this sphere is
neglected (in comparison with the radius $r_{q}$). Horizons are absent at the
initial instant (there is no BH). We are not interested in the matter beyond
the radius ${r_{i}=R_{1}}$.

At subsequent instants the dust sphere begins to fall to the center (towards
the charged sphere) under the action of gravity. At some instant the dust
sphere reaches its gravitational radius ${r_{h}^{+}}$ (see (53)) and finds
itself under the horizon of newly formed Reissner–Nordström BH at the next
instant.

When the dust sphere reaches the turning-point radius, according to (49), its
velocity $V$ decreases to zero, after which the dust sphere bounces to another
universe. This is due to the fact that the dust sphere is already located in
the inner R region (under both BH horizons). Therefore, the dust sphere cannot
escape backwards, and the solution indicates that its radius should increase
after the bounce. The only possibility is the escape of the dust sphere to
another universe. This event occurs at the instant of the formation of a new
white hole and new expanding universe. Apparently, the universe to which the
dust sphere escapes arises jointly with the formation of both horizons of
Reissner–Nordström BH.

One of the most interesting conclusions that can be drawn is that the falling
dust sphere in the model under consideration does not reach the charged sphere
of radius $r_{q}$. The reason is that the turning point for the dust sphere
with $r_{t}$ is located beyond the charged sphere. Therefore, a
Reissner–Nordström black-andwhite hole with all its horizons and R-T-regions
arises beyond the charged sphere. Apparently, nothing should change within the
charged sphere of radius $r_{q}$: the electric field and singularity are
absent even after the formation of a Reissner–Nordström BH outside the sphere.
The metric coefficients in the model should remain everywhere such as to
provide constant and continuous matching of two regions: within the charged
sphere and beyond it.

Thus, we proposed a model in which both the black-and-white hole and its
turning point are formed in already existing Universe; i.e., the black-and-
white hole in this model is not <<eternal>>. The analytical solution of a
similar self-consistent problem of accretion for a real rotating Kerr BH would
be practically impossible in view of its complexity. However, one can suggest
that a new Kerr black-and-white hole is formed according to a similar
mechanism during the rotating-star collapse. An important difference from the
dynamic Reissner–Nordström solution is that the Kerr solution does not contain
any limitations (see Section II) on the electric field magnitude and,
correspondingly, on the sizes and mass of black-and-white hole. The specific
features of the arising complex topology in a Kerr black-and-white hole were
considered, e.g., in Shatskiy (2020); B. Carter (1966, 1968).

ACKNOWLEDGMENTS

I am grateful to the organizers and all members of Zelmanov Memorial Seminars
in Gravitation and Cosmology at the Shternberg State Astronomical Institute of
Moscow State University for the questions, discussion, and fruitful remarks.

## References

  * Ch. Misner, K. Thorne, and J. Wheeler (1977) Ch. Misner, K. Thorne, and J. Wheeler, _Gravitation, vol. 3_ ((Freeman, San Francisco, 1973; Mir, Moscow, 1977).
  * D.A. Tretyakova, A.A. Shatskiy, I.D. Novikov, S. Alexeyev (2012) D.A. Tretyakova, A.A. Shatskiy, I.D. Novikov, S. Alexeyev, Phys. Rev. D 85, 124059 (2012).
  * A.A. Shatskiy and Novikov (2010) D. N. A.A. Shatskiy, A.G. Doroshkevich and I. Novikov, Journal of Experimental and Theoretical Physics 110, 268–276 (2010).
  * V.A.Berezin and V.I.Dokuchaev (2014) V.A.Berezin and V.I.Dokuchaev (2014), eprint gr-qc/1404.2726, gr-qc/1404.2727.
  * Shatskiy (2020) A. Shatskiy, Journal of Experimental and Theoretical Physics 130, 409–417 (2020).
  * (6) L.D. Landau, E.M. Lifshitz, _Course of Theoretical Physics, Volume 2: The Classical Theory of Fields_ (Nauka, Moscow; Oxford, 1990, 1988)
  * Tolman (1934) R. C. Tolman, Proc. Nat. Acad. Sci. US 20, 169 (1934).
  * Oppenheimer (1939) J. R. Oppenheimer, Phys. Rev. 56, 455 (1939).
  * A.A. Shatskiy, A.Yu. Andreev (1999) A.A. Shatskiy, A.Yu. Andreev, Journal of Experimental and Theoretical Physics 116, 353 (1999).
  * B. Carter (1966) B. Carter, Phys. Lett. 21, 423 (1966).
  * B. Carter (1968) B. Carter, Comm. Math. Phys. 10, 280 (1968).