# Robust quantum entanglement and teleportation in the tetrapartite spin-1/2
square clusters: Theoretical study on the effect of a cyclic four-spin
exchange

Hamid Arian Zad1,2, Azam Zoshki2, Nerses Ananikian1,3, Michal Jaščur2 1A.I.
Alikhanyan National Science Laboratory, 0036, Yerevan, Armenia  2Department of
Theoretical Physics and Astrophysics, Faculty of Science, P. J. Sǎfárik
University, Park Angelinum 9, 041 54 Košice, Slovak Republic  3CANDLE
Synchrotron Research Institute, Acharyan 31, 0040 Yerevan, Armenia

###### Abstract

The whole entanglement measure so-called geometric $\Pi_{4}$ average of
tangles and bipartite entanglement of the antiferromagnetic spin-1/2 XXX
Heisenberg model on a tetranuclear square cluster with cyclic four-spin
interaction are rigorously examined by the help of thermal negativities. The
model comprises two nearest-neighbor exchange couplings $J_{1}$ and $J_{2}$
such that $J_{1}\gg J_{2}$. When the cyclic four-spin exchange is zero, the
maximum value of whole entanglement $\Pi_{4}$ is achieved at low enough
temperatures and relatively high magnetic fields $(B\approx J_{1})$. Also,
maximum bipartite entanglement between pair spins with exchange coupling
$J_{1}$ is achievable at high temperature and high magnetic field. This
quantity remains alive for sufficiently high temperature and high magnetic
field values comparable with the relevant exchange coupling $J_{1}$. A nonzero
value of the cyclic four-spin exchange notably enhances the degree of the
whole entanglement, while it weakens the bipartite entanglement degree. We
demonstrate that the whole entanglement reaches an unconventional minimum at a
special parameter region of cyclic four-spin exchange almost ten order of
magnitude smaller than $J_{1}$, where the system is in a quantum
antiferromagnetic state. The real complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ as
a strong antiferromagnetic tetranuclear square compound provides us an
experimental representative to estimate the strength of the whole and
bipartite entanglements at high enough temperature. It is demonstrated that
the entanglement negativities of this complex are yet depend on the considered
cyclic four-spin interaction even though its magnitude is significantly
smaller than $J_{1}$. We realize that the possibility of quantum teleportation
through a couple of discrete $\text{Cu}_{4}^{\text{II}}$ complexes would
solely happen for the magnetic fields lower than a critical value. It is also
demonstrated that enhancing the cyclic four-spin exchange results in
diminishing the average fidelity. Moreover, we conclude that the quantum
Fisher information investigated for the pair of spins with exchange
interaction $J_{1}$ reaches its minimum close to the critical magnetic field
at which a ground-state phase transition occurs from highly entangled state to
separable one.

††preprint: APS/123-QED

##  I Introduction

Low-dimensional spin systems are of particular interest as magnetic materials
for recent developments in quantum information technologies and molecular
spintronics [1, 2, 3, 4, 5, 6]. Scientists have recently been interested in
the application of antiferromagnets by which the information could be encoded
in the direction of the Néel vector. This is due to the fact that
antiferromagnetic materials are not only able to write information much faster
but also are very stable against external fields [7].

The characterized tetranuclear copper(Cu${}^{\text{II}}$) complexes as
exchange-coupled systems have attracted a great deal of attention since past
decade owing to their potential applications in a number of branches such as
coordination polymers [8, 9], molecular magnetism [1], high-temperature
superconductors [10], string theory [11] and quantum information processing
[7]. Among them, strong antiferromagnetically tetranuclear Cu${}^{\text{II}}$
complexes have particularly studied due to bear magneto-structural properties
and reveal zero-temperature phase transitions between possible ground states.
For example, magnetic properties and crystal structure of a tetranuclear
cubane-type Cu${}^{\text{II}}$ complex based on a hexa-coordinating
$\text{N}_{2}\text{O}_{4}$ ligand have been reported in Ref. [12]. Next, S. C.
Manna et al. synthesized and characterized two tetranuclear closed-cubane like
complexes $[\text{Cu}_{4}(\text{L}^{1})_{4}]\cdot 3(\text{H}_{2}\text{O})$ and
$[\text{Cu}_{4}(\text{H}_{2}\text{L}^{2})_{4}(\text{H}_{2}\text{O})_{4}]$
[13], and based on the Tercero et al. [14] theoretical tasks, they determined
the nature and magnitude of exchange couplings between Cu atoms of these
complexes that are correlated with the $\text{Cu}-\text{O}-\text{Cu}$ angle.

The structural and magnetic properties of the ferromagnetically coupled
tetranuclear copper $\text{Cu}_{4}^{\text{II}}$ square complex
$[{\text{Cu}(\text{H}_{2}\text{L})}_{4}(\text{THF})]$ were experimentally
investigated in Ref. [15]. The spin arrangement of this special complex in the
crystal lattice leads to the formation of square structure. From the expected
values for the six coupling constants, four intracubane (square sides)
interactions and two intercubane (diagonals) interactions, intercubane
interactions are negligible and one can remove them from the magnetic coupling
pathway scheme. Owing to this fact, intracubane interactions have been often
observed as strong interactions in the fitting procedure, therefore the main
contribution being the coupling along the square sides, and thus the diagonal
contribution should be taken in to account with caution.

The crystal structure and magnetic properties of a tetranuclear
Cu${}^{\text{II}}$ complex containing the 2-(pyridine-2-yliminomethyl)-phenol
ligand have been studied in detail by S. Giri et al. in Ref. [16]. Indeed,
they synthesized a pure $\mathrm{Cu}^{\mathrm{II}}$ compound,
$[\mathrm{Cu}_{4}\mathrm{L}_{4}(\mathrm{H}_{2}\mathrm{O})_{4}](\mathrm{ClO}_{4})_{4}$,
that is a Schiff base ligand with the perchlorate counter ion. The structure
of this complex includes discrete planar tetranuclear
$[\mathrm{Cu}_{4}\mathrm{L}_{4}(\mathrm{H}_{2}\mathrm{O})_{4}]^{4+}$ cations
together with perchlorate anions. The four copper atoms that constitute a unit
tetrameric block are crystallographically equivalent and are five-coordinate
with square pyramidal environments. Fitting process of the magnetic
susceptibility data with the experiment indicated that the compound has a
strong antiferromagnetic, $J_{1}/k_{B}=-638\;\text{K}$, exchange interaction
between diphenoxo-bridged Cu${}^{\text{II}}$ centers, as well as, a moderate
antiferromagnetic, $J_{2}/k_{B}=-34\;\text{K}$, interaction between
$\text{N}-\text{C}-\text{N}$ bridged Cu${}^{\text{II}}$ centers. Such planar
cyclic $\mathrm{Cu}_{4}\mathrm{O}_{4}$ core structures involve a
$\mathrm{Cu}_{2}\mathrm{O}_{2}$ bridging moiety and have attracted a
considerable interest due to be momentous in seeking molecular precursor
compounds for synthesizing high-temperature superconductors [17, 18, 19].

Non-classical correlations between different parties of a quantum system can
be described by quantum entanglement where it plays as an important hallmark
that is able to separate quantum systems from their classical counterparts. To
quantify the entanglement degree of a given bipartite quantum state two
popular measures: Concurrence and negativity can be applied [20, 21, 22]. The
negativity [22] as a measure of either bipartite or multipartite entanglement
can be considered according to the necessary condition Peres–Horodecki
separability criterion [23, 24]. This criterion is called the PPT (positive
partial transpose) criterion and determines the separability of a two quantum
mechanical subsystems A and B of a global system (AB). Owing to this fact, the
quantity negativity has been widely employed to measure the degree of the
thermal entanglement of the multipartite Heisenberg spin models [25, 26, 27,
28, 29].

Motivated by the complex,
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$
[16] and its likely application in superconducting quantum computers at high
temperatures, in the present work, we consider this compound as an
experimental sample and proceed ourselves to theoretically investigate various
quantum protocols such as bipartite and tetrapartite entanglement negativities
[30, 33, 31, 32], quantum teleportation [34, 35, 36, 37, 38], and quantum
Fisher information (QFI) [39, 40, 41, 42, 43, 44, 45] over the room
temperature scales. The QFI is a useful tool that makes the possibility of
estimating the quantum metrology of a spin system and can be plausibly
connected to the quantum phase transition at higher temperatures. Hence, we
disclose anomalous behavior of these quantities nearby the critical points at
which a ground-state phase transition occurs.

The most important feature of spin systems is the type of interaction between
their ingredients. It has already been pronounced that in addition to the
bilinear exchange interaction, cyclic four-spin interactions known as the ring
exchange terms are momentous in the subject [46, 47, 48, 49, 50, 51, 52, 53,
54, 55]. In Ref. [56] the cyclic four-spin exchange term was detected in some
cuprates and calculated theoretically. The authors of Ref. [57] investigated
the susceptibility and specific heat of a two-leg $S-1/2$ ladder possessing
cyclic four-spin exchange and realized that a small cyclic exchange
interaction remarkably affects the thermodynamics of the model even at high
temperature. However, realizing the effects of ring exchange interaction on
the magnetic and quantum properties of small cluster magnets specially
tetranuclear Cu${}^{\text{II}}$ complexes is yet an open question and has not
been seriously attempted to answer. Therefore, in the present work we
theoretically address the question: How the bipartite and tetrapartite quantum
entanglements of the spin-1/2 antiferromagnetic Heisenberg model on a
tetranuclear square compound (i.g.,
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$)
are changed by cyclic four-spin exchange interaction?

This paper is organized as follows. In Sec. II we describe the spin-1/2
antiferromagnetic Heisenberg XXX model on a tetranuclear cluster magnet with
the cyclic four-spin interaction in the presence of an external magnetic
field. Next, we classify its eigenenergies with the corresponding eigenstates
and the most interesting results for the ground-state phase diagram is
figured, as well. In Sec. III, firstly we theoretically describe the
influences of a cyclic four-spin interaction on the whole entanglement
negativity of the model under consideration. Secondly, we consider the real
complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ as
an experimental example of the strong antiferromagnetic Heisenberg model on a
tetranuclear square compound that has a very strong exchange coupling
$J_{1}/k_{B}=-638\;\mathrm{K}$ and a moderate exchange coupling
$J_{2}/k_{B}=-34\;\mathrm{K}$, then we discuss the effects of a typical cyclic
four-spin exchange on the whole entanglement and bipartite quantum
entanglement between Cu ions with stronger antiferromagnetic interaction. In
the final part of this section, by quantifying average fidelity, we rigorously
examine the possibility of the quantum teleportation through this model at
high temperature. Furthermore, the behavior of QFI is compared with the
average fidelity and the most stimulating findings are reported, as well.
Finally, the concluding remarks are given in Sec. IV.

##  II The model

First let us define the Hamiltonian of the spin-1/2 Heisenberg model on a
tetranuclear square compound as

$\mathcal{H}=\sum\limits_{i=1}^{N}{H}_{i},$ (1)

where ${H}_{i}$ is the effective Hamiltonian of a tetrameric unit block,
indicating a discrete four-spin system (see supporting information of Ref.
[16]).

Figure 1: Structural arrangement of the spin-1/2 Heisenberg model on a square
compound with magnetic exchange interactions $J_{1}$, $J_{2}$ and cyclic four-
spin term $K$.

We suppose that the form of Hamiltonian of all unit blocks are identical (see
Fig. 1), therefore Hamiltonian $H={H}_{i}$ can be expressed as

$\begin{array}[]{lcl}{H}=-J_{1}\left[{\bm{S}}_{\mathrm{A}}\cdot{\bm{S}}_{\mathrm{B}}+{\bm{S}}_{\mathrm{C}}\cdot{\bm{S}}_{\mathrm{D}}\right]-J_{2}\left[{\bm{S}}_{\mathrm{B}}\cdot{\bm{S}}_{\mathrm{C}}+{\bm{S}}_{\mathrm{A}}\cdot{\bm{S}}_{\mathrm{D}}\right]\\\
\quad\quad\;+K[({\bm{S}}_{\mathrm{A}}\cdot{\bm{S}}_{\mathrm{B}})({\bm{S}}_{\mathrm{C}}\cdot{\bm{S}}_{\mathrm{D}})+({\bm{S}}_{\mathrm{A}}\cdot{\bm{S}}_{\mathrm{D}})({\bm{S}}_{\mathrm{B}}\cdot{\bm{S}}_{\mathrm{C}})\\\
\quad\quad\;-({\bm{S}}_{\mathrm{A}}\cdot{\bm{S}}_{\mathrm{C}})({\bm{S}}_{\mathrm{B}}\cdot{\bm{S}}_{\mathrm{D}})]-\mu_{B}gB(S_{\mathrm{A}}^{z}+S_{\mathrm{B}}^{z}+S_{\mathrm{C}}^{z}+S_{\mathrm{D}}^{z}).\end{array}$
(2)

