# VelCrys: Interactive web-based application to compute acoustic wave velocity
in crystals and its magnetic corrections

P. Nieves  I. Korniienko  A. Fraile  J.M. Fernández  R. Iglesias  D. Legut
Departamento de Física, Universidad de Oviedo, C. Leopoldo Calvo Sotelo, 18,
33007, Oviedo, Spain  IT4Innovations, VŠB - Technical University of Ostrava,
17. listopadu 2172/15, 70800 Ostrava, Czech Republic  Material Physics Center
(MPC) / Centro de Física de Materiales(CFM) CSIC-UPV/EHU, San Sebastián, Spain

###### Abstract

We present VelCrys, a web-based interactive tool, that allows to perform
further post-processing of the elastic tensor in order to compute and plot the
group velocity of the acoustic waves for any crystal symmetry. We also
implemented the calculation of effective magnetic corrections to the elastic
tensor and corresponding fractional change in group velocity under a magnetic
field. We apply it to dry sandstone, cubic CoPt and hcp Co to show some of the
program features. In the analysis of magnetic corrections, we find complex
landscapes of fractional change in group velocity as a function of ray
direction, as well as a field dependence consistent with Simon effect.

###### keywords:

Acoustic , Crystals , Magnetism

##  1 Introduction

The characterization of wave velocity is particularly useful in many fields
like geophysics[1, 2, 3, 4, 5], communication technology[6] or biosensors[7].
The way acoustic waves propagate in crystals can be highly anisotropic, and
thus it is different and more complicated than in isotropic media [8].
Providing the phase velocity direction ($\bm{k}$), one finds three possible
modes for an acoustic wave: one quasi-compressional wave (qP) and two quasi-
shear waves (qS1 and qS2)[1]. Over the last five decades, great efforts have
been made to develop theoretical techniques to compute wave velocity for a
general anisotropic media[9, 10, 11, 12, 1]. For instance, Zhou et al.[12]
proposed the eigenvalue and eigenvector methods to calculate the group
velocity. A more manageable approach has been recently proposed by Wang et
al.[1] based on analytical formulas for group velocity, which can deal with
the singular points and multivalue problems of the qS wave.

Magnetic effects on wave velocity in magnetostrictive materials are
interesting from a fundamental and applied point of view[13, 14, 15, 16]. The
effect of a magnetic field on wave velocity was originally discovered by Simon
in 1958 (Simon effect)[17]. Additionally, other magnetic effects have been
identified like isotropic exchange effects[18], anisotropic morphic
effects[19, 20, 13], rotational-magnetostrictive effects [13, 14], and $\Delta
s$-effect [14], which have been experimentally observed in ultrasonic elastic
waves[14, 21, 22, 23, 24, 25, 15, 26, 27, 16]. In 1979, a comprehensive
theoretical description of these effects were reported by Rouchy et al. for
cubic crystals [13]. Later, Rinaldi and Turilli showed that linear magnetic
effects on wave velocity can be equivalently taken into account as effective
corrections to the elastic constants [28]. The addition of these corrections
to the elastic tensor combined with Wang et al. method [1] to compute group
velocity may provide a practical way to efficiently compute the fractional
change in group velocity due to linear magnetic effects for arbitrary wave
propagation direction and applied magnetic field. Hence, it may facilitate the
analysis of this phenomenon beyond the cases with wave propagation and
magnetic field directions along high symmetry crystallographic directions
studied so far, as well as crystals with symmetry lower than cubic [13, 22,
18, 15]. Such kind of approach may be exploited in novel magneto-acoustic
phenomena like acoustic spin pumping [29], magnetization switching induced by
sound waves [30, 31, 32, 33, 34, 35, 36, 37] and picosecond ultrasonics [38].

The development advanced visualization tools in databases can help to
accelerate the design of new materials in different areas of condensed matter
physics[39, 40, 41]. For example, ELATE is an online user-friendly interface
that has been implemented to analyze the elastic tensor and visualize
anisotropic mechanical properties [42]. Similarly, MAELASviewer was developed
to visualize anisotropic magnetostrictive properties[43, 44, 45]. In this
work, we present VelCrys a web-based interactive tool that performs further
post-processing of the elastic tensor in order to compute and plot the sound
velocity of the three kinds of acoustic waves (qP, qS1, and qS2) for any
crystal symmetry based on the analytical formulas derived by Wang et al.[1].
Additionally, if the magnetoelastic constants are also provided, then it can
calculate the effective magnetic corrections to the elastic tensor and
corresponding fractional change in sound velocity under a magnetic field using
Rinaldi and Turilli approach[28].

##  2 Methodology

###  2.1 Software details

VelCrys is created with Dash which is an efficient Python framework for
building web applications [46]. Dash uses Flask as the web framework [47].
Visualization of 3D surfaces is done using the Python open source graphing
library Plotly [48]. We also make use of NumPy library for some mathematical
operations. The online application is currently available at http://www.md-
esg.eu/velcrys/, while the open source Python module is available at GitHub
repository [49]. VelCrys can be also used offline by executing the Python
module[49].

###  2.2 Calculation of sound velocity and magnetic corrections

The calculation of group velocity $v$ for any crystal symmetry is performed by
computing the analytical solution derived by Wang et al.[1] which makes use of
the elastic tensor $C_{ij}$ as input parameter

$v_{l}(C_{ij})=\frac{\partial c}{\partial
n_{l}}=\frac{\partial\omega}{\partial k_{l}},\quad l=x,y,z$ (1)

where $c=\omega/k$ is the phase velocity, and $\bm{n}=\bm{k}/k$ is the unit
vector of phase velocity direction (ray direction).

The effective magnetic corrections ($\Delta C_{ij}$) to elastic tensor are
obtained through the formulas derived by Rinaldi and Turilli which have the
form [28]

$\Delta
C_{ij}=\frac{1}{M_{s}^{2}}b_{pqi}b_{rsj}\alpha_{p}\alpha_{r}\chi_{qs},$ (2)

where $b$ are magnetoelastic constants, $M_{s}$ is saturation magnetization,
$\alpha$ is direction cosine of magnetization at equilibrium, and $\chi$ is
the susceptibility tensor, see Ref.[28] for details about subscript indexes.
The corresponding fractional change in group velocity due to magnetic effects
$(v-v_{0})/v_{0}$ is obtained by computing Wang et al. formulas[1] using the
elastic tensor including effective magnetic corrections ($C_{ij}+\Delta
C_{ij}$), that leads to $v$, and elastic tensor without effective magnetic
corrections ($C_{ij}$) which gives $v_{0}$, that is

$\frac{v-v_{0}}{v_{0}}=\frac{v(C_{ij}+\Delta C_{ij})-v(C_{ij})}{v(C_{ij})}.$
(3)

This option is currently supported only for Cubic I (point groups $432$,
$\bar{4}3m$, $m\bar{3}m$, space group numbers from 207 up to 230), and
Hexagonal I (point groups $6mm$, $622$, $\bar{6}2m$, $6/mmm$, space group
numbers from 177 up to 194) symmetries[28]. Note also that Rinaldi and Turilli
expressions does not include higher order corrections like morphic and
rotational-magnetostrictive effects [18]. These limitations are discussed in
Section 3.2.

##  3 Results

In this section, we perform some tests of the methodology implemented in
VelCrys.

###  3.1 Dry sandstone

Sandstone is sedimentary rock widely distributed in planet Earth mostly
composed of sand-sized minerals or rock grains cemented by siliceous,
argillaceous, and/or calcareous materials [50]. Its mechanical properties
depend highly upon the degree of cementation and grain composition [50] which
are affected by temperature and burial history [50, 51]. In this first
example, we verify the correct implementation in VelCrys of Wang et al.[1, 52]
methodology by computing the group velocity of complex dry sandstone and
compare with their results[1]. The elastic constants of this material are
shown in Table 1, while the mass density is $\rho=2080$ kg/m3 at room
temperature[1]. In Fig.1, we compare the group velocity for a qP wave obtained
with VelCrys and reported by Wang et al.[1], finding good agreement. A 3D plot
generated with VelCrys is shown in Fig.2, where one can see the group velocity
landscape of the qP wave as a function of the ray direction (phase velocity
direction), again consistent with Ref. [1].

Table 1: Elastic constants ($C_{ij}=C_{ji}$) of dry sandstone [1].

$C_{ij}$ | (GPa) | $C_{ij}$ | (GPa) | $C_{ij}$ | (GPa) | $C_{ij}$ | (GPa) | $C_{ij}$ | (GPa) | $C_{ij}$ | (GPa)
---|---|---|---|---|---|---|---|---|---|---|---
$C_{11}$ | 10.3 | $C_{12}$ | 0.9 | $C_{13}$ | 1.3 | $C_{14}$ | 1.4 | $C_{15}$ | 1.1 | $C_{16}$ | 0.8
|  | $C_{22}$ | 10.6 | $C_{23}$ | 2.1 | $C_{24}$ | 0.2 | $C_{25}$ | -0.2 | $C_{26}$ | -0.6
|  |  |  | $C_{33}$ | 14.1 | $C_{34}$ | 0.0 | $C_{35}$ | -0.5 | $C_{36}$ | -1.0
|  |  |  |  |  | $C_{44}$ | 5.1 | $C_{45}$ | 0.0 | $C_{46}$ | 0.2
|  |  |  |  |  |  |  | $C_{55}$ | 6.0 | $C_{56}$ | 0.0
|  |  |  |  |  |  |  |  |  | $C_{66}$ | 4.9

Figure 1: Isoline of the group velocity ($v=2.6$ km/s) of the qP wave for dry
sandstone with varied ray angle ($\phi,\theta$) obtained with VelCrys (red
open circles) and by Wang et al.[1] (blue solid line).  Figure 2: Group
velocity landscape as a function of the ray direction (phase velocity
direction)
$\bm{n}=\bm{k}/k=(\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta)$ of the qP
wave for dry sandstone generated with VelCrys. The color of the surface
corresponds to the magnitude of group velocity $|v|$. Table 2: Experimental
elastic constants ($C_{ij}$), anisotropic magnetoelastic constants ($b$),
anisotropic magnetostrictive coefficients ($\lambda$), magnetocrystalline
anisotropy constant ($K_{1}$), saturation magnetization ($M_{s}$) and mass
density ($\rho$) for cubic single crystals of disordered equiatomic Co-Pt
alloy at room temperature [14, 53].

$C_{ij}$ | (GPa) | $b$ | (MPa) | $\lambda$ | ($\times 10^{-6}$) |  $K_{1}$(kJ/m3) |  $\mu_{0}M_{s}$(T) |  $\rho$(kg/m3)
---|---|---|---|---|---|---|---|---
$C_{11}$ | 289.7 | $b_{1}$ | -35 | $\lambda_{100}$ | 210 | -40 | 0.942 | 16696.7
$C_{12}$ | 178.5 | $b_{2}$ | 11.9 | $\lambda_{111}$ | -53 |  |  |
$C_{44}$ | 124.1 |  |  |  |  |  |  |

Table 3: Theoretical expression for each magnetic contribution to fractional
change in group velocity $(v-v_{0})/v_{0}$ in cubic crystal [13, 14, 18].
Here, $\bm{u}$ corresponds to the displacement vector, while other quantities
are defined in the main text.

$\bm{k}\parallel$ | $\bm{u}\parallel$ | $\bm{H}\parallel$ | $G(m)$ | $R(\lambda)$ | $S(H)$ | $(l-l_{0})/l_{0}$ | Wave
---|---|---|---|---|---|---|---
$[001]$ | $[001]$ | $[001]$ | $\frac{2m_{1}^{\gamma,2}}{3C_{11}}$ | 0 | 0 | $\lambda_{100}$ | qP
$[001]$ | $[001]$ | $[100]$ | $-\frac{m_{1}^{\gamma,2}}{3C_{11}}$ | 0 | 0 | $-\frac{\lambda_{100}}{2}$ | qP
$[001]$ | $[100]$ | $[001]$ | $-\frac{m_{3}^{\gamma,2}}{12C_{44}}$ | $\frac{3}{2}(\lambda_{111}-\lambda_{100})$ | $-\frac{(b_{2})^{2}}{2C_{44}\mu_{0}M_{s}\left(H+H_{D}+\frac{2K_{1}}{\mu_{0}M_{s}}\right)}$ | $\lambda_{100}$ | qS1
$[001]$ | $[100]$ | $[100]$ | $-\frac{m_{3}^{\gamma,2}}{12C_{44}}$ | $\frac{3}{2}(\lambda_{100}-\lambda_{111})$ | $-\frac{(b_{2})^{2}}{2C_{44}\mu_{0}M_{s}\left(H+M_{s}+H_{D}+\frac{2K_{1}}{\mu_{0}M_{s}}\right)}$ | $-\frac{\lambda_{100}}{2}$ | qS1

Figure 3: Fractional change in ultrasound frequency for a longitudinal wave
(qP) under an applied magnetic field parallel or perpendicular to the wave
propagation vector for cubic single crystals of disordered Co-Pt alloy. Red
stars and orange diamonds show experimental data from Ref.[14] for the cases
$H\parallel[001]$ and $H\parallel[100]$, respectively. Red solid and orange
dashed lines stand for the theoretical values obtained from Eq.5 and Table 3.
Open blue squares and solid pink circles give the values derived from the
computed fractional change in group velocity $(v-v_{0})/v_{0}$ by VelCrys via
Eq.4. We added experimentally observed $\Delta s$-effect at saturated state
($\Delta^{sat}=0.00435$) both in the theoretical values and VelCrys data.
Figure 4: Fractional change in ultrasound frequency for a transverse wave
(qS1) under an applied magnetic field parallel or perpendicular to the wave
propagation vector for cubic single crystals of disordered equiatomic Co-Pt
alloy. Red stars and orange diamonds show experimental data from Ref.[14] for
the cases $H\parallel[100]$ and $H\parallel[001]$, respectively. Red solid and
orange dashed lines stand for the theoretical values obtained from Eq.5 and
Table 3. Open blue squares and solid pink circles give the values derived from
the computed fractional change in group velocity $(v-v_{0})/v_{0}$ by VelCrys
via Eq.4. We added experimentally observed $\Delta s$-effect at saturated
state ($\Delta^{sat}=0.00435$) both in the theoretical values and VelCrys
data. Figure 5: Landscape of the fractional change in group velocity
$(v-v_{0})/v_{0}$ for the 3 types of waves (qP, qS1 and qS2) as a function of
the ray direction (phase velocity direction) with applied magnetic field
$\mu_{0}\bm{H}=(0,0,5)$ T for cubic single crystals of disordered equiatomic
Co-Pt alloy generated with VelCrys.

###  3.2 Co-Pt alloy

Co-Pt (and Fe-Pt) are very interesting alloys since supported nanostructures
of these alloys present direct link between atom arrangements and magnetic
behavior. Additionally, both alloys are today model systems in the field of
nanoalloys, due to the diversity of atom arrangements present either in bulk
state or specific to the nanoscale (disordered fcc structures, chemically
ordered L10, L12, core-shell, five-fold structures - icosahedral or
decahedral, etc). Importantly, Co-Pt alloys that can be used as permanent-
magnet materials are found around the equiatomic composition where the
disordered A1 (fcc) phase transforms into the ordered L10 phase [54]. Here, we
further test VelCrys capabilities by evaluating magnetic corrections to group
velocity for cubic single crystals of disordered equiatomic Co-Pt alloy, and
comparing the results with theory and experiment [14, 53]. The fractional
change in group velocity $(v-v_{0})/v_{0}$ due to magnetic effects may be
written as [13, 18]

$\frac{v-v_{0}}{v_{0}}=\frac{F-F_{0}}{F_{0}}+\frac{l-l_{0}}{l_{0}},$ (4)

where $F$ is the frequency of sound wave measured in ultrasound experiments
(for instance, by the pulse echo method[14, 22]) and $l$ is the length in the
direction of wave propagation. The fractional change of frequency could be
further splitted into [13, 18, 15, 16, 14]

$\frac{F-F_{0}}{F_{0}}=G(m)+R(\lambda)+S(H)+\Delta,$ (5)

where $G(m)$ accounts for higher order correction arising from the morphic
coefficients $m$, $R(\lambda)$ comes from rotational and magnetostrictive
effects, $S(H)$ describes the Simon effect that depends on the external
magnetic field $\bm{H}$ [17], and $\Delta$ corresponds to other possible
effects like the $\Delta s$-effect associated to domain walls displacement in
the low magnetic field regime [14]. The anisotropic part of the fractional
change in length reads [18, 15]

$\displaystyle\frac{l-l_{0}}{l_{0}}\Bigg{|}_{\bm{\beta}}^{\bm{\alpha}}$
$\displaystyle=\lambda^{\gamma,2}\left(\alpha_{x}^{2}\beta_{x}^{2}+\alpha_{y}^{2}\beta_{y}^{2}+\alpha_{z}^{2}\beta_{z}^{2}-\frac{1}{3}\right)$
(6)
$\displaystyle+2\lambda^{\epsilon,2}(\alpha_{x}\alpha_{y}\beta_{x}\beta_{y}+\alpha_{y}\alpha_{z}\beta_{y}\beta_{z}+\alpha_{x}\alpha_{z}\beta_{x}\beta_{z}),$

where $\bm{\alpha}$ is the direction cosine of magnetization, $\bm{\beta}$ is
the normalized vector along the measuring length direction, and
$\lambda^{\gamma,2}=3\lambda_{100}/2$ and
$\lambda^{\epsilon,2}=3\lambda_{111}/2$ are the magnetostrictive coefficients.
In this example, we consider longitudinal (qP) and transverse (qS) elastic
waves propagating along the crystallographic direction $[001]$
($\bm{k}\parallel[001]$), where the magnetic field is applied either along
$\bm{H}\parallel[001]$ or $\bm{H}\parallel[100]$. In Table 3, we present the
theoretical expression for each magnetic contribution to fractional change in
group velocity $(v-v_{0})/v_{0}$ in cubic crystal [13, 14, 18] for these
cases. A comparison between experiment, theory and VelCrys is shown in Figs.3
and 4, finding consistent results. In function $G(m)$, we use experimental
values of the morphic coefficients $m_{1}^{\gamma,2}=-0.15$ GPa, and
$m_{3}^{\gamma,2}=0.77$ GPa, while for the other properties we use the
experimental values shown in Table 2[14]. The fractional change in frequency
$F$ is indirectly obtained in VelCrys by using the calculated fractional
change in group velocity in Eq.4. Here, the experimental observed $\Delta
s$-effect, about $\Delta=0.00435$ at saturated state[14], is included both in
theory and VelCrys. For the sake of simplicity, the demagnetizing field was
not included ($H_{D}=0$) in Simon effect $S(H)$. Note that VelCrys gives the
same value for the cases with $\bm{H}\parallel\bm{k}$ and $\bm{H}\perp\bm{k}$,
failing to reproduce the small observed shift $\delta$ between both cases. In
the case of transverse waves (Fig.4), this discrepancy is due to the fact that
Rinaldi and Turilli work[28] is based on infinitesimal strain theory, while
this shift arises from the rotational-magnetostrictive effect $(\delta=\delta
R=3[\lambda_{100}-\lambda_{111}]=0.0008)$ based on the finite strain theory
[13, 14, 15, 18]. Similarly, in the case of longitudinal wave (Fig.3), this
discrepancy is due to the lack of higher order corrections coming from morphic
coefficients ($\delta=\delta G=|m_{1}^{\gamma,2}|/C_{11}=0.0005$) [18, 15].
The low field regime below magnetic saturation is influenced by $\Delta
s$-effect due to domain walls displacement[14], and cannot be accounted by
Simon effect described by function $S(H)$ in Table 3, as well as VelCrys,
since they assume a magnetic saturated state[17, 28]. In Fig.5, we show the
landscape of the fractional change in group velocity for the 3 types of waves
(qP, qS1 and qS2) as a function of the ray direction (phase velocity
direction) with applied magnetic field $\mu_{0}\bm{H}=(0,0,5)$ T generated
with VelCrys. We see that for the qP wave we have $(v-v_{0})/v_{0}=0$ at
$\bm{k}\parallel[001]$ consistent with the Simon effect ($S(H)=0$), see Table
3. The largest value is observed in the qS1 wave about
$(v-v_{0})/v_{0}=-450\times 10^{-6}$. Note that the values shown in Fig.5 does
not include the $\Delta s$-effect. Theoretical expression for other high
symmetry ray and field directions were derived by Rouchy et al.[13, 14], which
can be also studied with VelCrys.

Figure 6: (Top) Fractional change in group velocity $(v-v_{0})/v_{0}$ for the
qP wave with ray direction $\bm{k}\parallel[001]$ as a function of the applied
magnetic field $\bm{H}\parallel[100]$ for hcp Co generated with VelCrys.
(Bottom) The $x$ and $z$ components of direction cosines of equilibrium
magnetization direction as function of the applied magnetic field. Figure 7:
Landscape of the fractional change in group velocity $(v-v_{0})/v_{0}$ for the
3 types of waves (qP, qS1 and qS2) as a function of the ray direction (phase
velocity direction) with applied magnetic field $\mu_{0}\bm{H}=(5,0,0)$ T for
hcp Co generated with VelCrys.

###  3.3 hcp Co

Most of the Co produced is used for special alloys [55]. A relatively large
percentage of the world’s production goes into magnetic alloys such as the
Alnicos for permanent magnets. Sizable quantities are utilized for alloys that
retain their properties at high temperatures and superalloys that are used
near their melting points (where steels would become too soft). Co is also
employed for hard-facing alloys, tool steels, low-expansion alloys (for glass-
to-metal seals), and constant-modulus (elastic) alloys (for precision
hairsprings). In this last example, we study magnetic corrections to group
velocity for Co with hexagonal symmetry (hcp Co). We point out that magnetic
field effect on group velocity is still rather unexplored in hexagonal
crystals [28, 56]. The used parameters in VelCrys for this material are given
in Table 4.

Table 4: Experimental elastic constants ($C_{ij}$), anisotropic magnetoelastic
constants ($b$), magnetocrystalline anisotropy constants ($K$), saturation
magnetization ($M_{s}$) and mass density ($\rho$) for hcp Co [57, 58].

$C_{ij}$ | (GPa) | $b$ | (MPa) | $K$ | (kJ/m3) |  $\mu_{0}M_{s}$(T) |  $\rho$(kg/m3)
---|---|---|---|---|---|---|---
$C_{11}$ | 307 | $b_{21}$ | -31.9 | $K_{1}$ | 410 | 1.76 | 8900
$C_{12}$ | 165 | $b_{22}$ | 25.5 | $K_{2}$ | 150 |  |
$C_{13}$ | 103 | $b_{3}$ | -8.1 |  |  |  |
$C_{33}$ | 358 | $b_{4}$ | 42.9 |  |  |  |
$C_{44}$ | 75 |  |  |  |  |  |

In Fig.6, we computed the fractional change in group velocity for the qP wave
with ray direction $\bm{k}\parallel[001]$ as a function of the applied
magnetic field $\bm{H}\parallel[100]$. At high fields ($\mu_{0}H>2$T), we
observe a field dependence $\sim 1/H$ consistent with Simon effect. In the low
field regime ($\mu_{0}H<2$T), a non-monotonic behaviour is observed probably
due to an interplay between magnetic anisotropy, applied magnetic field, and
rotation of magnetization. Similar peaks have been observed in Co-Pt
alloys[14]. In Fig.7, we present the landscape of the fractional change in
group velocity for the 3 types of waves (qP, qS1 and qS2) as a function of the
ray direction (phase velocity direction) with applied magnetic field
$\mu_{0}\bm{H}=(5,0,0)$ T generated with VelCrys. Largest values are found in
the qS2 wave ($\Delta v/v_{0}\sim-0.006$) along the direction of the magnetic
field $\bm{k}\parallel\bm{H}\parallel[100]$. More theoretical analysis and
experimental measurements in hexagonal crystals would be highly desirable to
have a more deep understanding of these results [56]. Note that the same
limitations of the implemented methodology as in the case of cubic crystals
(lack of higher order corrections and $\Delta s$-effect) are also expected for
the hexagonal symmetry.

##  4 Conclusion