where
${\bm{S}}_{a}\cdot{\bm{S}}_{b}=S_{a}^{x}S_{b}^{x}+S_{a}^{y}S_{b}^{y}+S_{a}^{z}S_{b}^{z}$
with $a,b=\\{\mathrm{A,B,C,D}\\}$. The coupling $J_{i}$ ($i=\\{1,2\\}$)
describes the strength of the spin interaction being antiferromagnetic
$J_{i}<0$ and ferromagnetic $J_{i}>0$. The local asymmetric formation of the
four different chemical environments of the copper centers results in two
different magnetic exchange interactions including $\\{J_{1},J_{2}\\}$
correspond to the square sides that play the main role of exchange
contribution [16]. $B$ is the external magnetic field which is considered in
the $z$-direction. $K$ is the cyclic four-spin exchange interaction. $\mu_{B}$
is Bohr magneton, $g$ is Landé g-factor that motivated by the experimental
representative considered in this work we set $g=2.19$. $S_{i}^{\alpha}$
($\alpha=x,y,z$) are the spin-1/2 operators with $\hbar=1$.

By implementing the exact analytical diagonalization one can directly solve
the model and derive the eigenvalues of the Hamiltonian $\mathcal{H}_{i}$ as

$\begin{array}[]{lcl}\left\\{\begin{array}[]{@{}l@{}}E_{1}=g\mu_{B}B-\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{2}=g\mu_{B}B-\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{3}=g\mu_{B}B+\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{4}=g\mu_{B}B+\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}-\dfrac{7}{16}K,\\\\[8.5359pt]
\end{array}\right.\quad\quad\left\\{\begin{array}[]{@{}l@{}}E_{5}=\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}+\dfrac{5}{16}K-\dfrac{1}{2}\sqrt{\Lambda},\\\\[8.5359pt]
E_{6}=\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}+\dfrac{5}{16}K+\dfrac{1}{2}\sqrt{\Lambda},\\\\[8.5359pt]
E_{7}=-\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{8}=-\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{9}=-g\mu_{B}B+\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{10}=\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\end{array}\right.\quad\left\\{\begin{array}[]{@{}l@{}}E_{11}=\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}-\dfrac{7}{16}K,\\\\[7.11317pt]
E_{12}=-g\mu_{B}B+\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}-\dfrac{7}{16}K,\\\\[7.11317pt]
E_{13}=-g\mu_{B}B-\dfrac{1}{2}J_{1}+\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{14}=-g\mu_{B}B+\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
\end{array}\right.\\\\[21.33955pt]
\left\\{\begin{array}[]{@{}l@{}}E_{15}=-2g\mu_{B}B-\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\\\\[8.5359pt]
E_{16}=2g\mu_{B}B-\dfrac{1}{2}J_{1}-\dfrac{1}{2}J_{2}+\dfrac{1}{16}K,\end{array}\right.\par\end{array}$

by adopting $\Lambda=4J_{1}^{2}-4J_{1}J_{2}+4J_{2}^{2}+2J_{1}K+2J_{2}K+K^{2}.$
The eigenvectors of the $\mathcal{H}_{i}$ in terms of standard basis are
classified according to the total quantum spin number
$S_{\text{T}}=S_{\mathrm{A}}+S_{\mathrm{B}}+S_{\mathrm{C}}+S_{\mathrm{D}}$ and
its $z-$component $S^{z}_{\text{T}}$ as follows

$\begin{array}[]{lcl}(S_{\text{T}}=1,S^{z}_{\text{T}}=-1):\left\\{\begin{array}[]{@{}l@{}}\left|\psi_{1}\right\rangle=\left|1000\right\rangle+\left|0100\right\rangle+\left|0010\right\rangle+\left|0001\right\rangle,\\\\[0.71114pt]
\left|\psi_{2}\right\rangle=-\left|1000\right\rangle-\left|0100\right\rangle+\left|0010\right\rangle+\left|0001\right\rangle,\\\\[0.71114pt]
\left|\psi_{3}\right\rangle=\left|1000\right\rangle-\left|0100\right\rangle-\left|0010\right\rangle+\left|0001\right\rangle,\\\\[0.71114pt]
\left|\psi_{4}\right\rangle=-\left|1000\right\rangle+\left|0100\right\rangle-\left|0010\right\rangle+\left|0001\right\rangle,\end{array}\right.\\\\[21.33955pt]
\quad\quad\quad\quad(S_{\text{T}}=0):\;\left\\{\begin{array}[]{@{}l@{}}\left|\psi_{5}\right\rangle=a(\left|1100\right\rangle+\left|0011\right\rangle)+b(\left|1010\right\rangle+\left|0101\right\rangle)-c(\left|1001\right\rangle+\left|0110\right\rangle),\\\
\left|\psi_{6}\right\rangle=d(\left|1100\right\rangle+\left|0011\right\rangle)-e(\left|1010\right\rangle+\left|0101\right\rangle)-f(\left|1001\right\rangle+\left|0110\right\rangle),\\\
\left|\psi_{7}\right\rangle=\left|1100\right\rangle+\left|1010\right\rangle+\left|1001\right\rangle+\left|0110\right\rangle+\left|0101\right\rangle+\left|0011\right\rangle,\\\
\left|\psi_{8}\right\rangle=-\left|1010\right\rangle+\left|0101\right\rangle,\\\
\left|\psi_{9}\right\rangle=-\left|1100\right\rangle+\left|0011\right\rangle,\\\
\left|\psi_{10}\right\rangle=-\left|1001\right\rangle+\left|0110\right\rangle,\end{array}\right.\\\\[21.33955pt]
\;\;(S_{\text{T}}=1,S^{z}_{\text{T}}=1):\;\left\\{\begin{array}[]{@{}l@{}}\left|\psi_{11}\right\rangle=\left|1110\right\rangle+\left|1101\right\rangle+\left|1011\right\rangle+\left|0111\right\rangle,\\\
\left|\psi_{12}\right\rangle=-\left|1110\right\rangle+\left|1101\right\rangle-\left|1011\right\rangle+\left|0111\right\rangle,\\\
\left|\psi_{13}\right\rangle=-\left|1110\right\rangle-\left|1101\right\rangle+\left|1011\right\rangle+\left|0111\right\rangle,\\\
\left|\psi_{14}\right\rangle=\left|1110\right\rangle-\left|1101\right\rangle-\left|1011\right\rangle+\left|0111\right\rangle,\\\
\end{array}\right.\\\\[21.33955pt] (S_{\text{T}}=2,S^{z}_{\text{T}}=\pm
2):\;\left\\{\begin{array}[]{@{}l@{}}\left|\psi_{15}\right\rangle=\left|1111\right\rangle,\\\
\left|\psi_{16}\right\rangle=\left|0000\right\rangle,\end{array}\right.\end{array}$

where 1 and 0 are nominated as spins-up and spins-down, respectively.
$\left|1\right\rangle$ and $\left|0\right\rangle$ are the eigenbasis of the 2
by 2 Pauli matrix $\sigma^{z}=\text{diag}(1,-1)$. Coefficients
$\\{a,b,c,d,e,f\\}$ are so long polynomial functions of $B$, $J_{1}$, $J_{2}$
and $K$ that for simplicity we numerically obtain them in order to examine
forthcoming quantities.

###  II.1 Ground-state phase transition

The spin-1/2 tetranuclear square complex with Hamiltonian (2) exhibits four
different ground states in the presence of an external magnetic field. The
ground-state phase diagram of the model involves two localized one spin-down
phases, $|\psi_{11}\rangle$ and $|\psi_{13}\rangle$, the localized two spin-up
phase $|\psi_{5}\rangle$ and the fully polarized phase $|\psi_{15}\rangle$.
These ground states can be represented by the form
$|S_{\text{T}},S_{\mathrm{AB}},S_{\mathrm{CD}}\rangle$ [58, 59] as

$\begin{array}[]{lcl}|\psi_{5}\rangle=|{0,1,1}\rangle\quad\text{with}\quad
M/M_{\text{s}}=0,\\\\[5.69046pt]
|\psi_{11}\rangle=|{1,1,1}\rangle_{\text{I}}\quad\text{with}\quad
M/M_{\text{s}}=\dfrac{1}{2},\\\\[5.69046pt]
|\psi_{13}\rangle=|{1,1,1}\rangle_{\text{II}}\quad\text{with}\quad
M/M_{\text{s}}=\dfrac{1}{2},\\\\[5.69046pt]
|\psi_{15}\rangle=|2,1,1\rangle\quad\text{with}\quad
M/M_{\text{s}}=1,\end{array}$ (3)

where $S_{\mathrm{AB}}=S_{\mathrm{A}}+S_{\mathrm{B}}$ and
$S_{\mathrm{CD}}=S_{\mathrm{C}}+S_{\mathrm{D}}$ are spin numbers of two faced
dimers. $M=-(\partial{\mathcal{G}}/\partial{B})_{T}$ represents the total
magnetization where $\mathcal{G}=-k_{B}T\ln Z$ is the Gibbs free energy with
$Z$ being the partition function of the model. $T$ is the temperature, $k_{B}$
is Boltzmann constant and $M_{\text{s}}$ is the saturation magnetization.

Figure 2: (a) Ground-state phase diagram in the $\mu_{B}B/|J_{1}|-K/|J_{1}|$
plane by assuming fixed $J_{2}/|J_{1}|=-0.053$ that indicates both exchange
couplings are antiferromagnetic such that $|J_{1}|$ is assumed to be the
energy unit and $J_{2}/|J_{1}|=-0.053$ being very weak antiferromagnetic with
respect to the $|J_{1}|$. The triple point at which three different ground
states coexist together is indicated by circle. Ground-state phase diagram in
the $\mu_{B}B/|J_{1}|-J_{2}/|J_{1}|$ plane by assuming; (b) $K=0$, (c)
$K/|J_{1}|=0.5$, (d) $K/J|J_{1}|=0.75$. Phase boundary between each two
different states is depicted by solid line with different color. Vertical
dotted line in lower panels shows the phase boundary between two ground states
$|{1,1,1}\rangle_{\text{I}}$ and $|{1,1,1}\rangle_{\text{II}}$ with the same
magnetization value $M/M_{\text{s}}=1/2$.

In Fig. 2 (a), we display the zero-temperature ground-state phase diagram in
the $\mu_{B}B/|J_{1}|-K/|J_{1}|$ plane where weak antiferromagnetic coupling
constant $J_{2}/|J_{1}|=-0.053$ is assumed. We choose this value according to
the experimental findings for the nature and strength of the exchange
couplings reported in Ref. [16]. It is noteworthy that here we evaluate
parameters $\\{\mu_{B}B,\;J_{2},\;K\\}$ with respect to the $|J_{1}|$ (where
$J_{1}<0$ being antiferromagnetic) in order to provide a dimensionless
parametric framework. This consideration not only enables one to use the
experimental data of exchange couplings, but also paves the way for readers to
apply our expressions for similar $\mathrm{Cu}_{4}^{\mathrm{II}}$ complexes
with different types of exchange interactions. Consequently, one finds three
different ground states: quantum ferromagnetic $|{1,1,1}\rangle_{\text{I}}$,
quantum antiferromagnetic $|{0,1,1}\rangle$ and classical ferromagnetic
$|{2,1,1}\rangle$ for $J_{2}/|J_{1}|<0$. The boundary between each two ground
states of Eq. (3) is given in table 1. For example, the boundary between
$|{1,1,1}\rangle_{\text{I}}$ and $|{0,1,1}\rangle$ phases is given by
$\mu_{B}B/|J_{1}|=-\frac{3}{4}K/|J_{1}|+\frac{1}{2|J_{1}|}\sqrt{\Lambda}$,
while between $|{1,1,1}\rangle_{\text{I}}$ and $|{2,1,1}\rangle$ phases is
given by $\mu_{B}B/|J_{1}|=-1-J_{2}/|J_{1}|+\frac{1}{2}K/|J_{1}|$. The phase
diagram pleasantly coincides the magnetization jumps and plateaus of the model
under consideration. All three ground states coexist together at the triple
point $K_{\text{c}}/|J_{1}|=2J_{2}/|J_{1}|=-0.106$ situated nearby the phase
boundary between two states $|{1,1,1}\rangle_{\text{I}}$ and
$|{1,1,1}\rangle_{\text{II}}$ with the same magnetization value. We discuss
this important term in what follows.

In Fig. 2 (b) the ground-state phase diagram of the model is shown in the
$\mu_{B}B/|J_{1}|-J_{2}/|J_{1}|$ plane when $K=0$. The same ground states to
panel Fig. 2 (a) are visible when antiferromagnetic case $J_{2}/|J_{1}|<0$ is
assumed, while for ferromagnetic range $J_{2}/|J_{1}|>0$ a different ground
state $|{1,1,1}\rangle_{\text{II}}$ with $M/M_{\text{s}}=1/2$ appears. Here, a
quadruple point is obvious at which four ground states coexist together.

The most interesting results of the ground-state phase diagram within the
$\mu_{B}B/|J_{1}|-J_{2}/|J_{1}|$ plane is achieved when the model involves
with a non zero cyclic four-spin exchange term. We plot in Figs. 2 (c) and 2
(d) the ground-state phase diagram for the fixed values $K/|J_{1}|=0.5$ and
$K/|J_{1}|=0.75$, respectively. These figures obviously represent the topology
of the first-order phase transitions between the ground states with different
magnetization magnitude (solid lines with different colors), and between the
ground states $|{1,1,1}\rangle_{\text{I}}$ and $|{1,1,1}\rangle_{\text{II}}$
with the same magnetization value (vertical dotted line) that comprise the
total magnetization at one-half of saturation magnetization.

Table 1:  The phase boundaries between four ground states of the considered spin-1/2 Heisenberg model on a tetranuclear square compound (see Fig. 1) with Hamiltonian (2). First column classifies each pair of ground states with the notation $|S_{\text{T}},S_{12},S_{34}\rangle$ that have mutual boundary. The second column accounts for the eigenvalues of the $S_{\text{T}}^{z}$ operator, while third column stands for the corresponding eigenenergies to the first column. And the last column figures out the relevant boundary curve formula by adopting $|J_{1}|$ as the energy unit. Ground states |  $S_{\text{T}}^{z}$ |  Energies |  Boundary curve
---|---|---|---
$|{1,1,1}\rangle_{\text{I}}$

$|0,1,1\rangle$ |  1

0 |  $E_{11}=-\mu_{B}B+\frac{1}{2}J_{1}+\frac{1}{2}J_{2}-\frac{7}{16}K$

$E_{16}=\frac{1}{2}J_{1}+\frac{1}{2}J_{2}+\frac{5}{16}K-\frac{1}{2}\sqrt{\Lambda}$ |  $\mu_{B}B/|J_{1}|=-\frac{3}{4}K/|J_{1}|+\frac{1}{2|J_{1}|}\sqrt{\Lambda}$
$|{1,1,1}\rangle_{\text{I}}$

$|2,1,1\rangle$ |  1

2 |  $E_{11}=-\mu_{B}B+\frac{1}{2}J_{1}+\frac{1}{2}J_{2}-\frac{7}{16}K$
$E_{14}=-2\mu_{B}B-\frac{1}{2}J_{1}-\frac{1}{2}J_{2}+\frac{1}{16}K$ |  $\mu_{B}B/|J_{1}|=-1-J_{2}/|J_{1}|+\frac{1}{2}K/|J_{1}|$
$|{1,1,1}\rangle_{\text{I}}$

$|{1,1,1}\rangle_{\text{II}}$ |  1

-1 |  $E_{11}=-\mu_{B}B+\frac{1}{2}J_{1}+\frac{1}{2}J_{2}-\frac{7}{16}K$
$E_{14}=-\mu_{B}B+\frac{1}{2}J_{1}-\frac{1}{2}J_{2}-\frac{1}{16}K$ |  $J_{2}/|J_{1}|=\frac{1}{2}K/|J_{1}|$
$|{1,1,1}\rangle_{\text{II}}$

$|2,1,1\rangle$ |  -1

2 |  $E_{11}=-\mu_{B}B+\frac{1}{2}J_{1}-\frac{1}{2}J_{2}-\frac{1}{16}K$
$E_{14}=-2\mu_{B}B-\frac{1}{2}J_{1}-\frac{1}{2}J_{2}+\frac{1}{16}K$ |  $\mu_{B}B/|J_{1}|=1$
$|{1,1,1}\rangle_{\text{II}}$

$|0,1,1\rangle$ |  -1

0 |  $E_{11}=-\mu_{B}B+\frac{1}{2}J_{1}-\frac{1}{2}J_{2}-\frac{1}{16}K$
$E_{16}=\frac{1}{2}J_{1}+\frac{1}{2}J_{2}+\frac{5}{16}K-\frac{1}{2}\sqrt{\Lambda}$ |  $\mu_{B}B/|J_{1}|=-J_{2}/|J_{1}|-\frac{1}{4}K/|J_{1}|+\frac{1}{2|J_{1}|}\sqrt{\Lambda}$

As a result, we observe that the ground-state phase diagram of the model does
not show phase transition between two ground states with the same
magnetization value when the cyclic four-spin exchange term is zero, while for
$K/|J_{1}|>0$ the system exhibits a typical first-order phase transition
between them in the ferromagnetic regime $J_{2}/|J_{1}|>0$ (see dotted line of
Fig. 2 (c)). With increase of the $K/|J_{1}|$ we witness that the boundary
between $|{1,1,1}\rangle_{\text{I}}$ and $|{1,1,1}\rangle_{\text{II}}$ becomes
wider and shifts towards higher values of the coupling constant
$J_{2}/|J_{1}|$.

##  III Results and discussion

In this section, we discuss the novel results obtained from our theoretical
investigations for the effects of the cyclic four-spin interaction on the
geometric $\Pi_{4}$ average of tangles as a measure of the tetrapartite
quantum entanglement. Next, we examine the whole entanglement $\Pi_{4}$,
bipartite negativity $N_{\text{AB}}$ as a measure of pairwise quantum
entanglement, average fidelity $F_{\text{A}}$ and QFI of the tetranuclear
$\mathrm{Cu}^{\mathrm{II}}$ square complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ at
temperatures as high as the room temperature scale and even more. This
compound exceptionally offers an experimental prototype with robust quantum
correlations allowed us to pleasantly elucidate a theoretical insight into the
quantum entanglement and the power of teleportation of a four-qubit quantum
system with strong antiferromagnetic Heisenberg interaction at sufficiently
high temperatures. Some interesting remarks will be delivered when a typical
cyclic four-spin interaction is considered for this exampled model.

###  III.1 Theoretical investigation of the tetrapartite entanglement

For our favorite system, there are different tetrapartite entanglement
quantifiers. Among them, $\pi-$tangle and four-tangle can be employed as the
most generic ones [60, 61, 62, 63]. To extend the definition of the
tetrapartite entanglement we invoke the whole entanglement measure
${\Pi}_{4}$. For a given four-qubit state we begin by defining below
negativities

$\begin{array}[]{lcl}N_{\mathrm{A(BCD)}}=||\rho^{T_{\mathrm{A}}}_{\mathrm{ABCD}}||-1,\\\\[2.84544pt]
N_{\mathrm{AB(CD)}}=||\rho^{T_{\mathrm{AB}}}_{\mathrm{ABCD}}||-1,\\\\[2.84544pt]
N_{\mathrm{A(BC)}}=||\rho^{T_{\mathrm{A}}}_{\mathrm{ABC}}||-1,\\\\[2.84544pt]
N_{\mathrm{AB}}=||\rho^{T_{\mathrm{A}}}_{\mathrm{AB}}||-1,\end{array}$ (4)

in which $\rho_{\text{ABCD}}(T)$ $=\exp\left(-\beta{H}\right)/Z$ denotes the
density matrix of four-spin system with Hamiltonian (2) in thermal equilibrium
at temperature $T$ where $Z=\mathrm{Tr}\left[\exp(-\beta{H})\right]$ denotes
the partition function of the four-spin system and $\beta={1/{{k_{B}}T}}$. The
density operator $\rho_{\text{ABCD}}(T)$ can be thus defined as

$\rho_{\mathrm{ABCD}}(T)=\frac{1}{Z}\sum\limits_{l=1}^{16}\exp\left(-\beta
E_{l}\right)\left|\psi_{l}\right\rangle\left\langle\psi_{l}\right|,$ (5)

where $E_{l}$ and $\left|\psi_{l}\right\rangle$ are, respectively, eigenvalues
and eigenvectors of the Hamiltonian $H$. The nonzero coefficients of the
density matrix $\rho_{\mathrm{ABCD}}$ are given in Appendix A.
${T_{\mathrm{A}}}$ and ${T_{\mathrm{AB}}}$ indicate the partial transpose over
the part $\mathrm{A}$ and pair $\mathrm{AB}$, respectively, and $||\cdot||$
stands for the trace norm of a matrix represented as
$||O||=\mathrm{Tr}{\sqrt{O^{\dagger}O}}$ [64].
$N_{\mathrm{A(BCD)}},N_{\mathrm{AB(CD)}},N_{\mathrm{A(BC)}}$, and
$N_{\mathrm{AB}}$ describe the entanglement between two parts, such as $1-3$
tangle, $2-2$ tangle, $1-2$ tangle, and $1-1$ tangle entanglement. For
instance, $N_{\mathrm{A(BCD)}}$ describes the entanglement between the part
$\mathrm{A}$ and the others $\\{\mathrm{B},\mathrm{C},\mathrm{D}\\}$.
Similarly, $N_{\mathrm{AB(CD)}}$ means the negativity between pair
$\\{\mathrm{A},\mathrm{B}\\}$ and remaining pair
$\\{\mathrm{C},\mathrm{D}\\}$, while $N_{\mathrm{AB}}$ introduces the
bipartite negativity between part $\mathrm{A}$ and part $\mathrm{B}$ where
$\rho_{\mathrm{AB}}=\mathrm{Tr}_{\mathrm{CD}}[\rho_{\mathrm{ABCD}}]$.
$N_{\mathrm{A(BC)}}$ demonstrates the entanglement negativity between part
$\mathrm{A}$ and pair $\\{\mathrm{B},\mathrm{C}\\}$ after tracing over
$\mathrm{D}$. Here, the $1-1$ and $1-3$ tangles satisfy the following
Coffman–Kundu–Wootters (CKW) monogamously inequality relation [20]

$\begin{array}[]{lcl}N^{2}_{\mathrm{A(BCD)}}\geq
N^{2}_{\mathrm{AB}}+N^{2}_{\mathrm{AC}}+N^{2}_{\mathrm{AD}}.\end{array}$ (6)

Accordingly, four-tangle entanglement (residual tangles) can be characterized
as

$\begin{array}[]{lcl}\pi_{\mathrm{A}}=N^{2}_{\mathrm{A(BCD)}}-N^{2}_{\mathrm{AB}}-N^{2}_{\mathrm{AC}}-N^{2}_{\mathrm{AD}},\\\
\pi_{\mathrm{B}}=N^{2}_{\mathrm{B(ACD)}}-N^{2}_{\mathrm{BA}}-N^{2}_{\mathrm{BC}}-N^{2}_{\mathrm{BD}},\\\
\pi_{\mathrm{C}}=N^{2}_{\mathrm{C(ABD)}}-N^{2}_{\mathrm{CA}}-N^{2}_{\mathrm{CB}}-N^{2}_{\mathrm{CD}},\\\
\pi_{\mathrm{D}}=N^{2}_{\mathrm{D(ABC)}}-N^{2}_{\mathrm{DA}}-N^{2}_{\mathrm{DB}}-N^{2}_{\mathrm{DC}}.\end{array}$
(7)

In the one hand, the existence of four different residual tangles
$\pi_{\mathrm{A}}$, $\pi_{\mathrm{B}}$, $\pi_{\mathrm{C}}$ and
$\pi_{\mathrm{D}}$ encourages one to use $\pi_{4}-$tangle;
$\pi_{4}=\frac{1}{4}(\pi_{\mathrm{A}}+\pi_{\mathrm{B}}+\pi_{\mathrm{C}}+\pi_{\mathrm{D}})$
for measuring the whole degree of entanglement of the tetrapartite spin
system. On the other hand, another measure of whole entanglement, described as
geometric mean ${\Pi}_{4}$ [65], can be defined as

$\begin{array}[]{lcl}{\Pi}_{4}=\sqrt[4]{\pi_{\mathrm{A}}\pi_{\mathrm{B}}\pi_{\mathrm{C}}\pi_{\mathrm{D}}}.\end{array}$
(8)

In what follows, we will examine the whole entanglement of the tetrameric
spin-1/2 Heisenberg system (2) by using ${\Pi}_{4}$.

Figure 3: Contour plot of the whole entanglement negativity $\Pi_{4}$ of the
spin-1/2 antiferromagnetic Heisenberg XXX model on a tetranuclear square
compound with Hamiltonian (2) consisting of the cyclic four-spin interaction
$K/|J_{1}|$. Panels (a) and (b) show the whole entanglement negativity
$\Pi_{4}$ in the $\mu_{B}B/|J_{1}|-k_{B}T/|J_{1}|$ plane for $K/|J_{1}|=0$ and
$K/|J_{1}|=0.1$, respectively. (c) Contour plot of $\Pi_{4}$ in the
$\mu_{B}B/|J_{1}|-K/|J_{1}|$ plane at low temperature $k_{B}T/|J_{1}|=0.02$.
(d) The whole entanglement $\Pi_{4}$ in the $k_{B}T/|J_{1}|-K/|J_{1}|$ plane
at $\mu_{B}B/|J_{1}|=0.5$. According to the coordinates of the triple point
appeared in the ground-state phase diagram of the model (see Fig. 2(a)), weak
antiferromagnetic exchange coupling $J_{2}/|J_{1}|=-0.053$ has been assumed.

We represent in Fig. 3 (a) the whole entanglement $\Pi_{4}$ of the four-qubit
state $\rho_{\mathrm{ABCD}}$ corresponds to the spin-1/2 Heisenberg model on a
square-shaped complex in the $\mu_{B}B/|J_{1}|-k_{B}T/|J_{1}|$ plane when
$K/|J_{1}|=0$. A weak antiferromagnetic exchange coupling
$J_{2}/|J_{1}|=-0.053$ is considered that coincides the coordinates of the
triple point at which three ground states coexist together (see Fig. 2(a)). It
is crystal clear that the tetrapartite entanglement $\Pi_{4}$ is in its
maximum value $\Pi_{4}\approx 0.621$ at sufficiently low temperature
$k_{B}T/|J_{1}|<0.01$ when the magnetic field is as large as the absolute
value $|J_{1}|$. In Fig. 3 (b) we show the whole entanglement in the
$\mu_{B}B/|J_{1}|-k_{B}T/|J_{1}|$ plane when a small nonzero value
$K/|J_{1}|=0.1$ is considered. One sees that a nonzero positive value of the
cyclic four-spin exchange $K/|J_{1}|$ remarkably broadens the region at which
the tetrapartite entanglement exists. Figure 3 (c) displays the tetrapartite
entanglement negativity in the $\mu_{B}B/|J_{1}|-K/|J_{1}|$ plane at low
temperature $k_{B}T/|J_{1}|=0.02$. This figure practically reproduces the
ground-state phase diagram plotted in Fig. 2 (a). The effect of cyclic four-
spin exchange on the whole entanglement negativity $\Pi_{4}$ is obvious in
this figure. When the ring exchange $K/|J_{1}|$ increases the maximum value of
the whole entanglement can be achieved at lower magnetic fields. We display in
Fig. 3 (d) the function $\Pi_{4}$ in the $k_{B}T/|J_{1}|-K/|J_{1}|$ plane when
the magnetic field is fixed at $\mu_{B}B/|J_{1}|=0.5$. From this figure one
finds that with increase of the temperature the both entangled states
$|{0,1,1}\rangle$ and $|{1,1,1}\rangle_{\text{I}}$ lose their entanglement
degree. The phase boundary between these states are evident by a narrow pale
stripe. This boundary separates two phases even at high enough temperature
region ($k_{B}T/|J_{1}|\approx 0.5$). It is also interesting to notify that
this separation starts from $K/|J_{1}|\approx\mu_{B}B/|J_{1}|$. It could be
concluded that the whole entanglement in the model under consideration is
still sizable at sufficiently high temperature when a special cyclic four-spin
interaction comes into play.

###  III.2 Experimental testing of theoretical results

From now on, we concentrate on the particular case
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$
that showcasts a good experimental platform of a strong antiferromagnetic
Heisenberg XXX model on a square-shaped compound. Then, we theoretically
investigate the effects of a deemed cyclic four-spin interaction on the whole
entanglement negativity, bipartite entanglement negativity, the ability of
quantum teleportation and QFI of this experimental example by considering
actual values of the exchange couplings $J_{1}/k_{B}=-638\;\text{K}$ and
$J_{2}/k_{B}=-34\;\text{K}$ reported in Ref. [16].

The spin density of the complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ is
mainly distributed between the $d_{x^{2}-y^{2}}$ orbital of the copper atoms
and the hybrid orbitals of the donor atoms in the equatorial plane. Hence, it
is expected that this compound does not consists of actual anisotropy neither
two-ion exchange anisotropy nor Dzyaloshinsky–Morya (DM) interaction due to
its geometry. For that reason, we proceeded to theoretically examine the
effect size value of the cyclic four-spin interaction $K/k_{B}$ on the
temperature dependence of $\chi_{\text{M}}T$ product of the complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$,
and we eventually realized that it is much smaller than other anisotropy
parameters specially for the temperature range $T<100\ \text{K}$. The results
concerning these findings have not been added to this work for brevity.
Hopefully, the theoretical studies reported in the current work will stimulate
future experimental testing of the entanglement negativity, quantum
teleportation and QFI of the tetranuclear square complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ at
high temperature through modern experimental procedures.

####  III.2.1 Tetrapartite entanglement

Figure 4:  Tetrapartite entanglement $\Pi_{4}$ of the complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ as
a function of the assumed cyclic four-spin interaction $K/k_{B}$ for a few
selected values of the magnetic field at moderate temperature $T=10\;\text{K}$
such that $J_{1}/k_{B}=-638\;\text{K}$, $J_{2}/k_{B}=-34\;\text{K}$.

The cyclic four-spin exchange dependence of the geometric $\Pi_{4}$ average is
manifested in Fig. 4 at moderate temperature $T=10\;\text{K}$ for different
values of the magnetic field where exchange couplings are taken up as actual
values $J_{1}/k_{B}=-638\;\text{K}$ and $J_{2}/k_{B}=-34\;\text{K}$. It is
evident from Fig. 4 that for the magnetic fields smaller than absolute value
$|J_{1}|$ the whole entanglement $\Pi_{4}$ of the tetranuclear
Cu${}^{\text{II}}_{4}$ complex exhibits an unconventional minimum at
$K/k_{B}=68\;\text{K}$ such that $\Pi_{4}=0$. It is worthwhile to note that a
very small value of the cyclic four-spin interaction ($K/k_{B}\ll|J_{1}|$) can
still alters the degree of the whole entanglement negativity $\Pi_{4}$ for
$B\ll|J_{1}|$ (see blue solid line of Fig. 4). This finding is in accordance
with that of discussed in Fig. 3(b). When $K$ increases further than the
observed minimum point ($K/k_{B}=68\;\text{K}$), the entanglement degree of
the system is significantly enhanced. We choose the intriguing value
$K/k_{B}=68\;\text{K}$ then examine the bipartite entanglement negativity
$N_{\mathrm{AB}}$ in what follows. In addition, a reentrant behavior can be
observed in $\Pi_{4}$ curve close to the critical points at which phase
transition occurs between $|0,1,1\rangle$ and $|1,1,1\rangle_{\text{I}}$
states (see Figs. 2 (a) and 3(c)). On the other hand, when magnetic fields
greater than the absolute value $|J_{1}|$ are considered the whole
entanglement $\Pi_{4}$ starts to increase from a critical point and rapidly
reaches its maximum value.

####  III.2.2 Bipartite entanglement

In this part we will put forward a theoretical inspection for the bipartite
entanglement between spins A and B of the complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$
with reduced density matrix $\rho_{\mathrm{AB}}$ (see Eq. (9) and Appendix A).
Since the system with Hamiltonian (2) is symmetric, all of the reduced density
matrices have the same formation, i.e., they are block diagonal. Therefore, In
the standard basis
$\\{\left|00\right\rangle,\left|01\right\rangle,\left|10\right\rangle,\left|11\right\rangle\\}$
the reduced density matrix
$\rho_{\mathrm{AB}}(T)=\mathrm{Tr}_{\mathrm{CD}}\left[\rho_{\mathrm{ABCD}}(T)\right]$
is given by

$\rho_{\mathrm{AB}}(T)=\begin{pmatrix}\varrho_{11}&0&0&0\\\
0&\varrho_{22}&\varrho_{23}&0\\\ 0&\varrho^{*}_{23}&\varrho_{33}&0\\\
0&0&0&\varrho_{44}\end{pmatrix}$ (9)

where $\varrho^{*}_{23}=\varrho_{23}$. Nonzero elements of the above density
operator are directly connected to the coefficients of the total density
operator $\rho_{\mathrm{ABCD}}(T)$ as we have composed in Appendix A. In order
to investigate bipartite thermal entanglement of the reduced density matrix
$\rho_{\mathrm{AB}}$, we employ negativity $N_{\mathrm{AB}}$ (4) that is one
the most widely accepted measures of entanglement for a two-qubit system. The
case $N_{\mathrm{AB}}=0$ happens when the system state is separable, whereas
$N_{\mathrm{AB}}=1$ indicates maximally entangled state.

Figure 5: (a) Bipartite entanglement negativity $N_{\mathrm{AB}}$ versus
temperature for several fixed values of the magnetic field when $K=0$, taking
$J_{1}/k_{B}=-638\;\text{K}$ and $J_{2}/k_{B}=-34\;\text{K}$. (b) Negativity
$N_{\mathrm{AB}}$ as a function of magnetic field for several selected
temperatures and the same set of exchange couplings to panel (a).

We have shown in Fig. 5 (a) the bipartite negativity $N_{\mathrm{AB}}$ as a
function of the temperature $T$, for several values of the magnetic field $B$
at fixed exchange interactions $J_{1}/k_{B}=-638\;\text{K}$ and
$J_{2}/k_{B}=-34\;\text{K}$, where the cyclic four-spin interaction $K$ is
zero. It is clear that the negativity $N_{\mathrm{AB}}$ decreases from its
maximum value $N_{\mathrm{AB}}=1$ when $T$ and $B$ increase. We find maximum
entanglement for notably high temperature ($T\approx 100\;\text{K}$) and high
magnetic fields ($B<300\;\text{T}$). This quantity remains robustly alive for
temperatures higher than room temperature. Interestingly, there is an
entanglement death at a critical temperature $T_{\text{c}}\approx
582\;\text{K}$ that is comparable with the absolute value of stronger
antiferromagnetic coupling constant $|J_{1}|$ between spins A and B.

The bipartite entanglement negativity $N_{\mathrm{AB}}$ as a function of the
magnetic field is displayed in Fig. 5 (b) for several temperatures and the
same exchange couplings to panel Fig. 5 (a). One observes maximum pairwise
entanglement at low temperature and remarkably high magnetic fields with
comparable magnitude to $|J_{1}|$. Close to $B=|J_{1}|$ there is a sharp
decrease in $N_{\mathrm{AB}}$ curve until a narrow plateau is observed at
$N_{\mathrm{AB}}\approx 0.207$, then it decreases and rapidly vanishes at a
critical magnetic field $B>|J_{1}|$. It is evident that nearby room
temperature the pairwise entanglement $N_{\mathrm{AB}}$ in this model is
sizable even for strong magnetic fields (see solid gray line and diamonds).

By inspecting Fig. 6 (a) we see that when a nonzero value of the cyclic four-
spin term, $K/k_{B}=68\;\text{K}$, is assumed the entanglement degree
decreases for the magnetic fields as large as $|J_{1}|$ (see dotted line).
Figure 6 (b) illustrates the field dependence of the bipartite negativity
$N_{\mathrm{AB}}$ for various selected values of the temperature and fixed
$K/k_{B}=68\;\text{K}$ that is very small compared with the $|J_{1}|$. By
Comparing solid gray lines and diamonds plotted in Figs. 5 (b) and 6 (b), one
quickly realizes that the interplay of even a small cyclic four-spin
interaction in the model leads to weaken the pairwise entanglement between Cu
atoms interacted through the stronger antiferromagnetic exchange coupling
$J_{1}$ nearby the room temperature. Also, the middle step appeared in the
negativity curve becomes wider when a typical cyclic exchange interaction
$K>0$ is considered. Referring to Fig. 2 (a), one sees that emerged step is
accompanied with the ground state $|1,1,1\rangle_{\text{I}}$ and the change in
the $N_{\mathrm{AB}}$ behavior is in a good agreement with the ground-state
phase transition of the model in $B-K$ plane. When the bipartite entanglement
$N_{\mathrm{AB}}$ is maximum, the system’s state is entangled state
$|0,1,1\rangle$, whereas for the case $N_{\mathrm{AB}}=0$ the system is in
fully polarized state $|2,1,1\rangle$.

Figure 6: (a) Bipartite entanglement negativity $N_{\mathrm{AB}}$ versus
temperature for several fixed values of the magnetic field, taking
$J_{1}/k_{B}=-638\;\text{K}$ and $J_{2}/k_{B}=-34\;\mathrm{K}$, when
$K/k_{B}=68\;\mathrm{K}$. (b) Magnetic field dependence of the bipartite
entanglement $N_{\mathrm{AB}}$ for various fixed values of the temperature,
assuming the same set of exchange couplings to panel (a).

The negativity $N_{\mathrm{AB}}$ behavior versus the ring exchange $K$ is
depicted in Fig. 7 for a few selected values of the magnetic field at
$T=10\;\text{K}$. Principally, this quantity behaves different from the
$\Pi_{4}$ against the cyclic four-spin exchange $K$. As can be seen in this
figure, $N_{\mathrm{AB}}$ remains almost steady at value
$N_{\mathrm{AB}}\approx 0.207$ for the magnetic field $B\approx|J_{1}|$ (see
blue solid line). For magnetic fields smaller than this value, when $K$
increases monotonically the negativity $N_{\mathrm{AB}}$ decreases from its
maximum value and drops down to approximately one-third of maximum value of
the whole entanglement $\Pi_{4}^{\mathrm{max}}\approx 0.621$, namely
$N_{\mathrm{AB}}\approx 0.207$, at a specific point that its position on the
$K-$axis can be controlled by the strength of the magnetic field. This
decrease is accompanied with the ground-state phase transition from
$|0,1,1\rangle$ to $|1,1,1\rangle_{\text{I}}$.

Figure 7:  Bipartite entanglement negativity as a function of the cyclic four-
spin interaction $K$ for a few selected values of the magnetic field and fixed
temperature $T=10\;\text{K}$ such that $J_{1}/k_{B}=-638\;\text{K}$,
$J_{2}/k_{B}=-34\;\text{K}$.

For the magnetic field larger that $|J_{1}|$, firstly, the negativity
$N_{\mathrm{AB}}$ arises from zero and then continues to increase upon the
increase of $K$ until reaches fixed value $N_{\mathrm{AB}}\approx 0.207$.
Under this conditions the ground-state phase transition from $|2,1,1\rangle$
to $|1,1,1\rangle_{\text{I}}$ happens. Here, any change in the behavior of
bipartite entanglement negativity reminisces us the position of triple point
at which three ground states $|0,1,1\rangle$, $|1,1,1\rangle_{\text{I}}$ and
$|2,1,1\rangle$ coexist together (see Fig. 2 (a)). In result, we understand
that the cyclic four-spin exchange (even for $K\ll|J_{1}|$) plays an important
role in the bipartite entanglement degree of the
$\mathrm{Cu}_{4}^{\mathrm{II}}$ compound, too.

####  III.2.3 Quantum teleportation and quantum Fisher information

Quantum teleportation is interpreted as transfer of quantum states from one
place to another. The teleportating of a quantum bit (qubit) is achieved using
quantum entanglement, in which two or more particles are indistinguishably
connected to each other. If an entangled pair of particles are located in two
spatially separate places, no matter the distance between them, the encoded
quantum information can be teleported through a channel. The teleportation is
the cornerstone of quantum networks, data storage, precision sensing and
computing for ushering in a new era of high-secure communication. In this
section, we study the quantum teleportation throughout a couple of
tetranuclear
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$
complexes as an optimal resource of a generalized depolarizing channels [66].
In fact, the main goal of this section is to investigate the effects of cyclic
four-spin exchange and the magnetic field on the possibility of teleportation
through a couple of independent tetrameric Cu${}^{\text{II}}_{4}$ systems. Let
us assume the input state being an arbitrary unknown two-qubit pure state
$\left|\psi_{\text{in}}\right\rangle=\cos(\frac{\theta}{2})\left|10\right\rangle+e^{i\phi}\sin(\frac{\theta}{2})\left|01\right\rangle$,
where $0\leq\theta\leq\pi$ and $0\leq\phi\leq 2\pi$. Angle $\theta$ describes
an arbitrary state with specific amplitude and $\phi$ is the phase corresponds
to this state. The concurrence of the input state
$\rho_{\text{in}}=\left|\psi_{\text{in}}\right\rangle\left\langle\psi_{\text{in}}\right|$
can be written as

${C}_{\text{in}}=2|e^{i\phi}\sin(\frac{\theta}{2})\cos(\frac{\theta}{2})|=|\sin(\theta)|.$
(10)

Owing to the fact that quantum channels are known as completely positive and
trace-preserving (CPTP) operators, an input density operator is mapped to an
output one [66].

Figure 8: The schematic representation for the teleportation process of the
input state $\rho_{\mathrm{in}}$ through a couple of discrete tetrameric
Cu${}^{\text{II}}_{4}$ unit cells as quantum channels with density operator
$\rho_{\mathrm{ch}}$. The teleported output state is demonstrated by
${\rho}_{\mathrm{out}}$.

Generally, when the quantum state is teleported through the mixed channel
${\rho}_{\text{ch}}$, the output replica state ${\rho}_{\text{out}}$ can be
obtained by applying joint measurements and the local unitary transformations
on the input state $\rho_{\text{in}}$, hence

${\rho}_{\text{out}}=\sum_{i,j\in\\{0,x,y,z\\}}p_{i}p_{j}(\sigma^{i}\otimes\sigma^{j})\rho_{\text{in}}(\sigma^{i}\otimes\sigma^{j}),$
(11)

where $\sigma_{0}=I$ with $I$ being the $2\times 2$ identity matrix.
$p_{i}=\mathrm{Tr}[E^{i}\mathcal{\rho}_{\text{ch}}]$ satisfies the condition
$\sum_{i}p_{i}=1$.
$E^{0}=\left|\Psi^{-}\right\rangle\left\langle\Psi^{-}\right|,$
$E^{1}=\left|\Phi^{-}\right\rangle\left\langle\Phi^{-}\right|,$
$E^{2}=\left|\Psi^{+}\right\rangle\left\langle\Psi^{+}\right|,$ and
$E^{3}=\left|\Phi^{+}\right\rangle\left\langle\Phi^{+}\right|,$ from which
$\\{\left|\Psi^{\pm}\right\rangle,$ $\left|\Phi^{\pm}\right\rangle\\}$ are
Bell states. In this paper, we consider the density operator channel as
${\rho}_{\text{ch}}=\rho_{\mathrm{AB}}(T)$ (see Appendix A). Therefore, the
output density operator ${\rho}_{\text{out}}$ takes the form

${\rho}_{\text{out}}=\begin{pmatrix}\widetilde{{\rho}}_{11}&0&0&0\\\
0&\widetilde{{\rho}}_{22}&\widetilde{{\rho}}_{23}&0\\\
0&\widetilde{{\rho}}_{23}&\widetilde{{\rho}}_{33}&0\\\
0&0&0&\widetilde{{\rho}}_{11}\end{pmatrix},$ (12)

where

$\displaystyle\widetilde{{\rho}}_{11}$ $\displaystyle=$
$\displaystyle({\varrho}_{11}+{\varrho}_{44})({\varrho}_{22}+{\varrho}_{33}),$
$\displaystyle\widetilde{{\rho}}_{22}$ $\displaystyle=$
$\displaystyle({\varrho}_{11}+\varrho_{44})^{2}\cos^{2}(\frac{\theta}{2})+({\varrho}_{22}+{\varrho}_{33})^{2}\sin^{2}(\frac{\theta}{2}),$
$\displaystyle\widetilde{{\rho}}_{33}$ $\displaystyle=$
$\displaystyle({\varrho}_{11}+{\varrho}_{44})^{2}\sin^{2}(\frac{\theta}{2})+({\varrho}_{22}+{\varrho}_{33})^{2}\cos^{2}(\frac{\theta}{2}),$
$\displaystyle\widetilde{{\rho}}_{23}$ $\displaystyle=$ $\displaystyle
2e^{i\phi}|{\varrho}_{23}|^{2}\sin\theta.$

Figure 9: Average fidelity ${F}_{\text{A}}$ versus temperature for parameter
set $J_{1}/k_{\text{B}}=-638\;\text{K}$, $J_{2}/k_{\text{B}}=-34\;\text{K}$
and fixed (a) $K=0$, (b) $K/k_{\text{B}}=68\;\text{K}$, where several values
of the magnetic field are considered.

To describe the quality of teleportation process, it is useful to study the
fidelity between $\rho_{\text{in}}$ and ${\rho}_{\text{out}}$ to characterize
the state that is teleported through a channel. When the input state is a pure
state, one can apply the concept of fidelity as a tool to prob the
teleportation performance of a quantum channel [67]. The fidelity is defined
by

${F}=\left(\mathrm{Tr}\sqrt{\sqrt{\rho_{\text{in}}}\mathcal{\rho}_{\text{out}}\sqrt{\rho_{\text{in}}}}\right)^{2},$
(13)

thus in our case $F$ is deduced from the straightforward calculation

${F}=\frac{\sin^{2}\theta}{2}[({\varrho}_{11}+{\varrho}_{44})^{2}+4{\varrho}_{23}^{2}-4{\varrho}_{22}^{2}]+4{\varrho}_{22}^{2}.$
(14)

The average fidelity $\mathcal{F}_{A}$ can be calculated from the formula

${F}_{\mathrm{A}}=\frac{1}{4\pi}\int_{0}^{2\pi}d\phi\int_{0}^{\pi}d\theta{F}\sin\theta.$
(15)

After integrating Eq. (15), the average fidelity ${F}_{\text{A}}$ for the
favorite bipartite subsystem can be given by

${F}_{\text{A}}=\frac{1}{3}[({\varrho}_{11}+{\varrho}_{44})^{2}+4{\varrho}_{23}^{2}-({\varrho}_{22}+{\varrho}_{33})^{2}]+({\varrho}_{22}+{\varrho}_{33})^{2}.$
(16)

By using quantum systems, scientists are seriously looking for a quantum state
transmitting better than the classical communication protocols. This
superiority is achieved if and only if ${F}_{\text{A}}>2/3$ (horizontal dotted
lines plotted in Fig. 9), which ${F}_{\text{A}}=2/3$ is the best fidelity
among the classical channels. Let us proceed to examine the ability of system
under study in teleporting information when the experimental data are
considered. In Fig. 9 (a) the average fidelity $F_{\mathrm{A}}$ is shown as a
function of the temperature for various magnetic fields such that exchange
couplings are being $J_{1}/k_{B}=-638\;\text{K}$, $J_{2}/k_{B}=-34\;\text{K}$
when $K=0$. It can be comprehended from this figure that quantity
${F}_{\mathrm{A}}$ does not reaches the limit of quantum fidelities for
$B\gtrsim 608\;\text{T}$ that is comparable with the absolute value $|J_{1}|$,
hence the teleportation of information happens for the magnetic field range
$B<608\;\text{T}$. By looking at Fig. 2 (b), one understands that the quantum
teleportation is solely possible for the area below than blue line at which a
ground-state phase transition occurs between $|1,1,1\rangle_{\text{I}}$ and
$|0,1,1\rangle$ states. The possibility of the teleportation decreases upon
increase of the temperature.

In panel Fig. 9 (b) we display the average fidelity versus the temperature for
several magnetic fields and the same coupling constants values
$J_{1}/k_{B}=-638\;\text{K}$, $J_{2}/k_{B}=-34\;\text{K}$ when the cyclic
four-spin exchange term has nonzero value $K/k_{B}=68\;\text{K}$ (at which
$\Pi_{4}=0$). Of course, for the magnetic fields remarkably lower than
absolute value $|J_{1}|$ the fidelity is in its maximum value
${F}_{\mathrm{A}}=1$ which monotonically decreases upon heating and tends to
the constant value ${F}_{\mathrm{A}}=1/3$. Considering the nonzero value of
$K$ results in decreasing the ability of teleportation for magnetic fields
close to $|J_{1}|$. In fact the critical magnetic field which determines the
frontier between classical transmitting and quantum teleportation of
information through this channel moves towards lower magnetic fields (see
long-dashed purple lines of Figs. 9 (a) and (b)).

Figure 10:  (a) Field dependence of the average fidelity ${F}_{\text{A}}$ for
various temperatures, where other parameters have been taken as Fig. 9 (a).
(b) The same function at $K/k_{\text{B}}=68\;\text{K}$ against the magnetic
field when different fixed values of the temperature are assumed. Horizontal
dotted lines denote the limit value 2/3 and identify the quantum teleportation
becomes possible for the case when the average fidelity is greater than this
value. Circles and diamonds illustrated in panels (a) and (b) represent the
QFI versus the magnetic field at high temperature $T=100\;\text{K}$ and
$T=200\;\text{K}$ for the same set of other parameters to the corresponding
$F_{\text{A}}$.

To uncover the influence of the cyclic four-spin interaction on the magnetic
field dependencies of average fidelity, a few plots of $F_{\text{A}}$ versus
the magnetic field are displayed in Figs. 10 (a) and 10 (b) for several
temperatures by assuming $K=0$ and $K/k_{B}=68\;\text{K}$, respectively. The
first important thing that can be viewed in these figures is the existence of
a sharp reentrance point in the average fidelity curve. This phenomenon
reminisces us the ground-state phase transition from $|2,1,1\rangle$ to
$|1,1,1\rangle_{\text{I}}$. At low temperatures this quantity remains maximum
(${F}_{\mathrm{A}}=1$) for notably high magnetic fields. Conventionally when
the temperature increases the average fidelity decreases and the possibility
of teleportation through this model is restricted. The interesting point to
emphasize is that, for the range $B<68\;\text{T}$ the average fidelity proves
that the system finds the possibility of quantum teleportation close to the
room temperature (see long-dashed purple line in Fig. 10 (a)). With increase
of the cyclic four-spin interaction $K$ (see Fig. 10 (b)) the average fidelity
decreases and tends to classical limit 2/3. This reduction is more perspicuous
close to the room temperature.

Finally, we examine the effect of ring exchange $K$ on the QFI. For
abbreviating the context, we have explained explicitly the QFI and its
corresponding formulae for a bipartite spin-1/2 system in Appendix B. In Figs.
10 (a) and 10 (b) the QFI alterations are demonstrated versus the magnetic
field for, respectively, $K=0$ and $K/k_{B}=68\;\text{K}$ at high temperatures
$T=100\;\text{K}$ (circles) and $T=200\;\text{K}$ (diamonds). Exchange
interactions have been considered antiferromagnetic as
$J_{1}/k_{B}=-638\;\text{K}$ and $J_{2}/k_{B}=-34\;\text{K}$. This quantity
gradually decreases from its maximum value with increase of the magnetic field
and reaches a global minimum nearby the critical field that the average
fidelity sharply decreases, then increases with further increase of the
magnetic field. This behaviour guarantees the ground-state phase transition
from $|2,1,1\rangle$ to $|1,1,1\rangle_{\text{I}}$. The change of term $K$
leads to change the magnetic position of the minimum appeared in QFI curve.

##  IV Conclusions

In this work, we have considered a cyclic four-spin interaction for a spin-1/2
tetranuclear square-shaped complex whose nearest-neighbor sites interact with
each other through the antiferromagnetic Heisenberg XXX exchange couplings
$J_{1}$ and $J_{2}$. An external magnetic field in the $z-$direction has been
hypothesized, as well. Next, we have investigated the ground-state phase
diagram of the model. The model exhibited four phases in its ground-state
phase diagram, such a way that several discontinuous phase transitions
together with a special triple point at which three different phases coexist
together were observed.

Next, we have theoretically investigated the effects of the cyclic four-spin
interaction on the whole entanglement degree of the spin-1/2 tetranuclear
square-shaped model. We found that the cyclic four-spin interaction could
notably enhance the whole entanglement negativity of the model. Afterwards, we
have paid attention to the special case tetranuclear square complex
$[\text{Cu}_{4}\text{L}_{4}(\text{H}_{2}\text{O})_{4}](\text{ClO}_{4})_{4}$ as
an experimental representative of the spin-1/2 antiferromagnetic Heisenberg
XXX model on a tetranuclear square compound [16]. This complex is symmetric
and possesses two different antiferromagnetic interactions among the copper
ions such that one of them is much more stronger than another one. Then, we
have theoretically examined for this system the whole entanglement intensity
and bipartite entanglement between pairs of CuII ions with stronger exchange
interaction. We have observed that considering a very small cyclic four-spin
interaction may change the size of whole entanglement in this complex at
moderate temperatures. We have also realized that for the magnetic fields
lower than the absolute value of the stronger coupling constant, the whole
entanglement of the system shows an unconventional minimum and a reentrant
behavior at higher values of the cyclic four-spin exchange where a ground-
state phase transition occurs. Increase of the cyclic four-spin term further
than minimum point remarkably enhances the tetrapartite entanglement degree,
where it does not vanish even at sufficiently high temperature.

Bipartite entanglement between pairs of CuII ions interacted through stronger
exchange interaction remains in its maximum value for high temperatures and
high magnetic fields. We understood that the bipartite entanglement remains
robust at sufficiently high temperatures as large as the room temperature
scale. At high temperature, bipartite entanglement is fragile against a cyclic
four-spin interaction at high magnetic field. Indeed, increasing the cyclic
four-spin exchange leads to decrease of the bipartite entanglement degree.

In addition, we have examined the possibility of teleportation through a
couple of discrete tetranuclear copper $\text{Cu}_{4}^{\text{II}}$ complexes.
We discovered that there is a critical magnetic field ($B\lesssim|J_{1}|$) at
which the average fidelity goes further than classical limit 2/3 and the
system gets the power of quantum teleportation of information. At weak
magnetic field the model preserves its ability of teleportation at high enough
temperatures. With increase of the cyclic four-spin term the average fidelity
decays and the quantum teleportation is achievable at lower magnetic fields.

Finally, the QFI has been studied for this system by considering the actual
values of the two-body exchange couplings. We have demonstrated that the QFI
decreases from its maximum value upon increasing the magnetic field and
reaches its minimum value close to the critical conditions at which the
ground-state phase transition occurs between separable state and an entangled
state. The increase of the cyclic four-spin interaction results in decaying
the QFI, as well.

## Acknowledgments

H. Arian Zad acknowledges the receipt of the grant from the Abdus Salam
International Centre for Theoretical Physics (ICTP), Trieste, Italy. N.
Ananikian has been supported by the grant from CSMES RA in the frame of the
research project No. SCS 21AG-1C0068. H. Arian Zad and A. Zoshki acknowledge
the receipt of the grant from the National scholarship Programme of the Slovak
Republic (SAIA). M. Jaščur has been supported by the grant APVV-16-0186. The
authors are also grateful to Prof. Jozef Strečka for his insightful
discussions.

##  Appendix A Density matrices

The partition function of the tetranuclear $\mathrm{Cu}_{4}^{\mathrm{II}}$
square complex on the spin-1/2 XXX Heisenberg model containing two coupling
constants $J_{1}$ and $J_{2}$ and a cyclic four-spin exchange $K$ in the
presence of an external magnetic field $B$ can be written as

$\begin{array}[]{lcl}Z=2e^{-\frac{1}{4}\beta(\frac{1}{4}K-2(J_{1}+J_{2}))}\cosh{(2\beta\mu_{B}gB)}+\\\\[2.84544pt]
2\big{[}2e^{-\frac{1}{16}\beta
K}\cosh(\frac{1}{2}\beta(J_{1}-J_{2}))+e^{-\frac{1}{4}\beta(-\frac{7}{4}K+2(J_{1}+J_{2}))}+\\\\[2.84544pt]
e^{-\frac{1}{4}\beta(\frac{1}{4}K-2(J_{1}+J_{2}))}\big{]}\cosh{(\beta\mu_{B}gB)}+\\\\[2.84544pt]
2e^{-\frac{1}{4}\beta(\frac{5}{4}K+2(J_{1}+J_{2}))}\cosh{(\beta\Lambda)}+\\\\[2.84544pt]
2e^{-\frac{1}{4}\beta(\frac{1}{4}K-2{J_{2}})}\cosh{(\frac{1}{2}\beta{J_{1}})}+e^{-\beta(\frac{1}{16}K-\frac{1}{2}J_{1}+\frac{1}{2}J_{2})}+\\\\[2.84544pt]
e^{-\beta(-\frac{7}{16}K+\frac{1}{2}J_{1}+\frac{1}{2}J_{2})}.\\\\[2.84544pt]
\end{array}$

Function $\Lambda$ has already been introduced in the Sec. II. The nonzero
coefficients of the tetramer density matrix $\rho_{\mathrm{ABCD}}$ are

$\begin{array}[]{lcl}\rho_{1,1},\rho_{2,2},\rho_{3,3},\rho_{4,4},\rho_{5,5},\rho_{6,6},\rho_{7,7},\rho_{8,8},\rho_{9,9},\rho_{10,10},\\\
\rho_{11,11},\rho_{12,12},\rho_{13,13},\rho_{14,14},\rho_{15,15},\rho_{16,16},\\\
\rho_{2,3},\rho_{2,5},\rho_{2,9},\rho_{3,5},\rho_{3,9},\rho_{4,6},\rho_{4,7},\rho_{4,10},\rho_{4,11},\rho_{4,13},\\\
\rho_{5,9},\rho_{6,7},\rho_{6,10},\rho_{6,11},\rho_{6,13},\rho_{7,10},\rho_{7,11},\rho_{7,13},\rho_{8,12},\\\
\rho_{8,14},\rho_{8,15},\rho_{10,11},\rho_{10,13},\rho_{11,13},\rho_{12,14},\rho_{12,15},\rho_{14,15}\end{array}$

where, $\rho_{i,j}$ are long polynomials of parameters $J_{1}$, $J_{2}$, $K$
and $B$. The nonzero coefficients of the density matrix
$\varrho^{\mathrm{A}}=\rho^{T_{\mathrm{A}}}_{\mathrm{ABCD}}$ are

$\begin{array}[]{lcl}\varrho^{\mathrm{A}}_{i,i}=\rho_{i,i}\quad i=1\cdots
16,\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{2,3}=\rho_{2,3},\varrho^{\mathrm{A}}_{2,5}=\rho_{2,5},\varrho^{\mathrm{A}}_{1,10}=\rho_{2,9},\varrho^{\mathrm{A}}_{3,5}=\rho_{3,5},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{1,11}=\rho_{3,9},\varrho^{\mathrm{A}}_{4,6}=\rho_{4,6},\varrho^{\mathrm{A}}_{4,7}=\rho_{4,7},\varrho^{\mathrm{A}}_{2,12}=\rho_{4,10},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{3,12}=\rho_{4,11},\varrho^{\mathrm{A}}_{5,12}=\rho_{4,13},\varrho^{\mathrm{A}}_{1,13}=\rho_{5,9},\varrho^{\mathrm{A}}_{6,7}=\rho_{6,7},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{2,14}=\rho_{6,10},\varrho^{\mathrm{A}}_{3,14}=\rho_{6,11},\varrho^{\mathrm{A}}_{5,14}=\rho_{6,13},\varrho^{\mathrm{A}}_{2,15}=\rho_{7,10},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{3,15}=\rho_{7,11},\varrho^{\mathrm{A}}_{5,15}=\rho_{7,13},\varrho^{\mathrm{A}}_{4,16}=\rho_{8,12},\varrho^{\mathrm{A}}_{6,16}=\rho_{8,14},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{7,16}=\rho_{8,15},\varrho^{\mathrm{A}}_{10,11}=\rho_{10,11},\varrho^{\mathrm{A}}_{10,13}=\rho_{10,13},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{11,13}=\rho_{11,13},\varrho^{\mathrm{A}}_{12,14}=\rho_{12,14},\varrho^{\mathrm{A}}_{12,15}=\rho_{12,15},\\\\[2.84544pt]
\varrho^{\mathrm{A}}_{14,15}=\rho_{14,15}\end{array}$

The nonzero coefficients of the density matrix
$\varrho^{\mathrm{B}}=\rho^{T_{\mathrm{B}}}_{\mathrm{ABCD}}$ are

$\begin{array}[]{lcl}\varrho^{\mathrm{B}}_{i,i}=\rho_{i,i}\quad i=1\cdots
16,\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{2,3}=\rho_{2,3},\varrho^{\mathrm{B}}_{1,6}=\rho_{2,5},\varrho^{\mathrm{B}}_{2,9}=\rho_{2,9},\varrho^{\mathrm{B}}_{1,7}=\rho_{3,5},\varrho^{\mathrm{B}}_{3,9}=\rho_{3,9},\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{4,6}=\rho_{4,6},\varrho^{\mathrm{B}}_{4,7}=\rho_{4,7},\varrho^{\mathrm{B}}_{4,10}=\rho_{4,10},\varrho^{\mathrm{B}}_{4,11}=\rho_{4,11},\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{4,13}=\rho_{4,13},\varrho^{\mathrm{B}}_{1,13}=\rho_{5,9},\varrho^{\mathrm{B}}_{6,7}=\rho_{6,7},\varrho^{\mathrm{B}}_{6,10}=\rho_{6,10},\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{6,11}=\rho_{6,11},\varrho^{\mathrm{B}}_{6,13}=\rho_{6,13},\varrho^{\mathrm{B}}_{7,10}=\rho_{7,10},\varrho^{\mathrm{B}}_{7,11}=\rho_{7,11},\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{7,13}=\rho_{7,13},\varrho^{\mathrm{B}}_{8,12}=\rho_{8,12},\varrho^{\mathrm{B}}_{8,14}=\rho_{8,14},\varrho^{\mathrm{B}}_{8,15}=\rho_{8,15},\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{10,11}=\rho_{10,11},\varrho^{\mathrm{B}}_{9,14}=\rho_{10,13},\varrho^{\mathrm{B}}_{9,15}=\rho_{11,13},\varrho^{\mathrm{B}}_{10,16}=\rho_{12,14},\\\\[2.84544pt]
\varrho^{\mathrm{B}}_{11,16}=\rho_{12,15},\varrho^{\mathrm{B}}_{14,15}=\rho_{14,15}\end{array}$

The nonzero coefficients of the density matrix
$\varrho^{\mathrm{C}}=\rho^{T_{\mathrm{C}}}_{\mathrm{ABCD}}$ are

$\begin{array}[]{lcl}\varrho^{\mathrm{C}}_{i,i}=\rho_{i,i}\quad i=1\cdots
16,\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{1,4}=\rho_{2,3},\varrho^{\mathrm{C}}_{2,5}=\rho_{2,5},\varrho^{\mathrm{C}}_{2,9}=\rho_{2,9},\varrho^{\mathrm{C}}_{1,7}=\rho_{3,5},\varrho^{\mathrm{C}}_{1,11}=\rho_{3,9},\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{2,8}=\rho_{4,6},\varrho^{\mathrm{C}}_{4,7}=\rho_{4,7},\varrho^{\mathrm{C}}_{2,12}=\rho_{4,10},\varrho^{\mathrm{C}}_{4,11}=\rho_{4,11},\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{2,15}=\rho_{4,13},\varrho^{\mathrm{C}}_{5,9}=\rho_{5,9},\varrho^{\mathrm{C}}_{5,8}=\rho_{6,7},\varrho^{\mathrm{C}}_{6,10}=\rho_{6,10},\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{8,9}=\rho_{6,11},\varrho^{\mathrm{C}}_{6,13}=\rho_{6,13},\varrho^{\mathrm{C}}_{5,12}=\rho_{7,10},\varrho^{\mathrm{C}}_{7,11}=\rho_{7,11},\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{5,15}=\rho_{7,13},\varrho^{\mathrm{C}}_{8,12}=\rho_{8,12},\varrho^{\mathrm{C}}_{6,16}=\rho_{8,14},\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{8,15}=\rho_{8,15},\varrho^{\mathrm{C}}_{9,12}=\rho_{10,11},\varrho^{\mathrm{C}}_{10,13}=\rho_{10,13},\varrho^{\mathrm{C}}_{9,15}=\rho_{11,13},\\\\[2.84544pt]
\varrho^{\mathrm{C}}_{10,16}=\rho_{12,14},\varrho^{\mathrm{C}}_{12,15}=\rho_{12,15},\varrho^{\mathrm{C}}_{13,16}=\rho_{14,15}\end{array}$

The nonzero coefficients of the density matrix
$\varrho^{\mathrm{D}}=\rho^{T_{\mathrm{D}}}_{\mathrm{ABCD}}$ are

$\begin{array}[]{lcl}\varrho^{\mathrm{D}}_{i,i}=\rho_{i,i}\quad i=1\cdots
16,\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{1,4}=\rho_{2,3},\varrho^{\mathrm{D}}_{1,6}=\rho_{2,5},\varrho^{\mathrm{D}}_{1,10}=\rho_{2,9},\varrho^{\mathrm{D}}_{3,5}=\rho_{3,5},\varrho^{\mathrm{D}}_{3,9}=\rho_{3,9},\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{2,8}=\rho_{4,6},\varrho^{\mathrm{D}}_{3,8}=\rho_{4,7},\varrho^{\mathrm{D}}_{4,10}=\rho_{4,10},\varrho^{\mathrm{D}}_{3,12}=\rho_{4,11},\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{3,14}=\rho_{4,13},\varrho^{\mathrm{D}}_{5,9}=\rho_{5,9},\varrho^{\mathrm{D}}_{5,8}=\rho_{6,7},\varrho^{\mathrm{D}}_{6,10}=\rho_{6,10},\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{5,12}=\rho_{6,11},\varrho^{\mathrm{D}}_{5,14}=\rho_{6,13},\varrho^{\mathrm{D}}_{8,9}=\rho_{7,10},\varrho^{\mathrm{D}}_{7,11}=\rho_{7,11},\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{7,13}=\rho_{7,13},\varrho^{\mathrm{D}}_{8,12}=\rho_{8,12},\varrho^{\mathrm{D}}_{8,14}=\rho_{8,14},\varrho^{\mathrm{D}}_{7,16}=\rho_{8,15},\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{9,12}=\rho_{10,11},\varrho^{\mathrm{D}}_{9,14}=\rho_{10,13},\varrho^{\mathrm{D}}_{11,13}=\rho_{11,13},\varrho^{\mathrm{D}}_{12,14}=\rho_{12,14},\\\\[2.84544pt]
\varrho^{\mathrm{D}}_{11,16}=\rho_{12,15},\varrho^{\mathrm{D}}_{13,16}=\rho_{14,15},\\\\[11.38092pt]
\end{array}$

The explicit nonzero elements of the reduced density matrix
$\rho_{\mathrm{AB}}$ (9) are

$\begin{array}[]{lcl}\varrho_{1,1}=\rho_{1,1}+\rho_{2,2}+\rho_{3,3}+\rho_{4,4},\\\\[2.84544pt]
\varrho_{2,2}=\rho_{5,5}+\rho_{6,6}+\rho_{7,7}+\rho_{8,8},\\\\[2.84544pt]
\varrho_{3,3}=\rho_{9,9}+\rho_{10,10}+\rho_{11,11}+\rho_{12,12},\\\\[2.84544pt]
\varrho_{4,4}=\rho_{13,13}+\rho_{14,14}+\rho_{15,15}+\rho_{16,16},\\\\[2.84544pt]
\varrho_{2,3}=\varrho_{2,3}^{*}=\rho_{5,9}+\rho_{6,10}+\rho_{7,11}+\rho_{8,12}.\end{array}$

The nonzero coefficients of the reduced density matrix
$\Omega^{\mathrm{A}}=\rho^{T_{\mathrm{A}}}_{\mathrm{AB}}$ are

$\begin{array}[]{lcl}\Omega^{\mathrm{A}}_{1,1}=\rho_{1,1}+\rho_{2,2}+\rho_{3,3}+\rho_{4,4},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{2,2}=\rho_{5,5}+\rho_{6,6}+\rho_{7,7}+\rho_{8,8},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{3,3}=\rho_{9,9}+\rho_{10,10}+\rho_{11,11}+\rho_{12,12},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{4,4}=\rho_{13,13}+\rho_{14,14}+\rho_{15,15}+\rho_{16,16},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{1,4}=\Omega^{\mathrm{A}}_{4,1}=\rho_{5,9}+\rho_{6,10}+\rho_{7,11}+\rho_{8,12}.\end{array}$

The nonzero coefficients of the reduced density matrix
$\Omega^{\mathrm{A}}=\rho^{T_{\mathrm{A}}}_{\mathrm{AC}}$ are

$\begin{array}[]{lcl}\Omega^{\mathrm{A}}_{1,1}=\rho_{1,1}+\rho_{2,2}+\rho_{5,5}+\rho_{6,6},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{2,2}=\rho_{3,3}+\rho_{4,4}+\rho_{7,7}+\rho_{8,8},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{3,3}=\rho_{9,9}+\rho_{10,10}+\rho_{13,13}+\rho_{14,14},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{4,4}=\rho_{11,11}+\rho_{12,12}+\rho_{15,15}+\rho_{16,16},\\\\[2.84544pt]
\Omega^{\mathrm{A}}_{1,4}=\Omega^{\mathrm{A}}_{4,1}=\rho_{3,9}+\rho_{4,10}+\rho_{7,13}+\rho_{8,14}.\end{array}$

The nonzero coefficients of the reduced density matrix
$\Omega^{\mathrm{B}}=\rho^{T_{\mathrm{B}}}_{\mathrm{BC}}$ are

$\begin{array}[]{lcl}\Omega^{\mathrm{B}}_{1,1}=\rho_{1,1}+\rho_{2,2}+\rho_{9,9}+\rho_{10,10},\\\\[2.84544pt]
\Omega^{\mathrm{B}}_{2,2}=\rho_{3,3}+\rho_{4,4}+\rho_{11,11}+\rho_{12,12},\\\\[2.84544pt]
\Omega^{\mathrm{B}}_{3,3}=\rho_{5,5}+\rho_{6,6}+\rho_{13,13}+\rho_{14,14},\\\\[2.84544pt]
\Omega^{\mathrm{B}}_{4,4}=\rho_{7,7}+\rho_{8,8}+\rho_{15,15}+\rho_{16,16},\\\\[2.84544pt]
\Omega^{\mathrm{B}}_{1,4}=\Omega^{\mathrm{B}}_{4,1}=\rho_{3,5}+\rho_{4,6}+\rho_{11,13}+\rho_{12,14}.\\\\[5.69046pt]
\end{array}$

##  Appendix B QFI

We here verify the QFI for the pair spins A and B. The alteration of the
thermal quantum state ${\rho}_{\mathrm{AB}}(T)$ (9) under a unitary
transformation can be described as
${\rho}(\Theta)=\exp[-i\mathcal{A}\Theta]{\rho}_{\mathrm{AB}}\exp[i\mathcal{A}\Theta]$,
where $\Theta$ is the phase shift and $\mathcal{A}$ is an operator. The
estimation accuracy for $\Theta$ is limited by the quantum Cramér-Rao
inequality
$\Delta\hat{\Theta}\geq\frac{1}{\sqrt{\nu\mathcal{F}({\rho}_{\Theta})}}$,
where $\hat{\Theta}$ is the unbiased estimator for $\Theta$, and $\nu$ is the
number of measurement repetitions. Therefore, the QFI of the given state can
be expressed as [68, 69]

$\begin{array}[]{cl}\mathcal{F}({\rho}_{\mathrm{AB}},\mathcal{A})=2\sum\limits_{i,j=1}^{4}\frac{(\tau_{i}-\tau_{j})^{2}}{\tau_{i}+\tau_{j}}\mid\langle\chi_{i}\mid\mathcal{A}\mid\chi_{j}\rangle\mid^{2},\end{array}$
(17)

where $|\chi_{i}\rangle$ and $\tau_{i}$ are, respectively, the eigenvectors
and their corresponding eigenvalues of the density matrix
${\rho}_{\mathrm{AB}}$. The complete sets of local orthonormal observables
$\\{\mathcal{A}_{\eta}\\}$ and $\\{\mathcal{B}_{\eta}\\}$ corresponding to the
subsystems, respectively, $\mathrm{A}$ and $\mathrm{B}$ might be employed to
quantify QFI of their states ${\rho}_{\mathrm{AB}}$. The QFI for a general
bipartite spin-1/2 system reads [70]

$\begin{array}[]{cl}\mathcal{F}=\sum\limits_{\eta}\mathcal{F}(\tilde{\rho},\mathcal{A}_{\eta}\otimes
I+I\otimes\mathcal{B}_{\eta}).\end{array}$ (18)

The local observables $\\{\mathcal{A}_{\eta}\\}$ and
$\\{\mathcal{B}_{\eta}\\}$ are orthonormal and can be written as

$\begin{array}[]{cl}\\{\mathcal{A}_{\eta}\\}=\\{\mathcal{B}_{\eta}\\}=\sqrt{2}\\{I,\;S^{x},\;S^{y},\;S^{z}\\},\end{array}$
(19)

which $I$ is $2\times 2$ identity matrix. Ultimately, by substituting Eq. (19)
to Eq. (17), QFI $\mathcal{F}$ is achieved.

##

## References

  * [1] O. Kahn, Molecular Magnetism; Wiley-VCH, (1993).
  * [2] C. Benelli, D. Gatteschi, Introduction to Molecular Magnetism: From Transition Metals to Lanthanides; Wiley-VCH Verlag GmbH & Co. KGaA: Weinheim, Germany, (2015).
  * [3] A. G. Ariño, F. Luis, S. Hill, and E. Coronado, Nat. Chem. 11, 301 (2019).
  * [4] S. V. Yurtaeva, I. F. Gilmutdinov and A. A. Rodionov, ACS Omega 4, 16450 (2019).
  * [5] B. Kintzel, M. Böhme, J. Liu et al., Chem. Commun. 54, 12934 (2018).
  * [6] J. Camarero and E. Coronado, J. Mater. Chem. 19, 1678 (2009).
  * [7] S.Y. Bodnar, L. Šmejkal, I. Turek et al., Nat. Commun. 9, 348 (2018).
  * [8] M. Sutradhar, E. C. B. A. Alegria, T. R. Barman, M. F. C. G. da Silva, C. -M. Liu and A. J. L. Pombeiro, Front. Chem. 8, 157 (2020).
  * [9] S. S. Tandon, S. D. Bunge , N. Patel, E. C. Wang and L. K. Thompson, Molecules 25, 5549 (2020)
  * [10] A. Burkhardt, E.T. Spielberg, H. Görls, W. Plass, Inorg. Chem. 47, 2485 (2008).
  * [11] L. Borsten, D. Dahanayake, M. J. Duff, A. Marrani and W. Rubens, Phys. Rev. Lett. 105, 100507 (2010).
  * [12] D. Maity, A. D. Jana, M. Debnath et al., Eur. J. Inorg. Chem. 3484 (2010).
  * [13] S. C. dra Manna, S. Manna, A. Paul et al., ChemSelec. 2, 3317 (2017).
  * [14] J. Tercero, E. Ruiz, S. Alvarez, A. R. Fortea, P. Alemany, J. Mater. Chem. 16, 2729 (2006).
  * [15] L. Salmon, P. Thuéry, E. Riviére, J. P. Costes, A. J . Motad and M. Ephritikhine,  New J. Chem. 38, 1306 (2014).
  * [16] S. Giri, S. Biswas, M. G.B. Drew, A. Ghosh and S. K. Saha, Inorgan. Chim. Act. 368, 152 (2011); Supporting Information of S. Giri, et al., Inorgan. Chim. Act. 368, 152 (2011).
  * [17] S. Wang, S.J. Trepanier, J.C. Zheng, Z. Pang, M.J. Wagner, Inorg. Chem. 31 2118 (1992).
  * [18] J. Fielden, J. Sprott, D. Long, P. Kögerler, L. Cronin, Inorg. Chem. 45, 2886 (2006).
  * [19] A. Burkhardt, E.T. Spielberg, H. Görls, W. Plass, Inorg. Chem. 47 2485 (2008).
  * [20] V. Coffman, J. Kundu, W.K. Wootters, Phys. Rev. A 61, 052306 (2000).
  * [21] M.B. Plenio, Phys. Rev. Lett. 95, 090503 (2005).
  * [22] G. Vidal, R.F. Werner, Phys. Rev. A 65, 032314 (2002).
  * [23] A. Peres, Phys. Rev. Lett. 77, 1413 (1996).
  * [24] M. Horodecki, P. Horodecki, and R. Horodecki, Phys. Lett. A  223, 1 (1996).
  * [25] J.-L. Guo, X.-L. Huang, and H.-S. Song, Phys. Scr. 76, 327 (2007).
  * [26] G.-H. Yang and L. Zhou, Phys. Scr. 78, 025703 (2008).
  * [27] S. Xu, X. Song, and L. Ye,  Quantum Inf. Process. 13, 1013 (2014).
  * [28] F. Wang, L.-P. Fu, and K.-T. Guo, Commun. Theor. Phys. 50, 341 (2008).
  * [29] S. D. Han, T. Tüfekçi, T. P. Spiller, and E. Aydiner,  Int. J. Theor. Phys. 56, 1474 (2017).
  * [30] W. K. Wootters, Phys. Rev. Lett. 80, 2245 (1998).
  * [31] V. Vedral, Rev. Mod. Phys. 74, 197 (2002).
  * [32] R. Horodecki, P. Horodecki, M. Horodecki and K. Horodecki, Rev. Mod. Phys. 81, 865 (2009).
  * [33] M.A. Nielsen and I.L. Chuang; Quantum Computation and Quantum Information; Cambridge University Press, Cambridge, (2000).
  * [34] C. H. Bennett, G. Brassard, C. Crepeau, R. Jozsa, A. Peres, and W. K. Wootters, Phys. Rev. Lett. 70, 1895 (1993).
  * [35] H. L. Huang, Int. J. Theor. Phys. 50, 70 (2011).
  * [36] H. Arian Zad and M. Rojas, Phys. E 126, 114455 (2021).
  * [37] M. Freitas, C. Filgueiras, and M. Rojas, Ann. Phys. 531, 1900261 (2019).
  * [38] G. C. Fouokeng, E. Tedong, A. G. Tene, M. Tchoffo and L. C. Fai, Phys. Lett. A 384, 126719 (2020).
  * [39] L. Bakmou, A. Slaoui, M. Daoud and R. A. Laamara,  Quant. Info.Process. 18, 163 (2019).
  * [40] A. Slaoui, L. Bakmou, M. Daoud and R. Ahl Laamara, Phys. Lett. A 383, 2241 (2019).
  * [41] R. Jafari and A. Akbari, Phys. Rev. A 101, 062105 (2020).
  * [42] U. Marzolino, T. Prosen, Phys. Rev. B 96, 104402 (2017).
  * [43] J. Liu, H. Yuan, X.-M. Lu and X. Wang,  J. Phys. A: Math. Theor. 53, 023001 (2020).
  * [44] B. L. Ye, L. Y. Xue, Y. L. Fang, S. Liu, Q. C. Wu, Y. H. Zhou and C. P. Yang,  Phys. E 115, 113690 (2020).
  * [45] B.-L. Ye, B. Li, Z.-X. Wang, X. Li-Jost and S.-M. Fei,  Sci. China Phys. Mech. Astron. 61, 110312 (2018).
  * [46] Y. Honda, Y. Kuramoto, and T. Watanabe, Phys. Rev. B 47, 11329 (1993).
  * [47] S. Brehmer, H.-J. Mikeska, M. Müller, N. Nagaosa, and S. Uchida, Phys. Rev. B 60, 329 (1999).
  * [48] T. Senthil and M. P. A. Fisher, Phys. Rev. Lett. 86, 292 (2001).
  * [49] A. Reischl and E. Müller-Hartmann, Eur. Phys. J. B 28, 173 (2002).
  * [50] L. Balents, M. P. A. Fisher, and S. M. Girvin, Phys. Rev. B 65, 224412 (2002).
  * [51] S. Capponi, P. Lecheminant and M. Moliner, Phys. Rev. B 88, 075132 (2013).
  * [52] A. Metavitsiadis and S. Eggert, Phys. Rev. B 95, 144415 (2017).
  * [53] H. Arian Zad and N. Ananikian, J. Phys.: Condens. Matter 29, 455402 (2017).;
  * [54] V. Štubňa and M. Jaščur, Acta Phys. Slov. 69, 75 (2019).
  * [55] H. Arian Zad and N. Ananikian, Eur. Phys. J. Plus 136, 597 (2021).
  * [56] C. J. Calzado, C. de Graaf, E. Bordas, R. Caballol and J. -P. Malrieu, Phys. Rev. B 67, 132409 (2003).
  * [57] A. Bühler, U. Löw, K. P. Schmidt and G. S. Uhrig, Phys. Rev. B 67, 134428 (2003).
  * [58] K. Karľová, J. Strečka, Act. Phys. Pol. A, 137, 595 (2020).
  * [59] K. Karľová, J. Strečka, J. Haniš and M. Hagiwara, Magnetochem. 6, 59 (2020).
  * [60] Y. Li, C. Liu, Q. Wang, H. Zhang and L. Hu, Optik 127, 9788 (2016).
  * [61] Q. Dong, A A. S. Manilla, I. L. Yáñez, G. H. Sun and S. H. Dong, Phys. Scr. 94, 105101 ( 2019).
  * [62] A. J. T. Arenas, Q. Dong, G. H. Sun, W.C. Qian and S. H. Dong, Phys. Lett. B 789, 93 (2019).
  * [63] Q. Dong, A. J. T. Arenas, G. H. Sun and S. H. Dong, Front. Phys. 15, 11602 (2020).
  * [64] R. A. Horn and C. R. Johnson, Matrix Analysis, (Cambridge University Press, New York) 415, pp. 205, (1985).
  * [65] C. Sabín and G. García-Alcaine, Eur. Phys. J. D 48, 435 (2008).
  * [66] G. Bowen, S. Bose, Phys. Rev. Lett. 87, 267901 (2001).; Y. Yeo, Phys. Rev. A 66, 062312 (2002).; Y. Zhou, G.F. Zhang, Eur. Phys. J. D 47, 227 (2008).
  * [67] R. Jozsa, J. Mod. Opt. 41, 12 (1994).
  * [68] A. S. Holevo, Probabilistic and Statistical Aspects of Quantum Theory, Math. Sci. Eng. (Edizioni della Normale, New York, NY, 1982).
  * [69] J. Liu, X. Jing, and X. Wang,  Phys. Rev. A 88, 042316 (2013).
  * [70] N. Li and S. Luo,  Phys. Rev. A 88, 014301 (2013).