In summary, we showed that the combination of Wang et al.[1] analytical
solution of group velocity for a general elastic tensor with Rinaldi and
Turilli[28] effective magnetic corrections could provide a convenient approach
to easily compute the group velocity and fractional change in group velocity
due to magnetic effects. Such strategy has been implemented in the open source
program VelCrys, facilitating the calculation of these quantities. VelCrys
also includes visualization features to generate 3D plots of the landscape of
these quantities as a function of the ray direction.

We have validated the program by recovering Wang et al. results of group
velocity for the dry sandstone [1]. The magnetic corrections[28] were tested
on cubic disordered Co-Pt alloy showing good agreement with theory and
experiment after adding the $\Delta$s-effect[14]. We also applied this method
to hcp Co finding a field dependence consistent with Simon effect[17], as well
as complex landscapes of fractional change in group velocity that have yet to
be fully explained. Current version of the program has some limitations. For
example, magnetic corrections to the elastic tensor has been implemented only
for cubic I and hexagonal I symmetries. Additionally, these corrections does
not include higher order terms associated to morphic coefficients and
rotational-magnetostrictive effect, as well as $\Delta s$-effect related to
domain walls displacement. Hopefully, some of these limitations might be
overcome in future releases of the program.

Potential applications of VelCrys could be in the research field of
geophysics[1, 2, 3, 4, 5], communication technology[6], biosensors[7],
acoustic spin pumping [29], magnetization switching induced by sound waves
[30, 31, 32, 33, 34, 35, 36, 37], and picosecond ultrasonics [38]. VelCrys
could also be useful for further post-processing of elastic tensor and
magnetoelastic constants obtained from automated schemes [42, 59, 60, 61, 43].
Moreover, its web-based interface might be exploited in databases of material
science [39, 40].

## Declaration of competing interest

The authors declare that they have no known competing financial interests or
personal relationships that could have appeared to influence the work reported
in this paper.

## Data availability

Source files of VelCrys tool are available on the GitHub repository [49].

## Acknowledgement

This work was supported by projects e-INFRA CZ (ID:90254)" and QM4ST
(CZ.02.01.01/00/22_008/0004572) by The Ministry of Education, Youth and Sports
of the Czech Republic and also by Czech Science Foundation of the Czech
Republic by grant No. 22-35410K. P. N. acknowledges support by grant
MU-23-BG22/00168 funded by The Ministry of Universities of Spain. A.F.
acknowledges funding from the Spanish Ministry of Science and Innovation
(grants nos. PID2022-139230NB-I00, and TED2021-132074B-C32) the Diputación
Foral de Gipuzkoa (Project No. 2023-CIEN-000077-01). Research conducted in the
scope of the Transnational Common Laboratory (LTC) Aquitaine-Euskadi Network
in Green Concrete and Cement-based Materials. R. I. acknowledges financial
support from the project MAGNES funded by the Principality of Asturias
Government (grant no. AYUD/2021/51822) and from and the project RADIAFUS V,
funded by the Spanish Ministry of Science and Innovation (grant no.
PID2019-105325RB-C32).

## References

  * [1] D. Wang, C.-Y. Bai, X.-W. Li, J. Hu, Analytic formula for the group velocity and its derivatives with respect to elastic moduli for a general anisotropic medium, Geophysics 88 (4) (2023) C111–C121.  arXiv:https://pubs.geoscienceworld.org/geophysics/article-pdf/88/4/C111/5894846/geo-2022-0566.1.pdf, doi:10.1190/geo2022-0566.1.
URL https://doi.org/10.1190/geo2022-0566.1

  * [2] J. Regan, D. L. Anderson, Anisotropic models of the upper mantle, Physics of the Earth and Planetary Interiors 35 (4) (1984) 227–263.  doi:https://doi.org/10.1016/0031-9201(84)90020-7.
URL https://www.sciencedirect.com/science/article/pii/0031920184900207

  * [3] M. Cara, J. Lévêque, Oriented olivine crystals in the upper mantle: a test from the inversion of multimode surface-wave data, Physics of the Earth and Planetary Interiors 47 (1987) 246–252.  doi:https://doi.org/10.1016/0031-9201(87)90080-X.
URL https://www.sciencedirect.com/science/article/pii/003192018790080X

  * [4] L. Thomsen, Weak elastic anisotropy, GEOPHYSICS 51 (10) (1986) 1954–1966.  arXiv:https://doi.org/10.1190/1.1442051, doi:10.1190/1.1442051.
URL https://doi.org/10.1190/1.1442051

  * [5] M. Schimmel, E. Stutzmann, S. Ventosa, Measuring group velocity in seismic noise correlation studies based on phase coherence and resampling strategies, IEEE Transactions on Geoscience and Remote Sensing (55) (2017) 1928–1935.  doi:10.1109/TGRS.2016.2631445.
  * [6] H. H. Yang Y, Dejous C, Trends and applications of surface and bulk acoustic wave devices: A review, Micromachines (Basel) (14(1):43) (2022).  doi:doi:10.3390/mi14010043.
  * [7] Y. Huang, P. K. Das, V. R. Bhethanabotla, Surface acoustic waves in biosensing applications, Sensors and Actuators Reports 3 (2021) 100041.  doi:https://doi.org/10.1016/j.snr.2021.100041.
URL https://www.sciencedirect.com/science/article/pii/S2666053921000175

  * [8] L. D. Landau, E. M. Lifshitz, Theory of elasticity / by L. D. Landau and E. M. Lifshitz ; translated from the Russian by J. B. Sykes and W. H. Reid, Pergamon London, 1959.
  * [9] J. R. Neighbours, G. E. Schacher, Determination of Elastic Constants from Sound-Velocity Measurements in Crystals of General Symmetry, Journal of Applied Physics 38 (13) (1967) 5366–5375.  arXiv:https://pubs.aip.org/aip/jap/article-pdf/38/13/5366/10553733/5366\\_1\\_online.pdf, doi:10.1063/1.1709328.
URL https://doi.org/10.1063/1.1709328

  * [10] V. Červený, Seismic Rays and Ray Intensities in Inhomogeneous Anisotropic Media, Geophysical Journal International 29 (1) (1972) 1–13.  arXiv:https://academic.oup.com/gji/article-pdf/29/1/1/1734163/29-1-1.pdf, doi:10.1111/j.1365-246X.1972.tb06147.x.
URL https://doi.org/10.1111/j.1365-246X.1972.tb06147.x

  * [11] S. Crampin, A review of wave motion in anisotropic and cracked elastic-media, Wave Motion 3 (4) (1981) 343–391.  doi:https://doi.org/10.1016/0165-2125(81)90026-3.
URL https://www.sciencedirect.com/science/article/pii/0165212581900263

  * [12] B. Zhou, S. Greenhalgh, On the computation of elastic wave group velocities for a general anisotropic medium, Journal of Geophysics and Engineering 1 (3) (2004) 205.  doi:10.1088/1742-2132/1/3/005.
URL https://dx.doi.org/10.1088/1742-2132/1/3/005

  * [13] J. Rouchy, E. du Tremolet de Lacheisserie, Elastic properties of ferromagnetic single crystal with cubic symmetry, Zeitschrift für Physik B Condensed Matter 36 (1) (1979) 67–80.  doi:10.1007/BF01333955.
URL https://doi.org/10.1007/BF01333955

  * [14] J. Rouchy, E. T. de Lacheisserie, J. Genna, A. Waintal, A 6 T ultrasonic study of the magnetoelastic coupling in ferromagnetic Co—Pt alloy, Journal of Magnetism and Magnetic Materials 21 (1) (1980) 69–79.  doi:https://doi.org/10.1016/0304-8853(80)90489-8.
URL https://www.sciencedirect.com/science/article/pii/0304885380904898

  * [15] E. D. T. de Lacheisserie, Magnetostriction: Theory and Application of Magnetoelasticity, CRC Press, Boca Raton, FL, 1993.
  * [16] A. Del Moral, Handbook of magnetostriction and magnetostrictive materials, Del Moral Zaragoza, Zaragoza, 2008.
  * [17] G. Simon, Das verhalten von ebenen ultraschallwellen in homogen magnetisierten einkristallen, Zeitschrift für Naturforschung A 13 (2) (1958) 84–89.  doi:doi:10.1515/zna-1958-0204.
URL https://doi.org/10.1515/zna-1958-0204

  * [18] P. Nieves, J. Tranchida, S. Nikolov, A. Fraile, D. Legut, Atomistic simulations of magnetoelastic effects on sound velocity, Phys. Rev. B 105 (2022) 134430.  doi:10.1103/PhysRevB.105.134430.
URL https://link.aps.org/doi/10.1103/PhysRevB.105.134430

  * [19] W. P. Mason, A phenomenological derivation of the first- and second-order magnetostriction and morphic effects for a nickel crystal, Phys. Rev. 82 (1951) 715–723.  doi:10.1103/PhysRev.82.715.
URL https://link.aps.org/doi/10.1103/PhysRev.82.715

  * [20] D. E. Eastman, Ultrasonic study of first-order and second-order magnetoelastic properties of yttrium iron garnet, Phys. Rev. 148 (1966) 530–542.  doi:10.1103/PhysRev.148.530.
URL https://link.aps.org/doi/10.1103/PhysRev.148.530

  * [21] E. du Tremolet de Lacheisserie, J. Rouchy, The magnetoelastic coupling in nickel, Journal of Magnetism and Magnetic Materials 28 (1) (1982) 77–87.  doi:https://doi.org/10.1016/0304-8853(82)90031-2.
URL https://www.sciencedirect.com/science/article/pii/0304885382900312

  * [22] J. Rouchy, P. Morin, E. T. de Lacheisserie, Magnetic and magnetoelastic properties of gdzn single crystals, Journal of Magnetism and Magnetic Materials 23 (1) (1981) 59–68.  doi:https://doi.org/10.1016/0304-8853(81)90068-8.
URL https://www.sciencedirect.com/science/article/pii/0304885381900688

  * [23] J. Sakurai, Ultrasonic propagation in nickel and mn-ferrite at high magnetic fields, Journal of the Physical Society of Japan 19 (3) (1964) 311–317.  arXiv:https://doi.org/10.1143/JPSJ.19.311, doi:10.1143/JPSJ.19.311.
URL https://doi.org/10.1143/JPSJ.19.311

  * [24] G. Alers, J. Neighbours, H. Sato, Dependence of sound velocity and attenuation on magnetization direction in nickel at high fields, Journal of Physics and Chemistry of Solids 9 (1) (1959) 21–27.  doi:https://doi.org/10.1016/0022-3697(59)90085-X.
URL https://www.sciencedirect.com/science/article/pii/002236975990085X

  * [25] G. Dietz, G. Wielpütz, Magnetoelastic coupling between sound waves in ferromagnetic crystals, Applied physics 11 (2) (1976) 131–136.  doi:10.1007/BF00920591.
URL https://doi.org/10.1007/BF00920591

  * [26] C. Kingner, J. Heil, B. Lüthi, Surface acoustic wave propagation in paramagnets and ferromagnets, Journal of Applied Physics 52 (3) (1981) 2270–2272.  arXiv:https://doi.org/10.1063/1.328902, doi:10.1063/1.328902.
URL https://doi.org/10.1063/1.328902

  * [27] M. Seavey, Acoustic resonance in the easy-plane weak ferromagnets alpha; fe2o3 and febo3, Solid State Communications 10 (2) (1972) 219–223.  doi:https://doi.org/10.1016/0038-1098(72)90385-7.
URL https://www.sciencedirect.com/science/article/pii/0038109872903857

  * [28] S. Rinaldi, G. Turilli, Theory of linear magnetoelastic effects, Phys. Rev. B 31 (1985) 3051–3058.  doi:10.1103/PhysRevB.31.3051.
URL https://link.aps.org/doi/10.1103/PhysRevB.31.3051

  * [29] K. Uchida, H. Adachi, T. An, T. Ota, M. Toda, B. Hillebrands, S. Maekawa, E. Saitoh, Long-range spin seebeck effect and acoustic spin pumping, Nature Materials 10 (10) (2011) 737–741.  doi:10.1038/nmat3099.
URL https://doi.org/10.1038/nmat3099

  * [30] I. S. Camara, J.-Y. Duquesne, A. Lemaître, C. Gourdon, L. Thevenard, Field-free magnetization switching by an acoustic wave, Phys. Rev. Applied 11 (2019) 014045.  doi:10.1103/PhysRevApplied.11.014045.
URL https://link.aps.org/doi/10.1103/PhysRevApplied.11.014045

  * [31] O. Kovalenko, T. Pezeril, V. V. Temnov, New concept for magnetization switching by ultrafast acoustic pulses, Phys. Rev. Lett. 110 (2013) 266602.  doi:10.1103/PhysRevLett.110.266602.
URL https://link.aps.org/doi/10.1103/PhysRevLett.110.266602

  * [32] L. Thevenard, J.-Y. Duquesne, E. Peronne, H. J. von Bardeleben, H. Jaffres, S. Ruttala, J.-M. George, A. Lemaître, C. Gourdon, Irreversible magnetization switching using surface acoustic waves, Phys. Rev. B 87 (2013) 144402.  doi:10.1103/PhysRevB.87.144402.
URL https://link.aps.org/doi/10.1103/PhysRevB.87.144402

  * [33] L. Thevenard, I. S. Camara, S. Majrab, M. Bernard, P. Rovillain, A. Lemaître, C. Gourdon, J.-Y. Duquesne, Precessional magnetization switching by a surface acoustic wave, Phys. Rev. B 93 (2016) 134430.  doi:10.1103/PhysRevB.93.134430.
URL https://link.aps.org/doi/10.1103/PhysRevB.93.134430

  * [34] V. S. Vlasov, A. M. Lomonosov, A. V. Golov, L. N. Kotov, V. Besse, A. Alekhin, D. A. Kuzmin, I. V. Bychkov, V. V. Temnov, Magnetization switching in bistable nanomagnets by picosecond pulses of surface acoustic waves, Phys. Rev. B 101 (2020) 024425.  doi:10.1103/PhysRevB.101.024425.
URL https://link.aps.org/doi/10.1103/PhysRevB.101.024425

  * [35] A. Stupakiewicz, C. S. Davies, K. Szerenos, D. Afanasiev, K. S. Rabinovich, A. V. Boris, A. Caviglia, A. V. Kimel, A. Kirilyuk, Ultrafast phononic switching of magnetization, Nature Physics 17 (4) (2021) 489–492.  doi:10.1038/s41567-020-01124-9.
URL https://doi.org/10.1038/s41567-020-01124-9

  * [36] W. Li, B. Buford, A. Jander, P. Dhagat, Acoustically assisted magnetic recording: A new paradigm in magnetic data storage, IEEE Transactions on Magnetics 50 (3) (2014) 37–40.  doi:10.1109/TMAG.2013.2285018.
  * [37] M. Strungaru, M. O. A. Ellis, S. Ruta, R. F. L. Evans, R. W. Chantrell, O. Chubykalo-Fesenko, A route to minimally dissipative switching in magnets via thz phonon pumping (2023).  arXiv:2209.04312.
  * [38] M. Mattern, A. von Reppert, S. P. Zeuschner, M. Herzog, J.-E. Pudell, M. Bargheer, Concepts and use cases for picosecond ultrasonics with x-rays, Photoacoustics 31 (2023) 100503.  doi:https://doi.org/10.1016/j.pacs.2023.100503.
URL https://www.sciencedirect.com/science/article/pii/S2213597923000563

  * [39] A. Jain, S. Ong, G. Hautier, W. Chen, W. Richards, S. Dacek, S. Cholia, D. Gunter, D. Skinner, G. Ceder, K. Persson, The materials project: A materials genome approach to accelerating materials innovation, APL Materials 1 (2013) 011002.
  * [40] S. Curtarolo, W. Setyawan, G. L. Hart, M. Jahnatek, R. V. Chepulskii, R. H. Taylor, S. Wang, J. Xue, K. Yang, O. Levy, M. J. Mehl, H. T. Stokes, D. O. Demchenko, D. Morgan, Aflow: An automatic framework for high-throughput materials discovery, Computational Materials Science 58 (2012) 218–226.  doi:https://doi.org/10.1016/j.commatsci.2012.02.005.
URL https://www.sciencedirect.com/science/article/pii/S0927025612000717

  * [41] M. Scheidgen, L. Himanen, A. N. Ladines, D. Sikter, M. Nakhaee, Ádám Fekete, T. Chang, A. Golparvar, J. A. Márquez, S. Brockhauser, S. Brückner, L. M. Ghiringhelli, F. Dietrich, D. Lehmberg, T. Denell, A. Albino, H. Näsström, S. Shabih, F. Dobener, M. Kühbach, R. Mozumder, J. F. Rudzinski, N. Daelman, J. M. Pizarro, M. Kuban, C. Salazar, P. Ondračka, H.-J. Bungartz, C. Draxl, Nomad: A distributed web-based platform for managing materials science research data, Journal of Open Source Software 8 (90) (2023) 5388.  doi:10.21105/joss.05388.
URL https://doi.org/10.21105/joss.05388

  * [42] R. Gaillac, P. Pullumbi, F.-X. Coudert, ELATE: an open-source online application for analysis and visualization of elastic tensors, Journal of Physics: Condensed Matter 28 (27) (2016) 275201.  doi:10.1088/0953-8984/28/27/275201.
URL https://doi.org/10.1088%2F0953-8984%2F28%2F27%2F275201

  * [43] P. Nieves, S. Arapan, A. P. Kądzielawa, D. Legut, MAELASviewer: An online tool to visualize magnetostriction, Sensors 20 (22) (2020).  doi:10.3390/s20226436.
URL https://www.mdpi.com/1424-8220/20/22/6436

  * [44] https://github.com/pnieves2019/MAELASwiewer.
  * [45] http://www.md-esg.eu/maelasviewer/.
  * [46] https://dash.plotly.com/.
  * [47] https://flask.palletsprojects.com/en/1.1.x/.
  * [48] https://plotly.com/python/.
  * [49] https://github.com/pnieves2019/VelCrys.
  * [50] H. Tian, T. Kempka, S. Yu, M. Ziegler, Mechanical properties of sandstones exposed to high temperature, Rock Mechanics and Rock Engineering 49 (1) (2016) 321–327.  doi:10.1007/s00603-015-0724-z.
  * [51] Q. Ben-dong, H. Jun, C. Lun-jian, Experimental research on mechanical properties of limestone and sandstone under high temperature (2009).
URL /article/id/615fae88ed73f804e4fe0df6

  * [52] https://github.com/li913/Velocity4GeneralAnisotropicMedium.
  * [53] P. Villars, K. Cenzual, R. Gladyshevskii, Handbook, De Gruyter, Berlin, Boston, 27 Nov. 2015.  doi:10.1515/9783110334678.
URL https://www.degruyter.com/view/title/304631

  * [54] A. S. Darling, Cobalt-platinum alloys, Platinum Metals Review 7 (3) (1963) 96–104.  doi:https://doi.org/10.1595/003214063X7396104.
URL https://technology.matthey.com/content/journals/10.1595/003214063X7396104

  * [55] Q. Dehaine, L. T. Tijsseling, H. J. Glass, T. Törmänen, A. R. Butcher, Geometallurgy of cobalt ores: A review, Minerals Engineering 160 (2021) 106656.  doi:https://doi.org/10.1016/j.mineng.2020.106656.
URL https://www.sciencedirect.com/science/article/pii/S0892687520304763

  * [56] V. G. Bar’yakhtar, A. G. Danilevich, Magnetoelastic waves in ferromagnets in the vicinity of lattice structural phase transitions, Ukrainian Journal of Physics 63 (9) (2018) 836.  doi:10.15407/ujpe63.9.836.
URL https://ujp.bitp.kiev.ua/index.php/ujp/article/view/2018150

  * [57] R. C. O’Handley, Modern magnetic materials, Wiley, 2000.
  * [58] P. Nieves, S. Arapan, S. Zhang, A. Kądzielawa, R. Zhang, D. Legut, Maelas: Magneto-elastic properties calculation via computational high-throughput approach, Computer Physics Communications 264 (2021) 107964.  doi:https://doi.org/10.1016/j.cpc.2021.107964.
URL https://www.sciencedirect.com/science/article/pii/S0010465521000801

  * [59] M. de Jong, W. Chen, T. Angsten, A. Jain, R. Notestine, A. Gamst, M. Sluiter, C. Krishna Ande, S. van der Zwaag, J. J. Plata, C. Toher, S. Curtarolo, G. Ceder, K. A. Persson, M. Asta, Charting the complete elastic properties of inorganic crystalline compounds, Scientific Data 2 (1) (2015) 150009.  doi:10.1038/sdata.2015.9.
URL https://doi.org/10.1038/sdata.2015.9

  * [60] S. Zhang, R. Zhang, AELAS: Automatic ELAStic property derivations via high-throughput first-principles computation, Computer Physics Communications 220 (2017) 403 – 416.  doi:10.1016/j.cpc.2017.07.020.
URL http://www.sciencedirect.com/science/article/pii/S0010465517302400

  * [61] https://github.com/pnieves2019/MAELAS.


Holographic Ricci DE as running vacuum with nonlinear interactions

Paxy George

Department of Physics,

Cochin University of Science and Technology, Kochi-22, India.

###### Abstract

The holographic Ricci dark energy can be treated as a running vacuum due to
its analogy in the energy density, which is a combination of H and $\dot{H}$,
the model can predict either eternal acceleration or eternal deceleration. In
the earlier works, we have shown that the presence of additive constant in the
energy density or by considering possible interaction between dark sectors
through a phenomenological term, the model can predict a transition from a
prior decelerated to a late accelerated epoch. This paper analyses the cosmic
evolution of holographic Ricci dark energy as running vacuum with a nonlinear
interaction between dark sectors in a flat FLRW universe. We consider three
possible nonlinear interaction forms which give analytically feasible
solutions. We have constrained the model using the Type1a
Supernova(Pantheon)+CMB(Planck 2018)+BAO(SDSS) data and evaluated the best-
estimated values of all the model parameters. We have analyzed the evolution
of the Hubble parameter and deceleration parameter of all three cases. We
perform state finder analysis of the model, which implies the quintessence
nature of the model and found that it is distinguishably different from the
standard $\Lambda$CDM model. The dynamical system analysis of all three cases
confirms the evolution of the universe from an unstable prior matter-dominated
era to a stable end de Sitter phase.

##  1 Introduction

The observational evidences from Type1a supernova (SN1a)[1, 2], large scale
structure (LSS)[3, 4], and cosmic microwave background (CMB)[5, 6] have
converged towards one of the focus of attention of modern cosmology- the
accelerated expansion of the universe. An exotic component called dark energy
(DE) with negative pressure was considered as the mysterious driving force
responsible for the accelerated expansion of the universe. The most prominent
and simplest model that explains the cosmic acceleration is the concordance
model (or $\Lambda$CDM model) that considers the cosmological constant as the
dark energy. The $\Lambda$CDM model describes a universe filled with cosmic
components such as cold dark matter (CDM) and the cosmological
constant($\Lambda$); the model accurately fits the current observational data.
This model, however, was plagued with two major issues such as (a) coincidence
problem, which is the fact that the present value of the matter density and
the vacuum energy (VE) density is of the same order, even though the former is
a rapidly decreasing function of time while the latter is just stationary, and
(b) fine-tuning problem for which the observed value of the VE density is 120
orders of magnitude below the value computed using the quantum field theory
(QFT). These drawbacks can be alleviated by assuming that CC can be a
dynamical quantity. The main objective of these classes of models is to
consider that the primary driving force that accelerates the universe is
dynamical VE and, therefore, a time variable CC. The possibility of a
dynamical DE is not excluded by current observations[7, 8, 9]. The recent
Planck observations indicate that the DE equation of state (EOS) parameter
$\omega\equiv P_{DE}/\rho_{DE}$ is close to $-1$ i.e., $\omega=-1.03\pm
0.03$[10]. Nowadays, models with time-dependent VE density $\Lambda(t)$, also
known as running vacuum energy (RVE) density attains much attention. The
primary feature in this ideology is that the equation of state parameter of
the VE is strictly equal to -1[11, 12].

The idea of RVE in an expanding universe was first proposed in the
literature[13, 14]. It is pretty natural to expect that the VE is a dynamical
quantity as the universe expands and thereby sensitive to time-evolving
functions such as the Hubble rate H = H(t). In these models, without invoking
the need for scalars, a phenomenologically viable description of the dynamical
nature of the VE can be achieved. By combining the phenomenological approach
with the renormalization group (RG) techniques of quantum field theory (QFT)
in curved spacetime, a large class of dynamical $\Lambda(H,\dot{H})$ models
described by an even power series expansion of the Hubble rate and its
derivatives were proposed by J.Sola and Shapiro[15, 16, 12]. However, for the
present universe, only the terms $\dot{H}$ and $H^{2}$ can be relevant,
together with an additive constant term. Recent studies remarkably shows that
these models that include the additive term plus one or both of the dynamical
components $\dot{H}$ and $H^{2}$ appeared to be more favoured than the
$\Lambda$CDM[17, 18, 19, 20, 21, 22]. The authors in the reference[23] shows
that these models can provide an excellent fit to the main cosmological data
SN1a+H(z)+BAO+LSS+CMB+BBN whereas $\Lambda$CDM model that supports a rigid
cosmological constant is currently disfavored at 4.2$\sigma$.
Phenomenologically, this model has a formal analogy to the holographic Ricci
dark energy (HRDE), where a combination of $\dot{H}$ and $H^{2}$ is present on
its energy density.

The holographic DE scenario is another fascinating approach to explain the
late-time accelerated expansion of the universe within the context of quantum
gravity. Inspired by applying the holographic principle to the universe,
Li[24, 25] first introduces holographic DE. The holographic DE density is
inversely proportional to the square of an appropriate length of the system
that characterizes the size of the universe. One of the natural choices of
this length is the inverse of the Hubble rate. However, this choice does not
induce acceleration in a homogeneous and isotropic universe[26, 27]. Another
choice for the length was suggested by Gao et al.[28, 29] in which the IR cut-
off was taken to be proportional to the Ricci scalar curvature, leading to the
introduction of the new holographic DE, the holographic Ricci DE. This model
can derive an accelerated universe that solves the cosmological constant
problem and fairly fits the observational data well. This model has been
studied widely in many literature[30, 31, 32, 33, 34]. The analogy between the
energy density of this model with the running vacuum model motivates us to
consider HRDE as a running vacuum.

In the previous works[35, 36], we have shown that considering HRDE as a
running vacuum can satisfy a conventional evolutionary status with a prior
deceleration epoch followed by an accelerated epoch. In the reference[35], we
have shown that a bare cosmological constant in the RVE density is essential
for the transition from a prior decelerated epoch to a late accelerated epoch.
Later, in the work[36], we have shown that, by considering the interaction
between the RVE density and the dark matter through the phenomenological term
$Q=3bH\rho_{m}$, the model can predict a transition from the prior decelerated
to a late accelerated phase even without an additive constant in the VE
density. The possibility to have a nonlinear interaction (NLI) between dark
sectors have been discussed in many literature[37, 38, 39, 40]. It is of great
interest to study the effects of nonlinear interaction between DE and DM on
the evolution of the universe, especially the influence of different forms of
nonlinear interaction terms on cosmological parameters. The authors in the
reference[41] have shown that nonlinear interaction between dark sectors is
more compatible with recent observational data than linear interaction
choices. In this paper, we have considered three different forms for the
nonlinear interaction terms, $Q$, which allow us to find analytically feasible
solutions. The paper is organized as follows. Section 2 discusses a brief
introduction of the RVE model and HRDE model. This section also consists of
our analysis of the evolution of the Hubble parameter corresponding to three
nonlinear interaction forms. In section 3, we study the evolution of the
deceleration parameter, and these models are compared with the standard
$\Lambda$CDM model by performing the state finder analysis of three cases.
Finally, the stability of all the nonlinear interaction terms was checked by
performing dynamical system analysis that predicts an early matter-dominated
and late accelerated epoch.

##  2 The model : Holographic Ricci DE as running vacuum with nonlinear
interactions

The running vacuum model proposed in the literature[42, 16, 15] is based on
the possibility that quantum effects in curved space time[43] can be
responsible for the renormalization group(RG) running of the VE, which is
mainly aimed to alleviate the cosmological constant problem[44]. The
corresponding VE density, $\rho_{\Lambda}$, should be time-dependent
quantities in cosmology. The general form of RVE density for the current
universe can be represented as a function of Hubble parameter and its
derivatives[45],

$\rho_{\Lambda}(H;\nu,\alpha)=\frac{3}{8\pi G}\left(c_{0}+\nu
H^{2}+\frac{2}{3}\alpha\dot{H}+\mathcal{O}(H^{4})\right),$ (1)

where, $c_{0}$ is an integration constant, $\nu$ and $\alpha$ are
dimensionless coefficients. The values of $\nu$ and $\alpha$ can be related to
the $\beta$-functions of the running vacuum and hence are assumed to be
naturally small in the order of $10^{-3}$. It is possible to recover
$\Lambda$CDM model for $\nu=\alpha=0$ but here we consider them as free
parameters hence we shall determine its values by fitting the model to
observations.The term $\mathcal{O}(H^{4})$ represent higher powers of the
Hubble rate which are relevant to describe inflation[46, 47, 48, 49], but
these terms has no significance for present universe and hence ignored. From
the expression(1), it is clear that the VE density is slowly evolving with the
Hubble function $H$ and its time derivative as the universe expands.

The holographic principle is believed to be a fundamental principle of quantum
gravity. Motivated by holographic principle, Cohen et al.[50] suggest that the
quantum zero-point energy of a system with size $L$ should not exceed the mass
of a black hole with the same size $L^{2}\rho_{\Lambda}\leq LM_{p}^{2}$. This
relation relates the UV cut-off of a system to its IR cut-off. The largest $L$
saturating this inequality will give rise to a holographic energy density,

$\rho_{\Lambda}=3c^{2}M_{p}^{2}L^{2},$ (2)

where c is a dimensionless constant and $M_{p}^{2}=\frac{1}{8\pi G}$ is the
reduced Planck mass. Gao et al[51] have proposed a holographic DE which is
proportional to the inverse of the Ricci scalar curvature, which is aimed to
avoid causality problems and is phenomenologically viable. The holographic DE
proportional to the Ricci scalar for a flat FLRW universe is,

$\rho_{de}(H,\dot{H})=3\beta M_{p}^{2}(\dot{H}+2H^{2}),$ (3)

where $\beta$ is a dimensionless coefficient determined by observations. The
major drawback of this model is that it will not allow a transition from a
decelerated epoch to a late accelerating epoch. To alleviate this issue, a
bare cosmological constant has been added to the initial form of the HRDE
density, and then it resembles the general form of RVE density as given in
equation(1). This motivates us to consider HRDE as a running vacuum in our
work. The model parameter $\beta$ can be related to the terms of the original
RG running parameter as $\beta\sim\frac{\nu}{2}\sim\frac{2\alpha}{3}$[36]. In
the work[36] we have shown that instead of considering an additive constant in
the energy density of HRDE, it is possible to attain a transition from
decelerated to accelerated epoch by accounting for a possible interaction
between dark sectors. The evolution equation for a spatially flat Friedmann
universe is,

$3H^{2}=\rho_{m}+\rho_{de},$ (4)

where $H=\frac{\dot{a}}{a},$ is the Hubble parameter, a is the scale factor,
$\rho_{m}$ is the matter density and $\rho_{de}$ is the DE density. Here we
also suppose that DM is pressureless and DE satisfies the equation of
state,$\omega=-1$ which gives $P_{de}=-\rho_{de}$. When a possible interaction
between dark sectors are taken into account, the conservation equation takes
the form,

$\begin{split}\dot{\rho_{m}}+3H(\rho_{m}+P_{m})=Q,\\\
\dot{\rho_{de}}+3H(\rho_{de}+P_{de})=-Q,\end{split}$ (5)

where Q is the interaction term which gives the rate of energy exchange
between dark sectors. The sign of Q determines the direction of flow of
energy. If Q is positive the DE decays to DM while if Q is negative DM decays
to DE. The interaction terms are introduced by phenomenological considerations
to avoid the coincidence problem[38]. The interaction term proportional to the
DE density or to the DM density or to a linear combination of both densities
was well studied in literature[52, 53, 54, 55]. But in some literature[39, 40,
56, 57, 37, 58] the non linear interaction i.e., an interaction proportional
to the product of DM density and DE energy density were studied and found to
be favored observationally than linear interaction between dark sectors[41].
The general form for the nonlinear interaction can be written as[39],

$Q=3bH\rho^{m+s}\rho_{m}^{n}\rho_{de}^{-s-n},$ (6)

where $b$ is the coupling constant, $\rho=\rho_{m}+\rho_{de}$ and the powers
$m,n$ and $s$ specify interaction. We have analyzed the case $Q=3bH\rho_{m}$
for $(m,n,s)=(1,1,-1)$ in reference[36]. In the present work we consider three
possible nonlinear forms for Q i.e., $Q=3bH\frac{\rho_{m}^{2}}{\rho}$ for
$(m,n,s)=(1,2,-2)$, $Q=3bH\frac{\rho_{de}^{2}}{\rho}$ for $(m,n,s)=(1,0,-2)$
and $Q=3bH\frac{\rho_{m}\rho_{de}}{\rho}$ for $(m,n,s)=(1,1,-2).$ These
particular combinations of the parameters$(m,n,s)$ gives rise to analytically
solvable models with nonlinear interaction terms.

###  2.1 Case 1: $Q=3bH\frac{\rho_{m}^{2}}{\rho}$

The model with this interaction term $Q=3bH\frac{\rho_{m}^{2}}{\rho}$ is named
as nonlinear interaction case1(NLI1). The conservation equation(5) for this
case takes the form,

$\begin{split}\frac{d\rho_{de}}{dx}=-b\frac{\rho_{m}^{2}}{H^{2}},\\\
\frac{d\rho_{m}}{dx}=\frac{b\rho_{m}^{2}}{H^{2}}-3\rho_{m}.\end{split}$ (7)

Here we have changed the variable from t to $x=\ln a.$ Combining the equations
(4),(7) and (3) we get

$H^{2}\left(\frac{d^{2}H^{2}}{dx^{2}}+3\frac{dH^{2}}{dx}\right)+9bH^{4}-6bH^{2}\rho_{de}+b\rho_{de}^{2}=0.$
(8)

The solution of the above differential equation give the Hubble parameter as

$H^{2}=c_{2}e^{-Ax}\cosh[B(x-4c_{1})]^{F},$ (9)

where $c_{1}$ and $c_{2}$ are integration constants, A, B and F are given by

$\begin{split}A=\frac{6(1+3b\beta(2\beta-1))}{4+9b\beta^{2}},\\\
B=\frac{3}{2}[\sqrt{1-2b(2-5\beta+2\beta^{2})}],\\\
F=\frac{4}{4+9b\beta^{2}}.\end{split}$ (10)

The values of the integration constants $c_{1}$ and $c_{2}$ can be obtained by
using the conditions

$H|_{x=0}=H_{0},\,\,\,\frac{dH^{2}}{dx}|_{x=0}=2H_{0}^{2}\left(\frac{\Omega_{de0}}{\beta}-2\right),$
(11)

where $\Omega_{de0}$ is the current value of the RVE density. The values of
the coefficients are found to be,

$\begin{split}c_{1}=\frac{1}{4B}\tanh^{-1}\left(\frac{-2}{BF}\left(\frac{\Omega_{de0}}{\beta}-2\right)-\frac{A}{BF}\right),\\\
c_{2}=\frac{H_{0}^{2}}{\cosh(4c_{1}B)^{F}}.\end{split}$ (12)

###  2.2 Case 2: $Q=3bH\frac{\rho_{de}^{2}}{\rho}$

We consider the interaction $Q=3bH\frac{\rho_{de}^{2}}{\rho}$, nonlinear
interaction case2 (NLI2) which reduce the conservation equation(5) as,

$\begin{split}\frac{d\rho_{de}}{dx}=-b\frac{\rho_{de}^{2}}{H^{2}},\\\
\frac{d\rho_{m}}{dx}=\frac{b\rho_{de}^{2}}{H^{2}}-3\rho_{m}.\end{split}$ (13)

Using equations(4),(13) and (3) we get a second order differential equation,

$H^{2}\left(\frac{d^{2}H^{2}}{dx^{2}}+3\frac{dH^{2}}{dx}\right)+9b\beta^{2}\left(4H^{4}+2H^{2}\frac{dH^{2}}{dx}+\frac{1}{4}(\frac{dH^{2}}{dx})^{2}\right)=0.$
(14)

The solution of the above equation is the Hubble parameter of this model NLI2

$H^{2}=c_{2}e^{-Ax}\cosh[B(x-4c_{1})]^{F}.$ (15)

The terms $A,B$ and $F$ are given by,

$\begin{split}A=\frac{6(1+6b\beta^{2})}{4+9b\beta^{2}},\\\
B=\frac{3}{2}\sqrt{1-4b\beta^{2}},\\\ F=\frac{4}{4+9b\beta^{2}}.\end{split}$
(16)

Here $c_{1}$ and $c_{2}$ are the integration constants whose values can be
obtained from the conditions(11),

$\begin{split}c_{1}=\frac{1}{4B}\tanh^{-1}\left(\frac{-2}{BF}\left(\frac{\Omega_{de0}}{\beta}-2\right)-\frac{A}{BF}\right),\\\
c_{2}=\frac{H_{0}^{2}}{\cosh(4c_{1}B)^{F}}.\end{split}$ (17)

###  2.3 Case 3: $Q=3bH\frac{\rho_{m}\rho_{de}}{\rho}$

The interaction considered in this case is
$Q=3bH\frac{\rho_{m}\rho_{de}}{\rho}$ which can be conveniently represented as
nonlinear interaction case3 (NLI3). The conservation equation(5) using this
type of interaction can be written as

$\begin{split}\frac{d\rho_{de}}{dx}=\frac{-b\rho_{de}\rho_{m}}{H^{2}},\\\
\frac{d\rho_{m}}{dx}=\frac{b\rho_{de}\rho_{m}}{H^{2}}-3\rho_{m}.\end{split}$
(18)

Combining the above conservation equation with (4) and (3) give a second order
differential equation,

$H^{2}\left(\frac{d^{2}H^{2}}{dx^{2}}+3\frac{dH^{2}}{dx}+9b\beta\left(2H^{2}+\frac{1}{2}\frac{dH^{2}}{dx}\right)\right)-9b\beta^{2}\left(4H^{4}+2H^{2}\frac{dH^{2}}{dx}+\frac{1}{4}\left(\frac{dH}{dx}\right)^{2}\right).$
(19)

The solution of the above differential equation gives Hubble parameter for
this case,

$H^{2}=c_{2}e^{-Ax}\cosh[B(x+4c_{1})]^{F}.$ (20)

The terms A, B and F in the Hubble parameter (20) are,

$\begin{split}A=\frac{3(3b\beta(4\beta-1)-2)}{9b\beta^{2}-4},\\\
B=\frac{3}{4}\sqrt{4+9b^{2}\beta^{2}+4b\beta(4\beta-5)},\\\
F=\frac{-4}{9b\beta^{2}-4}.\end{split}$ (21)

The integration $c_{1}$ and $c_{2}$ can be evaluated by using the conditions
(11) as,

$\begin{split}c_{1}=\frac{1}{4B}\tanh^{-1}\left(\frac{2}{BF}\left(\frac{\Omega_{de0}}{\beta}-2\right)+\frac{A}{BF}\right),\\\
c_{2}=\frac{H_{0}^{2}}{\cosh(4c_{1}B)^{F}}.\end{split}$ (22)

The values of model parameter $\beta$, interaction term $b$, $\Omega_{de0}$
and $H_{0}$ can be constrained by using $\chi^{2}$ minimization technique with
the following data sets. SN1a(supernova type 1a): the Pantheon sample,
composed of 1048 data points in the redshift range $0.01\leq z\leq 2.3$,
obtained from the compilation of 279 SN1a discovered by the Pan-STARRS1
medium-deep survey, the distance estimates from the Sloan Digital Sky
Survey(SDSS), Supernova Legacy Survey(SNLS) and from various low redshift and
Hubble Space Telescope(HST) samples[59]; CMB(cosmic microwave background)
data: CMB shift parameter from Planck 2018[60], $\mathcal{R}=1.7502\pm
0.0046$, the redshift at last scattering $z_{\ast}=1089.92$; BAO(baryon
acoustic oscillation)data: The acoustic parameter from SDSS(Solan Digital Sky
Survey), $\mathcal{A}=0.484\pm 0.016$, the redshift at signature of the peak
oscillation measured $z=0.35$[61, 62].

We have the distance modulus as a function of the model parameters,

$\mu_{t}(\beta,b,\Omega_{de0},H_{0},z_{i})=m-M\\\
=5\log_{10}\left[\frac{d_{L}(\beta,b,\Omega_{de0},H_{0},z_{i})}{Mpc}\right]+25,$
(23)

where $m$ and $M$ are the apparent and absolute magnitudes of the
supernovae(SN1a) respectively, $z_{i}$ is the red shift of the supernova and
the luminosity distance $d_{L}$ is defined as,

$d_{L}(\beta,b,\Omega_{de0},H_{0},z_{i})=c(1+z_{i})\int_{0}^{z_{i}}\frac{dz}{H(\beta,b,\Omega_{de0},H_{0},z)},$
(24)

where $H(\beta,b,H_{0},z)$ is the Hubble parameter as a function of the model
parameters and $c$ is the speed of light. Now the $\chi^{2}$ function which
include theoretical and observational magnitude can be written as,

$\chi^{2}(\beta,b,\Omega_{de0},H_{0})=\sum_{i=1}^{n}\frac{[\mu_{t}(\beta,b,\Omega_{de0},H_{0},z_{i})-\mu_{i}]^{2}}{\sigma_{i}^{2}},$
(25)

where $\mu_{i}$ is the observational distance moduli for the $i^{th}$
Supernova, $\sigma_{i}^{2}$ is the standard variation of the observation and
$n=1048,$ the total number of data points.

For obtaining the $\chi^{2}$ function, we also used The BAO signal has been
directly detected by the SDSS survey at a scale $\sim$100MPc. The BAO peak
parameter value was first proposed by D. J. Eisenstein, et al[63] and is
defined as,

$\mathcal{A}=\frac{\sqrt{\Omega_{m}}}{h(z_{1})^{\frac{1}{3}}}\left(\frac{1}{z_{1}}\int_{0}^{z_{1}}\frac{dz}{h(z)}\right)^{\frac{2}{3}},$
(26)

Here h(z) is the Hubble parameter, $z_{1}=0.35$ is the red shift of the SDSS
sample[64]. Using SDSS data from luminous red galaxies survey the value of the
parameter $\mathcal{A}$(for flat universe) is given by $\mathcal{A}=0.484\pm
0.016$[62]. The $\chi^{2}$ function for the BAO measurement takes the form,

$\chi^{2}_{BAO}=\frac{(\mathcal{A}-0.484)^{2}}{(0.016)^{2}}.$ (27)

The CMB shift parameter is the first peak of CMB power spectrum[65] can be
written as,

$\mathcal{R}=\sqrt{\Omega_{m}}\int_{0}^{z_{2}}\frac{dz}{h(z)},$ (28)

Here $z_{2}$ is the redshift at the last scattering surface.

From the Planck 2018 data, $z_{2}=1089.92$. At this redshift, $z_{2}$, the
value of shift parameter would be $\mathcal{R}=1.7502\pm 0.0046$[60]. The
$\chi^{2}$ function for the CMB measurement can be written as,

$\chi^{2}_{CMB}=\frac{(\mathcal{R}-1.7502)^{2}}{(0.0046)^{2}}.$ (29)

Considering three cosmological data sets together, i.e., (SN1a+BAO+CMB), the
total $\chi^{2}$ function is then given by,

$\chi_{total}^{2}(\beta,b)=\chi_{SN1a}^{2}+\chi_{CMB}^{2}+\chi_{BAO}^{2}.$
(30)

The best fit values for the model parameters corresponding to the 68.3$\%$
confidence interval for the three cases are given in Table 1. Along with this,
the best fit values of the absolute magnitude, M, of SN1a can be evaluated.
From the table, it is clear that values of the model parameter $\beta$ for
NLI2 and NLI3 are the same i.e.,$\beta\sim 0.456$ whereas NLI1 shows a slight
change in its value. The value of the interaction term is $\sim$ 0.003 for
NLI1 and NLI3, whereas for NLI2, its value is 0.001. For all the three cases,
the best fit values of $\Omega_{de0}$, $H_{0}$, and M are found to be 0.72,
68.7, 19.39, respectively. The chi-square degrees of freedom,
$\chi^{2}_{d.o.f}=\frac{\chi^{2}_{min}}{m-n}$, where $m$ is the number of data
points, and $n$ is the number of the model parameter, is found to be close to
one.

Table 1: The parameter values for the SN1a+CMB+BAO data. $Paramter$ | $NLI1$ | $NLI2$ | $NLI3$
---|---|---|---
$\beta$ | $0.4632_{-0.003}^{+0.003}$ | $0.4712_{-0.002}^{+0.004}$ | $0.4656_{-0.002}^{+0.002}$
$b$ | $0.003_{-0.006}^{+0.006}$ | $0.001_{-0.006}^{+0.005}$ | $0.003_{-0.004}^{+0.004}$
$\Omega_{de0}$ | $0.721_{-0.016}^{+0.015}$ | $0.719_{-0.016}^{+0.016}$ | $0.720_{-0.014}^{+0.013}$
$H_{0}$ | $68.74_{-0.42}^{+0.45}$ | $68.77_{-0.46}^{+0.47}$ | $68.76_{-0.46}^{+0.46}$
M | 19.39 | 19.39 | 19.39
$\chi^{2}_{min}$ | 1024.27 | 1024.23 | 1024.24

The evolution of the Hubble parameter with scale factor for the three cases
using the best fit values of the model parameters are shown in Figure (1), (2)
and (3). From the plot, it is clear that the density becomes infinity at the
origin, which predicts the presence of a big bang at the origin. This was
again confirmed by the asymptotic behavior of the Hubble parameter of three
cases given in equations (9), (15) and (20),i.e., when $a\to 0$ tends H $\to$
$\infty$.

Figure 1: The plot shows the evolution of Hubble parameter for the case NLI1.
Figure 2: The plot shows the evolution of Hubble parameter for the case NLI2.
Figure 3: The plot shows the evolution of Hubble parameter for the case NLI3.

##  3 Background analysis

###  3.1 Deceleration parameter

The decelerating/accelerating nature of the universe can be verified by
studying the evolution of the deceleration parameter, which can be defined as

$q=-1-\frac{\dot{H}}{H^{2}}.$ (31)

If the decelerating parameter is negative then the universe is in the
accelerating phase and vice versa. By substituting the values of the Hubble
parameter and its derivatives for three interaction cases we get,

$\begin{split}q_{(NLI1)}=-1-\frac{-Ac_{2}e^{-Ax}\cosh[B(-4c_{1}+x)]^{F}+Bc_{2}e^{-Ax}F\cosh[B(-4c_{1}+x)]^{-1+F}Sinh[B(-4c_{1}+x)]}{c_{2}e^{-Ax}\cosh[B(x-4c_{1})]^{F}}.\\\
q_{(NLI2)}=-1-\frac{-Ac_{2}e^{-Ax}\cosh[B(-4c_{1}+x)]^{F}+Bc_{2}e^{-Ax}F\cosh[B(-4c_{1}+x)]^{-1+F}Sinh[B(-4c_{1}+x)]}{c_{2}e^{-Ax}\cosh[B(x-4c_{1})]^{F}}.\\\
q_{(NLI3)}=-1-\frac{-Ac_{2}e^{-Ax}\cosh[B(4c_{1}+x)]^{F}+Bc_{2}e^{-Ax}F\cosh[B(4c_{1}+x)]^{-1+F}Sinh[B(4c_{1}+x)]}{c_{2}e^{-Ax}\cosh[B(x+4c_{1})]^{F}}.\end{split}$
(32)

The variation of the deceleration parameter $q$ with the red shift(z) for
three cases are shown in figure (4), (5) and (6). From the plots, the
transition red shift for the three cases are found to be $z_{T}\sim 0.72$. The
transition red shift of $\Lambda$CDM model in combination with SN1a+CMB data
is reported in the range $z_{T}=0.45-0.73$[66]. The transition red shift
predicted by the three cases of present model are found to be in agreement
with the upper limit of the corresponding $\Lambda$CDM range. The present
value of the deceleration parameter of the models corresponding to $z=0$ are
found to be $q_{0}\sim-0.577_{-0.037}^{+0.057}$ which is close to the WMAP
value of the deceleration parameter $q_{0}\sim-0.60$[67] and with the
Planck(2015) value $q_{0}\sim-0.55$[68]. From the q-z plot, it is clear that
the q-z curve stabilizes at $\sim-1$.

Figure 4: The q-z plot for the case NLI1. Figure 5: The q-z plot for the case
NLI2 Figure 6: The q-z plot for the case NLI3

###  3.2 Statefinder analysis

Statefinder analysis is a geometric diagnostic technique for contrasting
various DE models with the standard $\Lambda$CDM model introduced by Sahni et
al.[69]. This analysis allows one to characterize the properties of DE in a
model-independent manner. Usually, Hubble parameter, scale factor,
deceleration parameter, etc., are used for describing a model; it is better to
use quantities involving higher derivatives of H or scale factor. The state
finder parameter pair {r,s}depends on the third-order derivative of the scale
factor. The trajectories in the r-s plane for various DE models are generated
to bring about the discrimination between these models. For the concordance
$\Lambda$CDM model, the state finder parameter is a fixed point
$\\{r,s\\}=\\{1,0\\}$ but for quintessence and phantom models the trajectories
lie in the region $s>0;r<1$ whereas for Chaplygin gas models trajectories lie
in the region $s<0;r>1$. Thus the difference in the trajectory of a given
model from the standard $\Lambda$CDM model discriminates the model up to a
great extent. The state finder parameter $\\{r,s\\}$ can be defined as,

$\begin{split}\\\
r=\frac{\dddot{a}}{aH^{3}}=\frac{1}{2h^{2}}\frac{d^{2}h^{2}}{dx^{2}}+\frac{3}{2h^{2}}\frac{dh^{2}}{dx}+1,\\\
s=\frac{r-1}{3(q-\frac{1}{2})}=-\frac{\frac{1}{2h^{2}}\frac{d^{2}h^{2}}{dx^{2}}+\frac{3}{2h^{2}}\frac{dh^{2}}{dx}}{\frac{3}{2h^{2}}\frac{dh^{2}}{dx}+\frac{9}{2}}.\end{split}$
(33)

By substituting the value of Hubble parameter and its derivatives in the above
equation(33) for the three interaction cases, the evolution of the $\\{r,s\\}$
parameter was plotted and is given in the figure(7), (8) and (9).

Figure 7: The r-s plot for NLI1 Figure 8: The r-s plot for NLI2 Figure 9: The
r-s plot for NLI3

The value of the statefinder parameter corresponding to the current epoch for
the three cases is found to be

$\displaystyle(r_{0},s_{0})_{NLI1}=(0.999,0.003),\qquad(r_{0},s_{0})_{NLI2}=(0.997,0.007),$
(34) $\displaystyle(r_{0},s_{0})_{NLI3}=(0.997,0.008).$

This indicates that all three cases are distinguishably different from the
$\Lambda$CDM model. The present values of$\\{r,s\\}$ parameter of NLI2 and
NLI3 are almost the same, whereas there is a slight difference in its value
for the case of NLI1. The evolution of the trajectories of the three cases
shows that the state finder parameter lies in the region corresponding to
$r<1$ and $s>0$ until it reaches the fixed $\Lambda$CDM point resembles the
quintessence nature of DE[70]. The trajectory of NLI1 and NLI3 is approaching
the $\Lambda$CDM fixed point in the future, whereas NLI2 is diverging away
from the $\Lambda$CDM point.

##  4 Phase Space Analysis

We study the asymptotic behavior of the present model for the three nonlinear
interaction cases separately by using dynamical system tools. We define a new
set of dimensionless variables,

$u=\frac{\rho_{m}}{3H^{2}}\hskip 36.135ptv=\frac{\rho_{de}}{3H^{2}}$ (35)

We formulate the autonomous coupled differential equations by using these
variables, Friedmann equation and the corresponding conservation equations for
each case as,

For the case NLI1:

$\begin{split}u^{\prime}=\frac{du}{dx}=3(bu-1)u-6(1-v)\left(\frac{v}{2}-\frac{1}{2}\right)=f(u,v),\\\
v^{\prime}=\frac{dv}{dx}=-3bu^{2}-6v\left(\frac{v}{2}-\frac{1}{2}\right)=g(u,v).\end{split}$
(36)

For the case NLI2:

$\begin{split}u^{\prime}=\frac{du}{dx}=3(bv^{2}-u)-6u\left(\frac{v}{2}-\frac{1}{2}\right)=f(u,v),\\\
v^{\prime}=\frac{dv}{dx}=-3bv^{2}-6(1-u)\left(\frac{v}{2}-\frac{1}{2}\right)=g(u,v).\end{split}$
(37)

For the case NLI3:

$\begin{split}u^{\prime}=\frac{du}{dx}=3u(bv-1)-6u\left(\frac{v}{2}-\frac{1}{2}\right)=f(u,v),\\\
v^{\prime}=\frac{dv}{dx}=-3buv-6(1-u)\left(\frac{v}{2}-\frac{1}{2}\right)=g(u,v).\end{split}$
(38)

Here the prime refers to differentiation with the variable $x=\ln a.$ The
equilibrium solutions are obtained by equating $u^{\prime}=0$ and
$v^{\prime}=0.$ The corresponding critical points for the three interaction
cases are, $(u_{\ast},v_{\ast})=(0,1),(\frac{1}{1+b},\frac{b}{1+b})$ for NLI1,
$(u_{\ast},v_{\ast})=(1,0),(\frac{b}{1+b},\frac{1}{1+b})$ for NLI2,
$(u_{\ast},v_{\ast})=(0,1),(1,0)$ for NLI3. The state $(0,1)$ corresponds to
the end de Sitter phase where as the state $(1,0)$ represents early matter-
dominated phase. The state $(\frac{1}{1+b},\frac{b}{1+b})$ produce a matter-
dominated epoch and that $(\frac{b}{1+b},\frac{1}{1+b})$ shows a DE-dominated
era. The stability of the dynamic system in the neighborhood of the critical
point can be checked by by considering small perturbations around the critical
point as $u=u_{\ast}+\delta u$, $v=v_{\ast}+\delta v$. Linearize the system of
equations (36),(37) and (38) with respect to these perturbation gives a matrix
equation,

$\begin{bmatrix}\delta u^{\prime}\\\ \delta
v^{\prime}\end{bmatrix}=\begin{bmatrix}(\frac{\partial f}{\partial
u})_{\ast}&(\frac{\partial f}{\partial v})_{\ast}\\\ (\frac{\partial
g}{\partial u})_{\ast}&(\frac{\partial g}{\partial
v})_{\ast}\end{bmatrix}\begin{bmatrix}\delta u\\\ \delta v\end{bmatrix}$ (39)

The Jacobian matrix i.e., $2\times 2$ matrix in the right of equation(39) for
the three cases of the present model are,

for NLI1,

$\begin{bmatrix}3bu_{\ast}+3(bu_{\ast}-1)&-3(1-v_{\ast})+3(v_{\ast}-1)\\\
-6bu_{\ast}&-3(v_{\ast}-1)-3v_{\ast}\end{bmatrix}$ (40)

for NLI2,

$\begin{bmatrix}-3&-3(1-v_{\ast})+3(-1+v_{\ast})+6bv_{\ast}\\\
0&-3(-1+v_{\ast})-3v_{\ast}-6bv_{\ast}\end{bmatrix}$ (41)

for NLI3,

$\begin{bmatrix}-3(-1+v_{\ast})+3(-1+bv_{\ast})&-3u_{\ast}+3bu_{\ast}\\\
3(-1+v_{\ast})&-3(1-v_{\ast})-3bu_{\ast}\end{bmatrix}$ (42)

Diagonalizing the above matrices gives the eigenvalues. The stability of the
critical points are inferred from the sign of eigenvalues of the Jacobian
matrix. If all the eigenvalues are positive, then the trajectories from the
neighbourhood of these points diverge away from that critical point and then
it is an unstable node. If all the eigenvalues are negative, the trajectories
from the neighbourhood of the critical point converge to that point, and then
such critical points are stable attractors. If it consists of both positive
and negative eigenvalues, then the critical point is said to be a saddle point
for which the trajectories are diverging or converging from these points. The
critical points, eigenvalues and nature of the eigenvalues obtained for the
three different interaction cases are summarized in the table(2)

Table 2: Critical points and eigenvalue values $Models$ | $Criticalpoints$ | $Eigenvalues$ | $Nature$
---|---|---|---
NLI1 | (0.997,0.003) | (-3,3) | saddle
(0,1) | (-3,-3) | stable
NLI2 | (1,0) | (-3,3) | saddle
(0.001,0.999) | (-3,-3) | stable
NLI3 | (1,0) | (-3,2.99) | saddle
(0,1) | (-3,-2.99) | stable

For the case NLI1: the critical point $(0.979,0.021)$ is a saddle point since
one of the eigenvalues is positive and the other is negative. The critical
point $(0,1)$ is a stable point for which the eigenvalues are negative. For
the case, NLI2:the critical point $(1,0)$ is a saddle point seeing that the
eigenvalues consist of positive and negative values. The critical point
$(0.072,0.928)$ is stable since the eigenvalues are negative. For the case,
NLI3: the critical point $(1,0)$ is a saddle point for which the eigenvalues
consist of positive and negative values. The critical point $(0,1)$ is stable
since all the eigenvalues are negative. The asymptotic behaviour of all three
cases predicts the transition of the universe from the early decelerated epoch
to the stable later accelerated epoch.

Figure 10: Phase space plot for NLI1 Figure 11: Phase space plot for NLI2
Figure 12: Phase space plot for NLI3

##  5 Conclusion

RVE model has emerged as a new look into the cosmological constant problem.
These models are well motivated from the renormalization approach of quantum
field theory in curved spacetime indicates the possibility of slowly varying
CC[11, 15, 42]. These models fit the current observational data quite
successfully compared to the $\Lambda$CDM. The general form of RVE density
consists of a combination of $H^{2}$ and $\dot{H}$ terms. This is analogous to
the holographic Ricci DE density, where a combination of $H^{2}$ and $\dot{H}$
terms are also present in it. So we consider holographic Ricci DE as a running
vacuum. In the earlier works, we consider holographic Ricci DE with an
additive constant in its energy density. We have shown that this additive
constant in the energy density is responsible for causing a transition from
prior decelerated to a late accelerated phase of the universe[35]. Then we
have shown that by accounting for the interaction between the dark sectors
through a phenomenological term $Q=3bH\rho_{m}$ transition from decelerated to
accelerated epoch is possible without an additive constant in the DE
density[36]. The cosmology with the interaction between DM and DE has gained
massive attention because several cosmic puzzles such as the cosmological
constant problem, coincidence problem as well as some recent observational
issues, such as the $H_{0}$ tension, $\sigma$8 tension, and the crossing of
the phantom barrier without any scalar field theory can be solved by the
allowance of an interaction between the dark sectors[71, 72, 73]. The
observational data collected from various astrophysical sources have favoured
interaction between the dark sectors. Since the precise nature of the
interaction between the dark sectors is yet unknown, the possibility for the
cosmological models permitting nonlinear interaction is equally acceptable to
understand the dynamics of the universe. The nonlinear interaction with
holographic DE was discussed in the references[37, 38, 39]. The authors in
reference[74, 75] test the viability of the interacting holographic DE model
by using the combined observational data. The present work thus examines the
cosmological dynamics in the presence of various nonlinear interaction rates
between DM and DE in a spatially flat FLRW universe.

In this work, we have explored a model of HRDE as a running vacuum by
accounting for a nonlinear interaction between the DM and the DE. The
nonlinear interaction of the form
$Q=3bH\rho^{m+s}\rho_{m}^{n}\rho_{de}^{-s-n}$ where b is the coupling
constant, powers s,n, and m specify interactions, generalizing several models
studied so far in the literature[39]. The models corresponding to various
combinations of m, n, and s were found to be different from that of the
standard $\Lambda$CDM model[38]. In particular, we have considered three
distinct nonlinear interaction forms for the corresponding values of m,n,s as
(m,n,s)=(1,2,-2), (m,n,s)=(1,0,-2), (m,n,s)=(1,1,-2) in our work which gives
analytically feasible solutions. For each of these cases, we analytically
solve for the Hubble parameter. To examine the observational acceptance of the
models, we perform global fittings of the models using the latest
observational data SN1a(Pantheon)+CMB(Planck2018)+BAO(SDSS) data set. We have
extracted the value of model parameter $\beta$, interaction term $b$,
$\Omega_{de0}$, $H_{0}$, and $M$ for all the three interaction cases using the
$\chi^{2}$ minimization technique. The parameters of all the three cases are
found to behave similar values as $\beta\sim 0.45$, $\Omega_{de0}\sim 0.72$,
$H_{0}\sim 68.7$, $M\sim 19.39.$ The value of the interaction term b for NLI1
and NLI3 is $0.003$ whereas, for NLI2, its value is $0.001.$ The positive
value of the interaction parameter confirms the transfer of energy from DE to
the DM. The evolution of the Hubble parameter using these constraints of model
parameter for three interaction cases predicts the presence of the big bang at
the origin.

We also derive the deceleration parameter for the three cases of nonlinear
interaction. The present value of the deceleration parameter, $q_{0}$, agrees
with the observational results. The evolution of the q-z plot for the three
cases predicts a transition from decelerated to the accelerated epoch at the
redshift, $z_{T}\sim 0.72$, which agrees with the range of $\Lambda$CDM model.

The state finder diagnosis of the three cases shows that this model is
distinguishably different from the $\Lambda$CDM model. The trajectories of
$\\{r,s\\}$ parameter of NLI1 and NLI3 are approaching the $\Lambda$CDM fixed
point in the future, whereas that of NLI2 is diverging from the $\Lambda$CDM
fixed point. The trajectories of three cases lie in the region $r<1$, and
$s>0$ resembles the nature of quintessence DE.

The phase space analysis of the three cases was analyzed. The critical points
obtained having both positive and negative eigenvalues represent the prior
unstable matter-dominated phase. The critical points having both negative
eigenvalues correspond to the stable accelerated epoch. For the case NLI1: the
critical point $(0.997,0.003)$ represents an unstable matter-dominated phase,
and the critical point $(0,1)$ gives a stable de Sitter phase. For NLI2, the
critical point $(1,0)$ corresponds to an unstable matter-dominated epoch,
whereas the critical point $(0.001,0.999)$ shows an unstable DE-dominated
epoch. For NLI3, the critical point $(1,0)$ represents an unstable matter-
dominated early phase, and the critical point $(0,1)$ corresponds to the
stable DE-dominated late phase. Our analysis shows that all the three
interacting cases predict a prior matter-dominated universe with an unstable
node and a stable DE-dominated universe with acceleration, similar to the de
Sitter phase.

The authors in the reference[76] interpret the dark sector interaction model
in which the coupling form accounts for the interaction between the DM and the
DE $Q=\frac{\acute{\rho}}{3(\alpha+\beta\rho)}$ as a particular case of a
running vacuum model of the type $\Lambda(H)=n_{0}+n_{1}H^{2}+n_{2}H^{4}$ in
which VE decays to DM. In the reference[77] authors show that the nonlinear
interaction between dark sectors results in an ISW(Integrated Sachs–Wolfe)
effect. A more quantitative analysis of these interesting issues of non-
linearly interacting models will be the subject of future work.

## References

  * [1] Perlmutter, S. et al.,  Astrophys. J. 517 (1999) 565.
  * [2] Riess, A. G. et al.,  Astron. J. 116 (1998) 1009.
  * [3] Tegmark, M. et al.,  Phys. Rev. D69 (2004) 103501.
  * [4] Abazajian, K. et al.,  Astron. J. 128 (2004) 502.
  * [5] Peiris, H. V. et al.,  Astrophys. J. Suppl. 148 (2003) 213.
  * [6] Spergel, D. N. et al.,  Astrophys. J. Suppl. 148 (2003) 175.
  * [7] Spergel, D. N. et al.,  Astrophys. J. Suppl. 170 (2007) 377.
  * [8] Davis, T. M. et al.,  Astrophys. J. 666 (2007) 716.
  * [9] Komatsu, E. et al.,  Astrophys. J. Suppl. 180 (2009) 330.
  * [10] Aghanim, N. et al.,  ArXiv:1807.06209(2018) .
  * [11] Sola, J. and Stefancic, H.,  Phys. Lett. B624 (2005) 147.
  * [12] Basilakos, S., Plionis, M., and Solà, J.,  Phys. Rev. D80 (2009) 083511.
  * [13] Nelson, B. L. and Panangaden, P.,  Phys. Rev. D 25 (1982) 1019.
  * [14] Buchbinder, I. L.,  Fortsch. Phys. 34 (1986) 605.
  * [15] Sola, J.,  J. Phys. A41 (2008) 164066.
  * [16] Shapiro, I. L. and Sola, J.,  JHEP 02 (2002) 006.
  * [17] Basilakos, S., Plionis, M., and Sola, J.,  Phys. Rev. D82 (2010) 083512.
  * [18] Grande, J., Sola, J., Basilakos, S., and Plionis, M.,  JCAP 1108 (2011) 007.
  * [19] Basilakos, S., Polarski, D., and Sola, J.,  Phys. Rev. D86 (2012) 043010.
  * [20] Fabris, J. C., Shapiro, I. L., and Sola, J.,  JCAP 0702 (2007) 016.
  * [21] Rezaei, M., Malekjani, M., and Sola, J.,  Phys. Rev. D 100 (2019) 023539.
  * [22] Fritzsch, H., Solà, J., and Nunes, R. C.,  Eur. Phys. J. C 77 (2017) 193.
  * [23] Solà, J., Gómez-Valent, A., and de Cruz Pérez, J.,  Astrophys. J. 836 (2017) 43.
  * [24] Li, M.,  Phys. Lett. B603 (2004) 1.
  * [25] Hsu, S. D. H.,  Phys. Lett. B594 (2004) 13.
  * [26] Bouhmadi-Lopez, M., Errahmani, A., and Ouali, T.,  Phys. Rev. D84 (2011) 083508.
  * [27] Saridakis, E. N.,  Phys. Lett. B660 (2008) 138.
  * [28] Gao, C., Wu, F., Chen, X., and Shen, Y.-G.,  Phys. Rev. D79 (2009) 043511.
  * [29] Nojiri, S. and Odintsov, S. D.,  Gen. Rel. Grav. 38 (2006) 1285.
  * [30] Huang, Q.-G. and Li, M.,  JCAP 0408 (2004) 013.
  * [31] Jamil, M., Farooq, M. U., and Rashid, M. A.,  Eur. Phys. J. C61 (2009) 471.
  * [32] Chen, B., Li, M., and Wang, Y.,  Nucl. Phys. B774 (2007) 256.
  * [33] Feng, C.-J.,  Phys. Lett. B672 (2009) 94.
  * [34] Chang, Z., Wu, F.-Q., and Zhang, X.,  Phys. Lett. B633 (2006) 14.
  * [35] George, P. and Mathew, T. K.,  Mod. Phys. Lett. A31 (2016) 1650075.
  * [36] George, P., Shareef, V. M., and Mathew, T. K.,  Int. J. Mod. Phys. D28 (2018) 1950060.
  * [37] Oliveros, A. and Acero, M. A.,  Astrophys. Space Sci. 357 (2015) 12.
  * [38] Bolotin, Y. L., Kostenko, A., Lemets, O. A., and Yerokhin, D. A.,  Int. J. Mod. Phys. D 24 (2015) 1530007.
  * [39] Arevalo, F., Bacalhau, A. P. R., and Zimdahl, W.,  Class. Quant. Grav. 29 (2012) 235001.
  * [40] Baldi, M.,  Mon. Not. Roy. Astron. Soc. 411 (2011) 1077.
  * [41] He, J.-H. and Wang, B.,  JCAP 0806 (2008) 010.
  * [42] Sola, J.,  J. Phys. Conf. Ser. 283 (2011) 012033.
  * [43] Parker, L. E. and Toms, D.,  Quantum Field Theory in Curved Spacetime: Quantized Field and Gravity,  Cambridge Monographs on Mathematical Physics, Cambridge University Press, 2009.
  * [44] Lima, J. and Trodden, M.,  Phys. Rev. D 53 (1996) 4280.
  * [45] Sola, J.,  AIP Conf. Proc. 1606 (2015) 19.
  * [46] Lima, J., Basilakos, S., and Sola, J.,  Mon. Not. Roy. Astron. Soc. 431 (2013) 923.
  * [47] Lima, J. A. S., Basilakos, S., and Sola, J.,  Gen. Rel. Grav. 47 (2015) 40.
  * [48] Lima, J. A. S., Basilakos, S., and Solà, J.,  Eur. Phys. J. C76 (2016) 228.
  * [49] Solà, J. and Gómez-Valent, A.,  Int. J. Mod. Phys. D24 (2015) 1541003.
  * [50] Cohen, A. G., Kaplan, D. B., and Nelson, A. E.,  Phys. Rev. Lett. 82 (1999) 4971.
  * [51] Gao, C., Wu, F., Chen, X., and Shen, Y.-G.,  Phys. Rev. D 79 (2009) 043511.
  * [52] del Campo, S., Herrera, R., and Pavón, D.,  Phys. Rev. D 78 (2008) 021302.
  * [53] Pan, S., Bhattacharya, S., and Chakraborty, S.,  Mon. Not. Roy. Astron. Soc. 452 (2015) 3038.
  * [54] Pereira, S. H. and Jesus, J. F.,  Phys. Rev. D 79 (2009) 043517.
  * [55] Wang, B., Lin, C.-Y., and Abdalla, E.,  Physics Letters B 637 (2006) 357–361.
  * [56] Ma, Y.-Z., Gong, Y., and Chen, X.,  Eur. Phys. J. C69 (2010) 509.
  * [57] Mangano, G., Miele, G., and Pettorino, V.,  Mod. Phys. Lett. A18 (2003) 831.
  * [58] Sadri, E. and Khurshudyan, M.,  Int. J. Mod. Phys. D28 (2019) 1950152.
  * [59] Scolnic, D. M. et al.,  Astrophys. J. 859 (2018) 101.
  * [60] Chen, L., Huang, Q.-G., and Wang, K.,  JCAP 02 (2019) 028.
  * [61] Padmanabhan, N. et al.,  Mon. Not. Roy. Astron. Soc. 378 (2007) 852.
  * [62] Blake, C. et al.,  Mon. Not. Roy. Astron. Soc. 418 (2011) 1707.
  * [63] Eisenstein, D. J. et al.,  Astrophys. J. 633 (2005) 560.
  * [64] Tegmark, M. et al.,  Phys. Rev. D 69 (2004) 103501.
  * [65] Bond, J., Efstathiou, G., and Tegmark, M.,  Mon. Not. Roy. Astron. Soc. 291 (1997) L33.
  * [66] Alam, U., Sahni, V., and Starobinsky, A. A.,  JCAP 0406 (2004) 008.
  * [67] Hinshaw, G. et al.,  Astrophys. J. Suppl. 180 (2009) 225.
  * [68] Ade, P. A. R. et al.,  Astron. Astrophys. 594 (2016) A13.
  * [69] Sahni, V., Saini, T. D., Starobinsky, A. A., and Alam, U.,  JETP Lett. 77 (2003) 201.
  * [70] Zhang, X.,  Phys. Lett. B611 (2005) 1.
  * [71] Di Valentino, E., Melchiorri, A., Mena, O., and Vagnozzi, S.,  Phys. Rev. D 101 (2020) 063502.
  * [72] Kumar, S., Nunes, R. C., and Yadav, S. K.,  Eur. Phys. J. C 79 (2019) 576.
  * [73] Pan, S., Yang, W., and Paliathanasis, A.,  Mon. Not. Roy. Astron. Soc. 493 (2020) 3114.
  * [74] Feng, C., Wang, B., Gong, Y., and Su, R.-K.,  JCAP 09 (2007) 005.
  * [75] Sadri, E., Khurshudyan, M., and Zeng, D.-f.,  (2019).
  * [76] De-Santiago, J., Sánchez G., I. E., and Tamayo, D.,  Gen. Rel. Grav. 50 (2018) 101.
  * [77] Bouhmadi-López, M., Morais, J., and Zhuk, A.,  Phys. Dark Univ. 14 (2016) 11.