# Nonlinearity in spin dynamics of frustrated Kagomé lattice system under
harmonic perturbation

Saumen<EMAIL_ADDRESS>Arindam
<EMAIL_ADDRESS>Reeta<EMAIL_ADDRESS>and
Nimisha<EMAIL_ADDRESS>Department of Physics, Dibrugarh
University, Dibrugarh 786 004, Assam, India

###### Abstract

In this study, we investigate the spin dynamics of a frustrated Kagomé lattice
system, focusing on the nonlinearity of spin oscillations induced by a
harmonic magnetic field under varying strengths of Dzyaloshinskii-Moriya
interaction (DMI), exchange field, and anisotropy energy. We have utilized
Poincaré Surface Sections (PSS) and Power Spectra (PS) for different DMI and
anisotropy energy to study the spin dynamics. Our findings reveal that when
the DMI strength, external field, anisotropy, and applied magnetic field are
weak, the oscillations are quasi-periodic, mostly dominated by the exchange
field. With the increase in the DMI strength, the oscillation of the system
becomes highly aperiodic. Strong anisotropy tends to induce periodic
oscillations, but increasing DMI eventually leads to chaotic behaviour.
Additionally, the external magnetic field destabilizes the periodicity of
oscillations in systems with weak easy-axis anisotropy, but the systems with
strong anisotropy, the oscillations remain unaffected by the external field’s
strength. Our analysis of magnon dispersion and magnetic resonance (MR)
spectra reveals multiple resonance peaks at higher DMI strengths, indicating a
complex interplay between spin wave excitation and system parameters. These
results underscore the importance of understanding the inherent DMI and
anisotropy in the Kagomé lattice during fabrication for various applications.
Moreover, our comprehensive analysis of spin dynamics in a Kagomé lattice
system demonstrates a clear transition from quasi-periodic to chaotic
oscillations with the increase in DMI strength.

###### pacs:

67.30.hj, 85.75.-d, 74.90.+n

##  I Introduction

Frustrated magnetic systems have garnered significant attention in the
scientific community due to their intriguing and often unexpected properties
that challenge conventional paradigms of magnetic ordering. They exhibit
unconventional ground states and display exotic magnetic excitation and quasi-
particles ghimire ; mielke . Thus, it may behave unconventionally and
contribute to unique magnetic and thermal properties. In general, frustration
in a magnetic system arises when the interactions within a given system cannot
be simultaneously minimized, leading to a considerable multitude of low-energy
configurations characterized by an approximate degeneracy that scales
exponentially with the size of the system khatami ; saha ; farhan ; hog .
Moreover, these states can also include spin liquids han ; yan ; nayak , where
magnetic moments remain disordered even at very low temperatures. Among the
various frustrated geometries, the Kagomé lattice holds particular interest
because of its significantly smaller energy gaps than the other geometries,
with an exponentially large number of singlet states below the lowest triplet
state han ; mendels ; li ; zorko ; du ; yin . Conversely, experiments
conducted on materials featuring nearly decoupled Kagomé layers, exemplified
by Herbertsmithite, reveal the absence of magnetic order even at extremely low
temperatures and a lack of discernible evidence for a spin gap. Nevertheless,
when the study of these materials is delving into these materials, the
magnetic anisotropy and the Dzyaloshinskii - Moria interaction (DMI) dmi1 ;
dmi2  plays a significant role in dictating the nature of the ground state.

The DMI is an antisymmetric exchange interaction emerging in systems having no
inversion symmetry and promotes chiral magnetic phases such as skyrmions hog ;
liu1 ; chen ; petrova ; lee ; drozdov  in two dimensional (2D) lattices. These
unique magnetic textures have attracted significant attention due to their
potential applications in spintronics and information storage wiesendanger ;
zhang ; marrows ; dupe ; zhou ; fert ; ishida ; zhang2 ; tokura . DMI can
facilitate the controlled movement and manipulation of domain walls je ;
ajejas ; brandao ; wang ; dejong . Moreover, the nontrivial spin textures
induced by DMI leading to a transverse deflection of electrons as they move
through the material. Thus, it hosts the topological Hall effect and also
leads to a modification in the spin wave spectrum that can be probed
experimentally ren ; gao ; huang ; wang1 . So, the role of DMI in spin
dynamics and excitation spectrum of 2D lattice system is too inherent. The
Kagomé lattice, with its inherent frustration from competing interactions,
tends to resist conventional magnetic ordering. The DMI introduces non-
collinear spin arrangements, favouring the formation of chiral magnetic order.
Thus understanding the spin dynamics in the Kagomé lattice structure and its
interplay with DMI offers a rich playground for exploring novel magnetic
phenomena and opens avenues for potential applications in spintronics and
magnetic technology.

The study of spin dynamics is of paramount importance, as it forms the
cornerstone for understanding the behaviour of magnetic materials, providing
essential insights into their fundamental properties and enabling the design
of advanced technologies acharjee1 ; acharjee2 ; acharjee3 ; acharjee4 ;
acharjee5 ; liu . Over the years, the time evolution of magnetic moments
unveiled the intricate interplay of magnetic interactions by exploring diverse
magnetic phases in different magnetic systems. This understanding is crucial
for developing magnetic materials with tailored functionalities, including
applications in information storage, spintronics, and magnetic sensors
wiesendanger ; zhang ; marrows ; dupe ; zhou ; fert ; ishida ; zhang2 ; tokura
; gottscholl ; khan ; yan1 ; wong . The deformed Kagomé lattice, characterized
by its geometric frustration arising from competing interactions within an
array of corner-sharing triangles, serves as a fascinating platform for
investigating unconventional magnetic behaviour mendels ; li ; zorko ; du ;
yin . It is due to its ability to foster exotic magnetic ground states such as
spin liquids, thereby introducing a complex tapestry of interactions among
neighbouring spins. Moreover, the Kagomé system with anisotropy du ; yin ,
becomes even more nuanced, giving rise to intriguing phenomena that challenge
the conventional understanding of magnetic systems. This paper delves into the
intriguing dynamics of the frustrated spins within an anisotropic Kagomé
lattice system in the presence of DMI, induced by harmonic perturbation. The
basic objective of this work is to understand the role of magnetic anisotropy,
DMI in spin dynamics of Kagomé system under harmonic perturbation.

The paper is organized as follows. Section II outlines the minimal theoretical
model of the proposed frustrated Kagomé system and derives the explicit form
of the Hamiltonian. In Section III, we delve into the study of the Poincaré
Surface Section (PSS), power spectra, and magnetic resonance considering the
effect of anisotropy and DMI under harmonic perturbation. This analysis aims
to comprehend the nature of the oscillation and stability of spins within the
proposed system. Additionally, we explore the finite time Largest Lyapunov
exponent (LLE) to investigate the chaos in this system. Section V concludes
our work with a concise summary.

##  II Theoretical framework

###  II.1 Spin Dynamics

The schematic representation of two dimensional anisotropic Kagomé lattice
with exchange interactions, $J_{1}$, $J_{2}$, $J_{3}$ and $J_{4}$ is shown in
Fig. 1. We consider the lattice sites A, B and C for our reference with the
onsite magnetic moment $\mathbf{S}_{\text{A}}$, $\mathbf{S}_{\text{B}}$ and
$\mathbf{S}_{\text{C}}$ respectively. The Hamiltonian of the deformed Kagomé
lattice under harmonic perturbation can be expressed as

$\mathcal{H}=J_{\alpha}\sum_{\left<i,j\right>}\mathbf{S}_{i}.\mathbf{S}_{j}+K\sum_{i}(\mathbf{S}_{i}.\hat{e}_{z})^{2}+\sum_{ij}\mathbf{D}_{ij}^{\text{eff}}.\\{\mathbf{S}_{i}\times\mathbf{S}_{j}\\}\\\
-\mathbf{\mathcal{B}}_{\text{ext}}.\mathbf{S}_{i}$ (1)

Figure 1: Schematic representation of two-dimensional anisotropic Kagomé
lattice with the exchange interactions $J_{1}$, $J_{2}$, $J_{3}$ and $J_{4}$
between the spins at the sub-lattice ABCDE. The proposed unit cell with
lattice sites ABCDE is denoted by the dotted rectangle. The spin moments at
the lattice site A, B and C are denoted by $\mathbf{S}_{\text{A}}$,
$\mathbf{S}_{\text{B}}$ and $\mathbf{S}_{\text{C}}$ respectively while the
moments at the D and E are chosen to be (0,0,1) for simplicity.

where, the indices i and j denote the lattice sites and the sum
$\left<..\right>$ is over the nearest neighbours. The first term of Eq. (1)
represents the exchange interaction, and $J_{\alpha}$ are the exchange
coefficients between the spins $\mathbf{S}_{i}$ and $\mathbf{S}_{j}$ at the
atomic sites i and j. Here, we consider that the sublattice ABC for reference
with $\alpha=(1-4)$ corresponds to the nearest-neighbor coupling of the
lattice ABC as shown in Fig. 1. The parameters $J_{\alpha}>0$ for
antiferromagnetic Kagomé system; $K$ is the onsite anisotropy with easy axis
along $\hat{e}_{z}$ direction. The third term of Eq. (1) represents the
contribution due to the total DMI energy. The parameter
$\mathbf{D}_{ij}^{\text{eff}}$ incorporates the effective DMI vector
calculated by performing the sum over all nearest neighbours between the spins
obtained via the three-site model defined as

$\mathbf{D}_{ij}^{\text{eff}}=\sum_{k}\mathbf{\mathcal{D}}_{ijk}(\mathbf{R}_{ki},\mathbf{R}_{kj},\mathbf{R}_{ij})$
(2)

where, $\mathbf{R}_{ki}$, $\mathbf{R}_{kj}$ are the distance vector from the
site $k$ to the atomic sites $i$ and $j$ respectively. The last term of Eq.
(1) corresponds to the harmonic magnetic field defined as

$\mathbf{\mathcal{B}}_{\text{ext}}=B_{0}\sin(\Omega t)\hat{e}_{z}$ (3)

where, $B_{0}$ and $\Omega$ are the strength and frequency of the external
magnetic field applied along $z$-direction. The time evolution of the spins in
the i${}^{\text{th}}$ sublattice in Kagomé system can be studied by using the
Landau - Lifshitz - Gilbert (LLG) equation.

$\frac{d\mathbf{S}_{i}}{dt}=-\gamma(\mathbf{S}_{i}\times\mathbf{H}^{\text{eff}}_{i})+\alpha\mathbf{S}_{i}\times(\mathbf{S}_{i}\times\mathbf{H}^{\text{eff}}_{i})$
(4)

where $\gamma$ is the gyromagnetic ratio and $\alpha$ is the Gilbert damping
coefficient. The effective fields $\mathbf{H}^{\text{eff}}_{i}$ for the
respective sub lattices can be obtained by
$\mathbf{H}^{\text{eff}}_{i}=-\frac{\delta\mathcal{H}}{\delta\mathbf{S}_{i}}$.
Using Eq. (1) and (2), we obtain nine nonlinear coupled first-order
differential equations for the spins in the sublattice ABC. An explicit form
of the coupled time evolution equation is given in Appendix (A).

Figure 2: The PSS for the Kagomé system having
$(J_{1},J_{2},J_{3},J_{4})=(0.316,0.285,0.252,0.145)$ for different values of
$D$ and $K$ under magnetic field strength $B_{0}=0.1$ and frequency
$\Omega=1$. The plots in the top panel, middle and bottom panel are
respectively correspond to Spin at the lattice site A (green), B (blue) and C
(red) respectively.

Figure 3: The PS for Kagomé system for different values of $D$ and $K$ under
magnetic field strength $B_{0}=0.1$ and frequency $\Omega=1$.

Figure 4: The PSS for Kagomé system having
$(J_{1},J_{2},J_{3},J_{4})=(0.316,0.285,0.252,0.145)$ for different values of
$D$ and $K$ under magnetic field strength $B_{0}=1.0$ and frequency
$\Omega=1$. The plots in the top panel, middle and bottom panel are
respectively correspond to Spin at the lattice site A (green), B (blue) and C
(red) respectively.

Figure 5: The PS for Kagomé system for different values of $D$ and $K$ under
magnetic field strength $B_{0}=1.0$ and frequency $\Omega=1$.

##  III Simulations and Numerical Results

This section is divided into two parts. The first part describes the dynamical
behaviour of the Kagomé system using Poincaré Surface Section (PSS) and Power
Sectra (PS). The subsequent part showcases the nature of oscillations using
Magnetic resonance (MR) and Largest Lyapunov Exponents (LLE) and deliberates
on our numerical findings. For all our simulation we set the initial values
$S_{\text{A}}=(1,0,0)$, $S_{\text{B}}=(0,0,1)$ and $S_{\text{C}}=(1,0,0)$. The
simulations are performed by considering the exchange interaction values for
the deformed Kagomé system Rb2Cu3SnF12 with
$(J_{1},J_{2},J_{3},J_{4})=(0.316,0.285,0.252,0.145)$ khatami . Additionally,
we explore the impact of exchange coupling on spin dynamics in a subsequent
section of this work.

###  III.1 Poincaré Surface Section (PSS) and Power Spectra (PS)

Fig. 2 illustrates the PSS of the spins located at the lattice sites A, B, and
C for weak magnetic field strength, $B_{0}=0.1$ with $\Omega=1$ under weak and
strong anisotropy conditions with different DMI strength. It is found that the
system executes highly aperiodic and nonlinear oscillations for
$(K,D)=(0.1,0.1)$ as seen from Figs. 2(a) – 2(c). It is to be noted that the
mutual exchange coupling dominates in this condition suppressing the impact of
DMI and magnetic field. The corresponding PS in Fig. 3(a) indicates the
presence of multiple low-frequency modes present in the oscillations. The
presence of Kolmogorov-Arnold-Moser (KAM) islands indicates the non-chaotic
nature of the oscillations. With the increase in the strength of DMI to 0.5,
the number of frequency modes increases and thereby enhances the aperiodicity
of the oscillations as illustrated in the PSS Figs. 2(e) – 2(g) and
corresponding PS in Fig. 3(b). For a system with strong DMI and weak
anisotropy, i.e., $(K,D)=(0.1,1.0)$, the PSS is found to display random blotch
of points disappearing the KAM islands completely indicating the chaotic
nature of the oscillations as observed from Figs. 2(i) – 2(k). The
corresponding power spectra shown in Fig. 3(c), consist of a broad band of
frequency modes indicating the non-linear and chaotic nature of oscillations.
This is due to the reason that DMI dominates over the anisotropy and the
magnetic field. So, it is evident that the system undergoes a transition from
quasi-periodic to strongly aperiodic oscillation with the increase in DMI
strength for a system having low anisotropy and in the presence of a low
magnetic field.

Figs. 2(m) – 2(w), represent the PSS for oscillations of the coupled spin
oscillations in Kagomé system having strong uniaxial anisotropy $K=1$. The
coupled spin oscillations are found to be regular with the presence of
multiple KAM islands for a system with $(K,D)=(1.0,0.1)$ as seen from Figs.
2(m) – 2(o). The presence of only a few frequency modes in the corresponding
PS in Fig. 3(d) indicates the regular nature of the oscillation. It is due to
the stabilization of the spin oscillation due to the presence of strong
anisotropy energy in the system. With the increase in $D$ to $0.5$, the KAM
islands disappear with some splatter random points present in the PSS,
indicating that the motion undergoes a transition from periodic to aperiodic
oscillations as seen in Figs. 2(q) – 2(s). The corresponding PS in Fig. 3(e),
consists of multiple frequency modes and indicates similar characteristics of
the spin oscillations. A similar nature of oscillation in the PSS is also
observed in Figs. 2(u) – 2(w) for a system with strong anisotropy and strong
DMI strength, i.e., $(K,D)=(1.0,1.0)$. The corresponding PS consist of more
frequency bands as seen from Fig. 3(f). Thus, it is seen that the systems
having strong anisotropy execute nearly chaotic behaviour for moderate and
strong DMI strength. This is due to the reason that DMI introduces a torque
that disturbs the regularity of spin oscillations. So, for strong anisotropic
systems, the non-linearity is more evident as the non-linear oscillations are
more prominent along the anisotropic direction. For low anisotropic systems,
the nearly chaotic characteristics are seen only for the systems having strong
DMI strengths.

Figure 6: MR sprectra of Kagomé system for $B_{0}=1.0$ with DMI strength (a)
$D=0.1$, (b) $D=0.5$ and (c) $D=1.0$. The solid and dotted lines respectively
correspond to $K=0.1$ and $K=1.0$.

Figure 7: Density plot of MR with $D$ and frequency $\omega$ for $B_{0}=1$
with anisotropy energy $K=0.1$ (top panel) and $K=1.0$ (bottom panel).

Figure 8: Lyapunov exponents for Kagomé system for different choices of DMI
strength with $\Omega=1$ under (a) $B_{0}=0.1,K=0.1$, (b) $B_{0}=0.1,K=1.0$,
(c) $B_{0}=1.0,K=0.1$ and (d) $B_{0}=1.0,K=1.0$.

To understand the impact of a strong magnetic field, the Poincaré surface of
section (PSS) of the coupled spin system has been analyzed for a magnetic
field strength of $B_{0}=1.0$ in Fig. 4 under both weak and strong anisotropy
conditions, with varying Dzyaloshinskii-Moriya interaction (DMI) strengths.
The corresponding phase space (PS) is shown in Fig. 5. We observed that the
nature of the oscillations is quasi-periodic for $(K,D)=(0.1,0.1)$, as
evidenced by the PSS in Figs. 4(a) – 4(c) and the PS in Fig. 5(a). This
characteristic is similar to the system with a low magnetic field strength, as
previously discussed. The PSS exhibits random points for $D\geq 0.5$,
indicating aperiodic oscillations, as seen in Figs. 4(e) – 4(k) under low
anisotropy conditions. This behavior is attributed to the increased magnetic
field strength disrupting the spin oscillations, thus enhancing the
aperiodicity of the oscillations, as evident from the multiple frequency nodes
in Figs. 5(b) and 5(c). For materials with strong anisotropy, the system’s
characteristics in the PSS are similar to those of systems with strong
anisotropy and low magnetic field conditions, as shown in Figs. 2(m) – 2(w).
With an increase in anisotropy energy $(K=1)$, the oscillations of the coupled
spin system remain independent of the magnetic field strength, as demonstrated
in Figs. 4(m) – 4(w) and Figs. 2(m) – 2(w).

###  III.2 Magnon dispersion and Magnetic Resonance (MR) in Kagomé system

Over the years, Magnetic Resonance (MR) has been a standard tool for
researchers to understand various couplings like magnetic anisotropy, exchange
coupling, damping, and other key parameters that are essential for designing
and gaining a deeper understanding of the magnetic structures and devices
elkins ; bakuzis ; anisimov ; yuan3 ; Chen91 . Moreover, it is a standard tool
for probing spin wave excitation’s and spin dynamics bauer ; flovik ; bertelli
; crooker . Also, as already evident from Figs. 2 –5, the coupled spin
oscillations of our system are strongly dependent on DMI strength, anisotropy,
and applied magnetic field. Thus, to gain a better understanding of the nature
of oscillation and its correlation with the applied magnetic field it is
necessary to explore MR of our proposed system under weak and strong
anisotropy conditions for different DMI strengths.

To obtain the magnon dispersion relation we consider the sub-lattice ABC of
the Kagomé system defined by the order parameters $\mathbf{S}_{\text{A}}$,
$\mathbf{S}_{\text{B}}$ and $\mathbf{S}_{\text{C}}$. In view of Eq. (4), the
equation of motion of the sub-lattice ABC of the Kagomé system can be written
as

$\displaystyle\dot{\mathbf{S}}_{\text{A}}$
$\displaystyle=-\gamma(\mathbf{S}_{\text{A}}\times\mathbf{H}^{\text{eff}}_{\text{A}})+\alpha\mathbf{S}_{\text{A}}\times(\mathbf{S}_{\text{A}}\times\mathbf{H}^{\text{eff}}_{\text{A}})$
(5) $\displaystyle\dot{\mathbf{S}}_{\text{B}}$
$\displaystyle=-\gamma(\mathbf{S}_{\text{B}}\times\mathbf{H}^{\text{eff}}_{\text{B}})+\alpha\mathbf{S}_{\text{B}}\times(\mathbf{S}_{\text{B}}\times\mathbf{H}^{\text{eff}}_{\text{B}})$
(6) $\displaystyle\dot{\mathbf{S}}_{\text{C}}$
$\displaystyle=-\gamma(\mathbf{S}_{\text{C}}\times\mathbf{H}^{\text{eff}}_{\text{C}})+\alpha\mathbf{S}_{\text{C}}\times(\mathbf{S}_{\text{C}}\times\mathbf{H}^{\text{eff}}_{\text{C}})$
(7)

where, $\mathbf{H}^{\text{eff}}_{\text{A}}$,
$\mathbf{H}^{\text{eff}}_{\text{B}}$ and $\mathbf{H}^{\text{eff}}_{\text{C}}$
are the effective field of sub-lattice A, B and C respectively. In order to
obtain the magnon dispersion we consider a spin wave of frequency $\omega$ and
propagation vector $k$ defined as
$\mathbf{S}_{j}=\mathbf{S}_{j}^{0}+\delta\mathbf{S}_{j}e^{\\{i(\mathbf{k}.\mathbf{r}-\omega
t)\\}}$ with, $j=1,2,3$ corresponding to the sub-lattice A, B and C
respectively, with $\mathbf{S}_{j}^{0}$ being the ground state spin moment of
the sub-lattice $j$. Here, we consider a small deviation perpendicular to
$\mathbf{S}_{j}^{0}$ as $\delta\mathbf{S}_{j}$ is. The Eqs. (5) - (7) can be
linearized by using Kittel’s approach. In view of this the coefficient of
$\delta\mathbf{S}_{\text{A}}$, $\delta\mathbf{S}_{\text{B}}$ and
$\delta\mathbf{S}_{\text{C}}$ can be expressed in the form of a determinant as

$\begin{matrix}\mathcal{S}==\left(\begin{array}[]{cccccc}-i\omega&\mathcal{R}_{1}&D^{z}&-J_{1}&-D^{z}&-J_{4}\\\
-\mathcal{R}_{1}&-i\omega&J_{1}&D^{z}&J_{4}&-D^{z}\\\
-D^{z}&-J_{1}&-i\omega&\mathcal{R}_{2}&D^{z}&-J_{2}\\\
J_{1}&-D^{z}&-\mathcal{R}_{2}&-i\omega&J_{2}&D^{z}\\\
D^{z}&-J_{4}&-D^{z}&-J_{2}&-i\omega&\mathcal{R}_{3}\\\
J_{4}&D^{z}&J_{2}&-D^{z}&-\mathcal{R}_{3}&-i\omega\\\
\end{array}\right)\end{matrix}$ (8)

where, $\mathcal{R}_{1}=-i\alpha\omega-\gamma(J_{1}-J_{4}-2K)$,
$\mathcal{R}_{2}=-i\alpha\omega+\gamma(J_{1}+J_{2}+2K)$ and
$\mathcal{R}_{3}=i\alpha\omega+\gamma(J_{2}-J_{4}-2K)$. The magnon dispersion
relation can be obtained by solving the secular determinant
$det(\mathcal{S})=0$ for the frequency $\omega$, can be expressed as

$\mathcal{A}_{1}\omega+\mathcal{A}_{2}\omega^{2}+\mathcal{A}_{3}\omega^{3}+\mathcal{A}_{4}\omega^{4}+\mathcal{A}_{5}\omega^{5}+\mathcal{A}_{6}\omega^{6}+\xi=0$
(9)

where, the $\mathcal{A}_{i}$ with $i=1-6$ and $\xi$ are the coefficients whose
explicit form is displayed in the Appendix (B).

The numerically simulated results of Eq. (9) has been presented in Figs. 6 and
7 under weak and strong anisotropy conditions for different DMI strengths. We
observed that under low DMI conditions, the system displays a single resonance
peak as shown in Fig. 6(a) for low and high anisotropy conditions. This is due
to the reason that the oscillations of the coupled system are nearly periodic
in this condition. With the rise in DMI strength to $0.5$ and $1$, the system
displays multipeaked resonance spectra as seen in Figs. 6(b) and 6(c). This is
due to the reason that: (1) As the DMI strength increases, the oscillations of
the coupled spin system undergo a transition from regular to aperiodic
oscillations as observed earlier in Figs. 2 –5 and (2) DMI can cause a change
in the excitation of standing spin wave modes resulting in multipeaked MR
spectra. The presence of multiple frequencies in the PS in this regime
indicates the possibility of multiple resonance frequencies in the system. It
is to be noted that systems with strong magnetic anisotropy can significantly
influence the excitation of standing spin wave modes, resulting in the
increase of MR peaks in the resonance spectra as depicted by the dotted lines
in Figs. 6(b) and 6(c). Additionally, due to the aperiodic nature of the
oscillations, it is essential to explore different initial conditions to
understand the system’s sensitive dependence.

Figure 9: Variation of $S_{z}(t)$ with $t$ and $J_{1}$ considering
$J_{3}=0.25$, $J_{4}=0.15$, $B_{0}=0.1$, $D=0.5$, $K=0.5$ and $\Omega=1$. The
plots in the left panel (a)-(c) are for $J_{2}=0.25$, middle panel (d)-(f) are
for $J_{2}=0.5$ and (g)-(i) are for $J_{2}=0.25$.

Figs. 7 displays the density plot of the MR for different DMI and anisotropy
strength for $B_{0}=1.0$. In this case, we set the initial conditions as
$\\{\textbf{S}_{\text{A}},\textbf{S}_{\text{B}},\textbf{S}_{\text{C}}\\}=\\{(1,0,0),(0,1,0),(1,0,0)\\}$.
Under low anisotropy conditions, the oscillations of the first spin
($\textbf{S}_{\text{A}}$) of the system are predominantly in the low DMI
profile while it remains nearly independent with a maximum $\sim 3$ GHz for
moderate and high DMI strength. Moreover, the presence of multiple frequency
bands at low DMI conditions indicates the sensitive dependence on initial
conditions and the chaotic nature of oscillations. A similar characteristic is
observed for ($\textbf{S}_{\text{B}}$) and ($\textbf{S}_{\text{C}}$) in the
low DMI regime. However, for $\textbf{S}_{\text{B}}$ the MR spectra are
predominated by broad peaks around low, mid, and high frequency for moderate
and high DMI profiles. For $\textbf{S}_{\text{C}}$, MR peaks are observed at
moderate and high-frequency conditions for strong DMI strength as seen in Fig.
7(c). For systems having high anisotropy, multiple resonance frequencies are
observed as the oscillations of the coupled spin system are highly aperiodic
under low and moderate DMI strength while a quasi-periodic nature of resonance
spectra is found for strong DMI systems as seen in Figs. 7(d) - 7(f). Thus,
from Fig. 7, we observe that the oscillations of the coupled spin system not
only depend on the DMI strength but also have a significant dependence on the
anisotropic energy of the system. Also, by understanding the MR spectroscopy,
it may be possible to probe the complete magnetization dynamics of the
proposed system can be utilized to fabricate microwave devices.

###  III.3 Largest Lyapunov Exponents $(\Lambda)$

To comprehend the nature of spin oscillations of Kagomé spectra we analyze the
time evolution of the order parameters in Figs. 2 \- 5. However, to gain
insight into the stability of the system we have calculated the Lyapunov
Exponents (LE) sandri ; wolff  of the system. The LEs are the indicators that
measure the rate of convergence or divergence of the nearby trajectories. In
the context of chaotic systems, a substantial increase in the rate of
divergence is observed, leading to the emergence of positive values of LE.
Conversely, for regular systems, the LE exhibits a negative value.

The LE can be obtained by computing the average sequence of distances
$\Pi_{i}$ (where, $i=0,1,....N$) between the nearest neighbour trajectories
over a finite time. The underlying concept involves initiating two closely
situated points with initial displacements of $\mathcal{P}|_{(0,0)}$ and
$\mathcal{Q}|_{(0,0)}$ at time $t=0$. Thus, we can express
$\mathcal{Q}|_{(0,0)}=\mathcal{P}|_{(0,0)}+\bf{\Pi}_{0}$ initially at $t=0$,
where $\bf{\Pi}_{0}$ signify the initial displacement between the two points.
The points $\mathcal{P}|_{(0,0)}$ and $\mathcal{Q}|_{(0,0)}$ evolve after a
time interval $\tau$, resulting in new positions $\mathcal{P}|_{(0,\tau)}$ and
$\mathcal{Q}|_{(0,\tau)}$. So, we have
$\mathcal{Q}|_{(0,\tau)}=\mathcal{P}|_{(0,\tau)}+\Pi_{1}$, where
$\Pi_{1}=|\Pi_{1}|$ is the new distance between the two points. Considering it
to be the fresh reference point, we can write,
$\mathcal{Q}|_{(1,0)}=\mathcal{P}|_{(1,0)}+\frac{\Pi_{0}}{\Pi_{1}}\mathbf{\Pi_{1}}$,
with $\mathcal{P}|_{(1,0)}=\mathcal{P}|_{(0,\tau)}$. Following this iterative
approach one can generate the second set of points ($\mathcal{P}|_{(1,\tau)}$,
$\mathcal{Q}|_{(1,\tau)}$), where,
$\mathcal{Q}|_{(1,\tau)}=\mathcal{P}|_{(1,\tau)}+\mathbf{\Pi}_{2}$, Resetting
the reference points again to
$\mathcal{Q}|_{(2,0)}=\mathcal{P}|_{(2,0)}+\frac{\Pi_{0}}{\Pi_{2}}\mathbf{\Pi_{2}}$,
where, $\mathcal{P}|_{(2,0)}=\mathcal{P}|_{(1,\tau)}$ and continuing the
iteration further one can obtain new set of points for a single time step of
duration $\tau$. Repeating the steps for N times and generating a sequence of
distances $\Pi_{i}=|\bf{\Pi}_{i}|$, where i = 1, 2,… N we obtain the finite
time Largest Lyapunov Exponent (LLE), defined as acharjee5 ; liu

$\Lambda=\frac{1}{N\tau}\sum_{i=1}^{N}\log_{e}\left(\frac{\bf{\Pi}_{i}}{\bf{\Pi}_{0}}\right)$
(10)

Here, $\Lambda$ is independent of $\tau$ as it is considered $\tau\ll 0$. For
normal trajectories $\lim_{N\rightarrow\infty}\Lambda\leq 0$ while
$\lim_{N\rightarrow\infty}\Lambda>0$ for chaotic trajectories liu .

The LLE of the Kagomé system corresponding to the PSS in Fig. 2 and Fig. 4, is
displayed in Fig. 8. The inset pictures corroborate the precise DMI value
above which the system exhibits chaotic oscillations. The spin oscillations
are found to be mostly stable and non-chaotic for a system with weak
anisotropy energy and weak DMI strength i.e., $(K,D)=(0.1,0.1)$ in the
presence of a weak harmonic field of frequency $\Omega=1$. With the increase
in DMi strength to $0.5$ the system displays highly unstable characteristics.
For a system with strong DMI strength and weak anisotropy energy i.e.,
$(K,D)=(0.1,1.0)$, the LLE is found to be positive suggesting that spin
oscillations are strongly chaotic as seen from Fig. 8(a). The fig inset shows
that the LLE, ($\Lambda<0$) for $D\leq 0.5$ indicate the normal to chaotic
crossover region in this condition. A similar characteristic is also seen in
Fig. 8(b) for $K=1.0$ and $B_{0}=0.1$. In this case, the system also shows
chaotic behaviour near the region $D\geq 0.5$. Although the system’s
characteristics are identical, in this condition the system becomes strongly
chaotic for $D\sim 0.5$ and then execute unstable characteristics for $0.5\leq
D<0.7$. With the further increase in DMI strength, the system becomes strongly
chaotic as seen from the inset of Fig. 8(b). The spin oscillations are stable
for a system with weak anisotropy energy $K=0.1$ and weak DMI strength $D=0.1$
under strong magnetic field strength $B_{0}=1.0$. However, it becomes highly
chaotic for $D>0.5$ as seen from Fig. 8(c). It is to be noted that a stiff
rise in the LLE is observed at $D\sim 0.5$ indicating a sudden transit from
regular to chaotic regime in this condition. The nature of the oscillations is
found to be drastically different in Fig. 8(d) for a system with
$(K,D)=(1.0,1.0)$ with $B_{0}=1.0$. In this case, the LLE is positive for
$D\geq 0.3$ indicating the early transit from regular to chaotic region in
this regime.

###  III.4 Coupled Spin Oscillations for different exchange strength

While our primary analyses are confined to the exchange interaction of
Rb2Cu3SnF12, we are also interested in examining the impact of various
exchange interactions on spin dynamics. In Fig. 9, we present the study of the
z-component of spin for different Kagomé lattices, characterized by different
values of $J_{1}$ and $J_{2}$. For this analysis, we set, $J_{3}=0.25$,
$J_{4}=0.15$, $B_{0}=0.1$, $D=0.5$, $K=0.5$ and $\Omega=1$ for our analysis.
The oscillations of the system are found to be nearly periodic for low values
of $J_{1}$ and die out over time. A signature of aperiodic oscillations is
noticed with the rise in $J_{1}$. However, the effect of damping torque is
enhanced in this condition resulting in the spin oscillations at site A dying
out too rapidly in this scenario as seen from Fig. 9(a). For spins at the
lattice site B and spin C, the oscillations are found to have low periodicity
and are not much affected by the increase in $J_{1}$. Although some
aperiodicity in spin oscillations is noticed, it is found to have periodicity
and display weak damping with the increase in time $(t)$ as depicted in Fig.
9(b) and Fig. 9(c). The rapid damping of the oscillations at the lattice site
A is due to its coupling with the spin at the site B which tends to stabilize
the oscillations.

To explore the effect of exchange coupling between spins at sites B and C on
the coupled spin dynamics, we increased $J_{2}$ to $0.5$. The damping in the
spin oscillations is reduced and a nonlinear characteristic is observed for
low values of $J_{1}$ in this condition. However, significant damping is
observed for higher values of $J_{1}$. It is to be noted that the frequency of
oscillations also increases with the increase in the value of $J_{2}$ to
$0.5$. A similar characteristic is also observed for spins at the lattice
sites B and C, characterized by a significant damping noticed for higher
values of $J_{1}$. This is due to the reason that with the increase in
exchange interaction, the spin moments at sites A and B tend to align and
hence enhance the damping. However, as spin at site B is connected to A and C
via the exchange coupling $J_{1}$ and $J_{2}$, so for low values of $J_{1}$
and $J_{2}$ the damping is comparatively lower as compared to the higher
values of it as observed from Figs. 9(d) – 9(f). However, the frequency of
oscillations increases with the increase in the strength of $J_{1}$ for the
spins at sites A and B. Although the characteristics are found to be similar
in Figs. 9(g) – 9(i), the further increase of $J_{2}$ to $0.75$, a higher
damping is observed. In summary, oscillations in the Kagomé system exhibit
periodicity with low exchange interaction, moderate DMI strength, under
moderate or strong anisotropy energy. Nonlinear characteristics emerge for the
Kagomé systems with moderate exchange coupling, while strong exchange coupling
results in highly damped oscillations under moderate DMI strength and
anisotropy conditions.

##  IV Conclusions

In summary, we have investigated the spin dynamics of frustrated Kagomé
lattice system, characterizing the dependence of the nonlinearity of the spin
oscillations under different DMI strengths, exchange field, and anisotropy
energy induced via harmonic magnetic field. We showed that when the strength
of all the interacting parameters DMI, external field, exchange field,
anisotropy and applied magnetic field are low, the oscillations are quasi-
periodic as the exchange field dominates the spin oscillations. The increase
in DMI strength adds to the aperiodicity thereby making the oscillations
highly aperiodic. A strong anisotropy tends to make the oscillations periodic
but an increase in DMI eventually makes the oscillations chaotic. It was also
found that the strength of the external magnetic tends to destabilize the
periodicity of oscillations for a Kagomé system with weak easy-axis anisotropy
energy while the nature of the oscillations remains independent of the
strength of the external magnetic field for a Kagomé system with strong
anisotropy. In conclusion, our comprehensive analysis of spin dynamics in a
Kagomé lattice system demonstrates a clear transition from quasi-periodic to
chaotic oscillations as DMI strength increases. The study of magnon dispersion
and MR spectra reveals multiple resonance peaks at higher DMI strengths,
suggesting a complex interplay between spin wave excitation and system
parameters. The results of our analysis shed light on the importance of
understanding the inherent DMI and anisotropy of the Kagomé lattice structure
during the fabrication process for different applications. Though, Kagomé
lattice structure boasts of unique and exotic characteristics, the inherent
DMI and anisotropy of the material can greatly impact its overall behaviour
and performance. Another encouraging route for future study is to explore the
different means by which the DMI and the anisotropy of the system can be
controlled. As a concluding remark, this study provides valuable insights into
the intricate dynamics of spin systems, with potential implications for
designing advanced magnetic materials and devices. Such studies can add a new
dimension to the Kagomé lattice structure from the research perspective as
well as in potential technological applications.

##  Appendix A Time evolution equations for the Spin at the sub-lattice ABC

The time evolution equation for the spin at the lattice site A, B and C are
given by

$\dot{S}_{\text{A}}^{x}(t)=\Gamma[-D^{z}\\{\alpha
S_{\text{A}}^{x}S_{\text{A}}^{y}(S_{\text{B}}^{x}+S_{\text{C}}^{x})+S_{\text{A}}^{z}(S_{\text{B}}^{x}+S_{\text{C}}^{x})+\alpha\\{(S_{\text{A}}^{y})^{2}+(S_{\text{A}}^{z})^{2})(S_{\text{B}}^{y}+S_{\text{C}}^{y})\\}+D^{y}\\{\alpha((S_{\text{A}}^{y})^{2}+(S_{\text{A}}^{z})^{2})(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}\\\
+S_{\text{E}}^{z})-S_{\text{A}}^{y}(S_{\text{B}}^{x}+S_{\text{C}}^{x})+\alpha
S_{\text{A}}^{x}S_{\text{A}}^{z}(S_{\text{B}}^{x}+S_{\text{C}}^{x})\\}+D^{x}\\{S_{\text{A}}^{y}(\alpha
S_{\text{A}}^{x}(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})+S_{\text{B}}^{y}+S_{\text{C}}^{y})+S_{\text{A}}^{z}(\alpha
S_{\text{A}}^{x}\\{-(S_{\text{B}}^{y}+S_{\text{C}}^{y})\\}\\\
+S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})\\}-(S_{\text{A}}^{y}-\alpha
S_{\text{A}}^{x}S_{\text{A}}^{z})\\{2K_{0}S_{\text{A}}^{z}-B_{0}\sin(\Omega
t)+J_{3}(S_{\text{D}}^{z}+S_{\text{E}}^{z})\\}+J_{1}\\{\alpha
S_{\text{A}}^{x}S_{\text{A}}^{y}S_{\text{B}}^{y}-S_{\text{A}}^{y}(\alpha
S_{\text{A}}^{y}S_{\text{B}}^{x}+S_{\text{B}}^{z})\\\ +\alpha
S_{\text{A}}^{x}S_{\text{A}}^{z}S_{\text{B}}^{z}-\alpha(S_{\text{A}}^{z})^{2}S_{\text{B}}^{x}+S_{\text{A}}^{z}S_{\text{B}}^{y}\\}+J_{4}\\{\alpha
S_{\text{A}}^{x}S_{\text{A}}^{y}S_{\text{C}}^{y}-S_{\text{A}}^{y}(\alpha
S_{\text{A}}^{y}S_{\text{C}}^{x}+S_{\text{C}}^{z})+\alpha
S_{\text{A}}^{x}S_{\text{A}}^{z}S_{\text{C}}^{z}-\alpha(S_{\text{A}}^{z})^{2}S_{\text{C}}^{x}+S_{\text{A}}^{z}S_{\text{C}}^{y}\\}]$
(11) $\dot{S}_{\text{A}}^{y}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{A}}^{x}\\{S_{\text{A}}^{x}(S_{\text{B}}^{x}+S_{\text{C}}^{x})+S_{\text{A}}^{y}(S_{\text{B}}^{y}+S_{\text{C}}^{y})\\}+\alpha(S_{\text{A}}^{z}){}^{2}(S_{\text{B}}^{x}+S_{\text{C}}^{x})-S_{\text{A}}^{z}(S_{\text{B}}^{y}+S_{\text{C}}^{y})\\}+D^{y}\\{S_{\text{A}}^{x}(-\alpha
S_{\text{A}}^{y}(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})\\\
+S_{\text{B}}^{x}+S_{\text{C}}^{x})+S_{\text{A}}^{z}(\alpha
S_{\text{A}}^{y}(S_{\text{B}}^{x}+S_{\text{C}}^{x})+S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})\\}+D^{x}(-\\{\alpha((S_{\text{A}}^{x}){}^{2}+(S_{\text{A}}^{z}){}^{2})(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})+S_{\text{A}}^{x}(S_{\text{B}}^{y}+S_{\text{C}}^{y})+\alpha
S_{\text{A}}^{y}S_{\text{A}}^{z}(S_{\text{B}}^{y}\\\
+S_{\text{C}}^{y})\\})+(S_{\text{A}}^{x}+\alpha
S_{\text{A}}^{y}S_{\text{A}}^{z})\\{2K_{0}S_{\text{A}}^{z}-B_{0}\sin(\Omega
t)+J_{3}(S_{\text{D}}^{z}+S_{\text{E}}^{z})\\}+J_{1}\\{\alpha(S_{\text{A}}^{x}S_{\text{A}}^{y}S_{\text{B}}^{x}-(S_{\text{A}}^{x}){}^{2}S_{\text{B}}^{y}-(S_{\text{A}}^{z}){}^{2}S_{\text{B}}^{y}+S_{\text{A}}^{y}S_{\text{A}}^{z}S_{\text{B}}^{z})+S_{\text{A}}^{x}S_{\text{B}}^{z}\\\
-S_{\text{A}}^{z}S_{\text{B}}^{x}\\}+J_{4}\\{\alpha(S_{\text{A}}^{x}S_{\text{A}}^{y}S_{\text{C}}^{x}-(S_{\text{A}}^{x}){}^{2}S_{\text{C}}^{y}-(S_{\text{A}}^{z}){}^{2}S_{\text{C}}^{y}+S_{\text{A}}^{y}S_{\text{A}}^{z}S_{\text{C}}^{z})+S_{\text{A}}^{x}S_{\text{C}}^{z}-S_{\text{A}}^{z}S_{\text{C}}^{x}\\}]$
(12) $\dot{S}_{\text{A}}^{z}(t)=\Gamma[D^{z}\\{S_{\text{A}}^{x}\\{\alpha
S_{\text{A}}^{z}(S_{\text{B}}^{y}+S_{\text{C}}^{y})+S_{\text{B}}^{x}+S_{\text{C}}^{x}\\}+S_{\text{A}}^{y}\\{-\alpha
S_{\text{A}}^{z}(S_{\text{B}}^{x}+S_{\text{C}}^{x})+S_{\text{B}}^{y}+S_{\text{C}}^{y}\\}\\}-D^{y}\\{\alpha
S_{\text{A}}^{x}\\{S_{\text{A}}^{z}(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})+S_{\text{A}}^{x}(S_{\text{B}}^{x}\\\
+S_{\text{C}}^{x})\\}+S_{\text{A}}^{y}(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})+\alpha(S_{\text{A}}^{y}){}^{2}(S_{\text{B}}^{x}+S_{\text{C}}^{x})\\}+D^{x}\\{\alpha
S_{\text{A}}^{y}\\{S_{\text{A}}^{z}(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})+S_{\text{A}}^{y}(S_{\text{B}}^{y}+S_{\text{C}}^{y})\\}-S_{\text{A}}^{x}(S_{\text{B}}^{z}+S_{\text{C}}^{z}+S_{\text{D}}^{z}+S_{\text{E}}^{z})\\\
+\alpha(S_{\text{A}}^{x}){}^{2}(S_{\text{B}}^{y}+S_{\text{C}}^{y})\\}+\alpha((S_{\text{A}}^{x}){}^{2}+(S_{\text{A}}^{y}){}^{2})\\{-2K_{0}S_{\text{A}}^{z}+B_{0}\sin(\Omega
t)-J_{3}(S_{\text{D}}^{z}+S_{\text{E}}^{z})\\}+J_{1}\\{\alpha
S_{\text{A}}^{z}(S_{\text{A}}^{x}S_{\text{B}}^{x}+S_{\text{A}}^{y}S_{\text{B}}^{y})-\alpha
S_{\text{B}}^{z}((S_{\text{A}}^{x}){}^{2}+(S_{\text{A}}^{y}){}^{2})\\\
+S_{\text{A}}^{x}(-S_{\text{B}}^{y})+S_{\text{A}}^{y}S_{\text{B}}^{x}\\}+J_{4}\\{\alpha
S_{\text{A}}^{z}(S_{\text{A}}^{x}S_{\text{C}}^{x}+S_{\text{A}}^{y}S_{\text{C}}^{y})-\alpha
S_{\text{C}}^{z}((S_{\text{A}}^{x}){}^{2}+(S_{\text{A}}^{y}){}^{2})+S_{\text{A}}^{x}(-S_{\text{C}}^{y})+S_{\text{A}}^{y}S_{\text{C}}^{x}\\}]$
(13) $\dot{S}_{\text{B}}^{x}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{B}}^{y}(S_{\text{A}}^{y}S_{\text{B}}^{y}-S_{\text{B}}^{x}S_{\text{C}}^{x}-S_{\text{B}}^{y}S_{\text{C}}^{y})+\alpha(S_{\text{B}}^{z}){}^{2}(S_{\text{A}}^{y}-S_{\text{C}}^{y})+S_{\text{A}}^{x}(\alpha
S_{\text{B}}^{x}S_{\text{B}}^{y}+S_{\text{B}}^{z})-S_{\text{B}}^{z}S_{\text{C}}^{x}\\}+D^{y}\\{\alpha
S_{\text{B}}^{z}(-S_{\text{A}}^{z}S_{\text{B}}^{z}+S_{\text{B}}^{x}S_{\text{C}}^{x}+S_{\text{B}}^{z}S_{\text{C}}^{z})\\\
+\alpha(S_{\text{B}}^{y}){}^{2}(S_{\text{C}}^{z}-S_{\text{A}}^{z})+S_{\text{A}}^{x}(S_{\text{B}}^{y}-\alpha
S_{\text{B}}^{x}S_{\text{B}}^{z})-S_{\text{B}}^{y}S_{\text{C}}^{x}\\}+D^{x}\\{\alpha
S_{\text{B}}^{x}(S_{\text{A}}^{y}S_{\text{B}}^{z}-S_{\text{A}}^{z}S_{\text{B}}^{y}-S_{\text{B}}^{z}S_{\text{C}}^{y}+S_{\text{B}}^{y}S_{\text{C}}^{z})-S_{\text{A}}^{y}S_{\text{B}}^{y}-S_{\text{A}}^{z}S_{\text{B}}^{z}+S_{\text{B}}^{y}S_{\text{C}}^{y}+S_{\text{B}}^{z}S_{\text{C}}^{z}\\}\\\
+J_{1}\\{S_{\text{A}}^{y}(\alpha
S_{\text{B}}^{x}S_{\text{B}}^{y}+S_{\text{B}}^{z})+S_{\text{A}}^{z}(\alpha
S_{\text{B}}^{x}S_{\text{B}}^{z}-S_{\text{B}}^{y})+\alpha(-S_{\text{A}}^{x})((S_{\text{B}}^{y}){}^{2}+(S_{\text{B}}^{z}){}^{2})\\}+J_{2}\\{\alpha
S_{\text{B}}^{x}S_{\text{B}}^{y}S_{\text{C}}^{y}-S_{\text{B}}^{y}(\alpha
S_{\text{B}}^{y}S_{\text{C}}^{x}+S_{\text{C}}^{z})+\alpha
S_{\text{B}}^{x}S_{\text{B}}^{z}S_{\text{C}}^{z}-\alpha(S_{\text{B}}^{z}){}^{2}S_{\text{C}}^{x}\\\
+S_{\text{B}}^{z}S_{\text{C}}^{y}\\}+(S_{\text{B}}^{y}-\alpha
S_{\text{B}}^{x}S_{\text{B}}^{z})(B_{0}\sin(\Omega
t)-2K_{0}S_{\text{B}}^{z})]$ (14)
$\dot{S}_{\text{B}}^{y}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{B}}^{x}(-S_{\text{A}}^{x}S_{\text{B}}^{x}+S_{\text{B}}^{x}S_{\text{C}}^{x}+S_{\text{B}}^{y}S_{\text{C}}^{y})+\alpha(S_{\text{B}}^{z}){}^{2}(S_{\text{C}}^{x}-S_{\text{A}}^{x})+S_{\text{A}}^{y}(S_{\text{B}}^{z}-\alpha
S_{\text{B}}^{x}S_{\text{B}}^{y})-S_{\text{B}}^{z}S_{\text{C}}^{y}\\}+D^{y}\\{\alpha
S_{\text{B}}^{y}(S_{\text{A}}^{z}S_{\text{B}}^{x}-S_{\text{B}}^{x}S_{\text{C}}^{z}+S_{\text{B}}^{z}S_{\text{C}}^{x})\\\
+S_{\text{A}}^{x}(-(S_{\text{B}}^{x}+\alpha
S_{\text{B}}^{y}S_{\text{B}}^{z}))-S_{\text{A}}^{z}S_{\text{B}}^{z}+S_{\text{B}}^{x}S_{\text{C}}^{x}+S_{\text{B}}^{z}S_{\text{C}}^{z}\\}+D^{x}\\{\alpha\\{S_{\text{B}}^{y}S_{\text{B}}^{z}(S_{\text{A}}^{y}-S_{\text{C}}^{y})+S_{\text{A}}^{z}((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{z}){}^{2})-S_{\text{C}}^{z}((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{z}){}^{2})\\}+S_{\text{B}}^{x}(S_{\text{A}}^{y}-S_{\text{C}}^{y})\\}\\\
+J_{1}\\{\alpha S_{\text{A}}^{x}S_{\text{B}}^{x}S_{\text{B}}^{y}-\alpha
S_{\text{A}}^{y}((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{z}){}^{2})+S_{\text{A}}^{z}(S_{\text{B}}^{x}+\alpha
S_{\text{B}}^{y}S_{\text{B}}^{z})-S_{\text{A}}^{x}S_{\text{B}}^{z}\\}+J_{2}\\{S_{\text{B}}^{x}(\alpha
S_{\text{B}}^{x}(-S_{\text{C}}^{y})+\alpha
S_{\text{B}}^{y}S_{\text{C}}^{x}+S_{\text{C}}^{z})+S_{\text{B}}^{z}(\alpha
S_{\text{B}}^{y}S_{\text{C}}^{z}-S_{\text{C}}^{x})\\\
+\alpha(-(S_{\text{B}}^{z}){}^{2})S_{\text{C}}^{y}\\}+(S_{\text{B}}^{x}+\alpha
S_{\text{B}}^{y}S_{\text{B}}^{z})(-(B_{0}\sin(\Omega
t)-2K_{0}S_{\text{B}}^{z}))]$ (15)
$\dot{S}_{\text{B}}^{z}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{B}}^{z}(S_{\text{A}}^{x}S_{\text{B}}^{y}-S_{\text{A}}^{y}S_{\text{B}}^{x}-S_{\text{B}}^{y}S_{\text{C}}^{x}+S_{\text{B}}^{x}S_{\text{C}}^{y})-S_{\text{A}}^{x}S_{\text{B}}^{x}-S_{\text{A}}^{y}S_{\text{B}}^{y}+S_{\text{B}}^{x}S_{\text{C}}^{x}+S_{\text{B}}^{y}S_{\text{C}}^{y}\\}+D^{y}\\{\alpha(S_{\text{A}}^{x}-S_{\text{C}}^{x})((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{y}){}^{2})\\\
+S_{\text{A}}^{z}(\alpha
S_{\text{B}}^{x}S_{\text{B}}^{z}+S_{\text{B}}^{y})+\alpha
S_{\text{B}}^{x}(-S_{\text{B}}^{z})S_{\text{C}}^{z}-S_{\text{B}}^{y}S_{\text{C}}^{z}\\}+D^{x}\\{\alpha\\{S_{\text{B}}^{y}(-S_{\text{A}}^{z}S_{\text{B}}^{z}+S_{\text{B}}^{y}S_{\text{C}}^{y}+S_{\text{B}}^{z}S_{\text{C}}^{z})-S_{\text{A}}^{y}((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{y}){}^{2})+(S_{\text{B}}^{x}){}^{2}S_{\text{C}}^{y}\\}+S_{\text{B}}^{x}(S_{\text{A}}^{z}-S_{\text{C}}^{z})\\}\\\
+J_{1}\\{S_{\text{A}}^{y}(\alpha
S_{\text{B}}^{y}S_{\text{B}}^{z}-S_{\text{B}}^{x})-\alpha
S_{\text{A}}^{z}((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{y}){}^{2})+S_{\text{A}}^{x}S_{\text{B}}^{y}+\alpha
S_{\text{A}}^{x}S_{\text{B}}^{x}S_{\text{B}}^{z}\\}+J_{2}\\{\alpha
S_{\text{B}}^{z}(S_{\text{B}}^{x}S_{\text{C}}^{x}+S_{\text{B}}^{y}S_{\text{C}}^{y})-\alpha
S_{\text{C}}^{z}((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{y}){}^{2})+S_{\text{B}}^{x}(-S_{\text{C}}^{y})\\\
+S_{\text{B}}^{y}S_{\text{C}}^{x}\\}+\alpha((S_{\text{B}}^{x}){}^{2}+(S_{\text{B}}^{y}){}^{2})(B_{0}\sin(\Omega
t)-2K_{0}S_{\text{B}}^{z})]$ (16)
$\dot{S}_{\text{C}}^{x}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{C}}^{y}(S_{\text{B}}^{y}S_{\text{C}}^{y}-S_{\text{C}}^{x}S_{\text{A}}^{x}-S_{\text{C}}^{y}S_{\text{A}}^{y})+\alpha(S_{\text{C}}^{z}){}^{2}(S_{\text{B}}^{y}-S_{\text{A}}^{y})+S_{\text{B}}^{x}(\alpha
S_{\text{C}}^{x}S_{\text{C}}^{y}+S_{\text{C}}^{z})-S_{\text{C}}^{z}S_{\text{A}}^{x}\\}+D^{y}\\{\alpha
S_{\text{C}}^{z}(-S_{\text{B}}^{z}S_{\text{C}}^{z}+S_{\text{C}}^{x}S_{\text{A}}^{x}+S_{\text{C}}^{z}S_{\text{A}}^{z})\\\
+\alpha(S_{\text{C}}^{y}){}^{2}(S_{\text{A}}^{z}-S_{\text{B}}^{z})+S_{\text{B}}^{x}(S_{\text{C}}^{y}-\alpha
S_{\text{C}}^{x}S_{\text{C}}^{z})-S_{\text{C}}^{y}S_{\text{A}}^{x}\\}+D^{x}\\{\alpha
S_{\text{C}}^{x}(S_{\text{B}}^{y}S_{\text{C}}^{z}-S_{\text{B}}^{z}S_{\text{C}}^{y}-S_{\text{C}}^{z}S_{\text{A}}^{y}+S_{\text{C}}^{y}S_{\text{A}}^{z})-S_{\text{B}}^{y}S_{\text{C}}^{y}-S_{\text{B}}^{z}S_{\text{C}}^{z}+S_{\text{C}}^{y}S_{\text{A}}^{y}+S_{\text{C}}^{z}S_{\text{A}}^{z}\\}\\\
+J_{1}\\{S_{\text{B}}^{y}(\alpha
S_{\text{C}}^{x}S_{\text{C}}^{y}+S_{\text{C}}^{z})+S_{\text{B}}^{z}(\alpha
S_{\text{C}}^{x}S_{\text{C}}^{z}-S_{\text{C}}^{y})+\alpha(-S_{\text{B}}^{x})((S_{\text{C}}^{y}){}^{2}+(S_{\text{C}}^{z}){}^{2})\\}+J_{1}\\{\alpha
S_{\text{C}}^{x}S_{\text{C}}^{y}S_{\text{A}}^{y}-S_{\text{C}}^{y}(\alpha
S_{\text{C}}^{y}S_{\text{A}}^{x}+S_{\text{A}}^{z})+\alpha
S_{\text{C}}^{x}S_{\text{C}}^{z}S_{\text{A}}^{z}-\alpha(S_{\text{C}}^{z}){}^{2}S_{\text{A}}^{x}\\\
+S_{\text{C}}^{z}S_{\text{A}}^{y}\\}+(S_{\text{C}}^{y}-\alpha
S_{\text{C}}^{x}S_{\text{C}}^{z})(B_{0}\sin(\Omega
t)-2K_{0}S_{\text{C}}^{z})]$ (17)
$\dot{S}_{\text{C}}^{y}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{C}}^{x}(-S_{\text{B}}^{x}S_{\text{C}}^{x}+S_{\text{C}}^{x}S_{\text{A}}^{x}+S_{\text{C}}^{y}S_{\text{A}}^{y})+\alpha(S_{\text{C}}^{z}){}^{2}(S_{\text{A}}^{x}-S_{\text{B}}^{x})+S_{\text{B}}^{y}(S_{\text{C}}^{z}-\alpha
S_{\text{C}}^{x}S_{\text{C}}^{y})-S_{\text{C}}^{z}S_{\text{A}}^{y}\\}+D^{y}\\{\alpha
S_{\text{C}}^{y}(S_{\text{B}}^{z}S_{\text{C}}^{x}-S_{\text{C}}^{x}S_{\text{A}}^{z}+S_{\text{C}}^{z}S_{\text{A}}^{x})\\\
+S_{\text{B}}^{x}(-(S_{\text{C}}^{x}+\alpha
S_{\text{C}}^{y}S_{\text{C}}^{z}))-S_{\text{B}}^{z}S_{\text{C}}^{z}+S_{\text{C}}^{x}S_{\text{A}}^{x}+S_{\text{C}}^{z}S_{\text{A}}^{z}\\}+D^{x}\\{\alpha\\{S_{\text{C}}^{y}S_{\text{C}}^{z}(S_{\text{B}}^{y}-S_{\text{A}}^{y})+S_{\text{B}}^{z}((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{z}){}^{2})-S_{\text{A}}^{z}((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{z}){}^{2})\\}+S_{\text{C}}^{x}(S_{\text{B}}^{y}-S_{\text{A}}^{y})\\}\\\
+J_{1}\\{\alpha S_{\text{B}}^{x}S_{\text{C}}^{x}S_{\text{C}}^{y}-\alpha
S_{\text{B}}^{y}((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{z}){}^{2})+S_{\text{B}}^{z}(S_{\text{C}}^{x}+\alpha
S_{\text{C}}^{y}S_{\text{C}}^{z})-S_{\text{B}}^{x}S_{\text{C}}^{z}\\}+J_{1}\\{S_{\text{C}}^{x}(\alpha
S_{\text{C}}^{x}(-S_{\text{A}}^{y})+\alpha
S_{\text{C}}^{y}S_{\text{A}}^{x}+S_{\text{A}}^{z})+S_{\text{C}}^{z}(\alpha
S_{\text{C}}^{y}S_{\text{A}}^{z}-S_{\text{A}}^{x})\\\
+\alpha(-(S_{\text{C}}^{z}){}^{2})S_{\text{A}}^{y}\\}+(S_{\text{C}}^{x}+\alpha
S_{\text{C}}^{y}S_{\text{C}}^{z})(-(B_{0}\sin(\Omega
t)-2K_{0}S_{\text{C}}^{z}))]$ (18)
$\dot{S}_{\text{C}}^{z}(t)=\Gamma[D^{z}\\{\alpha
S_{\text{C}}^{z}(S_{\text{B}}^{x}S_{\text{C}}^{y}-S_{\text{B}}^{y}S_{\text{C}}^{x}-S_{\text{C}}^{y}S_{\text{A}}^{x}+S_{\text{C}}^{x}S_{\text{A}}^{y})-S_{\text{B}}^{x}S_{\text{C}}^{x}-S_{\text{B}}^{y}S_{\text{C}}^{y}+S_{\text{C}}^{x}S_{\text{A}}^{x}+S_{\text{C}}^{y}S_{\text{A}}^{y}\\}+D^{y}\\{\alpha(S_{\text{B}}^{x}-S_{\text{A}}^{x})((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{y}){}^{2})\\\
+S_{\text{B}}^{z}(\alpha
S_{\text{C}}^{x}S_{\text{C}}^{z}+S_{\text{C}}^{y})+\alpha
S_{\text{C}}^{x}(-S_{\text{C}}^{z})S_{\text{A}}^{z}-S_{\text{C}}^{y}S_{\text{A}}^{z}\\}+D^{x}\\{\alpha\\{S_{\text{C}}^{y}(-S_{\text{B}}^{z}S_{\text{C}}^{z}+S_{\text{C}}^{y}S_{\text{A}}^{y}+S_{\text{C}}^{z}S_{\text{A}}^{z})-S_{\text{B}}^{y}((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{y}){}^{2})+(S_{\text{C}}^{x}){}^{2}S_{\text{A}}^{y}\\}+S_{\text{C}}^{x}(S_{\text{B}}^{z}\\\
-S_{\text{A}}^{z})\\}+J_{1}\\{S_{\text{B}}^{y}(\alpha
S_{\text{C}}^{y}S_{\text{C}}^{z}-S_{\text{C}}^{x})-\alpha
S_{\text{B}}^{z}((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{y}){}^{2})+S_{\text{B}}^{x}S_{\text{C}}^{y}+\alpha
S_{\text{B}}^{x}S_{\text{C}}^{x}S_{\text{C}}^{z}\\}+J_{1}\\{\alpha
S_{\text{C}}^{z}(S_{\text{C}}^{x}S_{\text{A}}^{x}+S_{\text{C}}^{y}S_{\text{A}}^{y})-\alpha
S_{\text{A}}^{z}((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{y}){}^{2})\\\
+S_{\text{C}}^{x}(-S_{\text{A}}^{y})+S_{\text{C}}^{y}S_{\text{A}}^{x}\\}+\alpha((S_{\text{C}}^{x}){}^{2}+(S_{\text{C}}^{y}){}^{2})(B_{0}\sin(\Omega
t)-2K_{0}S_{\text{C}}^{z})]$ (19)

where, $\Gamma=\frac{\gamma}{1+\alpha^{2}}$ and $\\{D^{x},D^{y},D^{z}\\}$ are
the components of DMI vector.

##  Appendix B Explicit form of the coefficients $\mathcal{A}_{i}$ and $\xi$

An explicit form of the coefficients
$\mathcal{A}_{1},\mathcal{A}_{2},\mathcal{A}_{3},\mathcal{A}_{4},\mathcal{A}_{5},\mathcal{A}_{6}$
and $\xi$ are expressed below

$\mathcal{A}_{1}=-1-3\alpha^{2}-3\alpha^{4}-\alpha^{6}$ (20)
$\mathcal{A}_{2}=-4i\alpha\gamma\left(1+\alpha^{2}\right)^{2}(J_{4}+3K)$ (21)
$\mathcal{A}_{3}=2(\alpha^{2}+1)\gamma^{2}(\alpha^{2}((D^{z})^{2}-2J_{1}^{2}+J_{1}(J_{4}-J_{2})+J_{4}(J_{2}\\\
+4J_{4})+20J_{4}K+30K^{2})+3(D^{z})^{2}+2J_{1}^{2}+J_{1}J_{2}\\\
-J_{1}J_{4}+2J_{2}^{2}-J_{2}J_{4}+2J_{4}^{2}+4J_{4}K+6K^{2})$ (22)
$\mathcal{A}_{4}=-4i\alpha\gamma^{3}\\{(D^{z})^{2}(\alpha^{2}(J_{1}-4K)+J_{1}+2J_{2}-2J_{4}-8K\\}-J_{1}^{2}\\{\alpha^{2}(J_{2}-3J_{4}-8K)-3J_{2}+J_{4}\\}+\alpha^{2}J_{1}\\{J_{2}(J_{4}\\\
+4K)-2J_{4}(J_{4}+2K)\\}+J_{1}J_{2}(2J_{2}-3J_{4})-2(J_{2}^{2}(J_{4}+2K)+J_{4}^{3}+6J_{4}^{2}K+12J_{4}K^{2}+12K^{3})-2\alpha^{2}\\{J_{4}^{2}(J_{2}+J_{4})\\\
+2J_{4}K(J_{2}+4J_{4})+20J_{4}K^{2}+20K^{3})\\}$ (23)
$\mathcal{A}_{5}=\gamma^{4}[-(\alpha^{2}+9)(D^{z})^{4}+2(D^{z})^{2}(2(\alpha^{2}-3)J_{1}^{2}+J_{1}(\alpha^{2}(J_{2}+3J_{4}+12K)+J_{2}+3J_{4}+4K)-\alpha^{2}\\{J_{4}(J_{2}-2J_{4})+24K^{2}\\}\\\
+7J_{2}J_{4}+8J_{2}K-6J_{4}^{2}-8J_{4}K-16K^{2})-4(\alpha^{2}+1)J_{1}^{4}-4(\alpha^{2}+1)J_{1}^{3}(J_{2}-J_{4})+J_{1}^{2}(-\alpha^{2}(J_{2}^{2}+24K(J_{2}-3J_{4})+2J_{2}J_{4}-15J_{4}^{2}\\\
-96K^{2})-9J_{2}^{2}+14J_{2}J_{4}+24J_{2}K-9J_{4}^{2}-8J_{4}K)+2J_{1}(J_{2}^{2}(\alpha^{2}J_{4}+J_{4}+8K)+\alpha^{2}J_{2}(J_{4}^{2}+12J_{4}K+24K^{2})-J_{2}J_{4}(7J_{4}+12K)\\\
+2J_{4}^{3}-6\alpha^{2}J_{4}(J_{4}+2K)^{2})-9J_{2}^{2}J_{4}^{2}-\alpha^{2}\\{J_{2}^{2}J_{4}^{2}+12J_{2}J_{4}(J_{4}+2K)^{2}+4(J_{4}^{4}+12J_{4}^{3}K+48J_{4}^{2}K^{2}+80J_{4}K^{3}+60K^{4})\\}-16J_{2}^{2}J_{4}K\\\
-16J_{2}^{2}K^{2}+4J_{2}J_{4}^{3}-4J_{4}^{4}-16J_{4}^{3}K-48J_{4}^{2}K^{2}-64J_{4}K^{3}-48K^{4}]$
(24)
$\mathcal{A}_{6}=-4i\alpha\gamma^{5}\\{(D^{z})^{2}-2J_{1}^{2}+J_{1}(J_{4}-J_{2})+J_{4}(J_{2}+2J_{4})+8J_{4}K+12K^{2}\\}(-(D^{z})^{2}(J_{1}+J_{4}-K)+J_{1}^{2}(J_{2}-J_{4}-2K)-J_{1}J_{2}K\\\
+J_{1}J_{4}(J_{4}+K)+J_{2}J_{4}^{2}+J_{4}K(J_{2}+2J_{4})+4J_{4}K^{2}+4K^{3})$
(25)
$\xi=4\gamma^{6}[-(D^{z})^{2}(J_{1}+J_{4}-K)+J_{1}^{2}(J_{2}-J_{4}-2K)-J_{1}J_{2}K+J_{1}J_{4}(J_{4}+K)+J_{2}J_{4}^{2}+J_{4}K(J_{2}+2J_{4})+4J_{4}K^{2}+4K^{3}]^{2}$
(26)

## References

  * (1) N. J. Ghimire and I. I. Mazin, Nat. Mater., 19, 137–138 (2020).
  * (2) A. Mielke, J. Phys. A: Math. Gen., 25 4335 (1992).
  * (3) E. Khatami, R. R. P. Singh, and M. Rigol, Phys. Rev. B, 84, 224411 (2011).
  * (4) P. Saha et. al. Phys. Rev. B, 103, 224402 (2021).
  * (5) A. Farhan et. al. Phys. Rev. B, 96, 064409 (2017).
  * (6) S. E. Hog et. al., 563, 169920 (2022).
  * (7) T. H. Han, et al. Nature, 492, 406–410 (2012).
  * (8) S. Yan, D. A. Huse, and S. R. White, Science, 332, 1173–1176 (2011).
  * (9) A. K. Nayak, et al. Sci. Adv., 2, e1501870 (2016).
  * (10) P. Mendels et. al. Phys. Rev. Lett., 98, 077204 (2007).
  * (11) Y. Li et. al. Phys. Rev. B, 101, 161115(R) (2020).
  * (12) A. Zorko et. al. Phys. Rev. B, 99 214449 (2019).
  * (13) Q. Du et. al., Phys. Rev. Lett. 129, 236601 (2022).
  * (14) J. X. Yin et. al., Nature 562, 91–95 (2018).
  * (15) I. Dzyaloshinsky, Phys. Chem. Solids, 4 (4), 241 (1958).
  * (16) T. Moria Phys. Rev., 120 (1), 91 (1960).
  * (17) Y. H. Liu and Y.Q. Li Chinese Phys. B, 24 017506 (2015).
  * (18) Z. Chen et. al. Phys. Rev. Appl., 17 L011002 (2022).
  * (19) O. Petrova and O. Tchernyshyov, Phys. Rev. B 84, 214433 (2011).
  * (20) M. K. Lee and M. Mochizuki, Phys. Rev. Appl. 18, 014074 (2022).
  * (21) I. K. Drozdov et. al. Phys. Rev. Lett. 104, 136804 (2010).
  * (22) R. Wiesendanger, Nat. Rev. Mater., 1, 16044 (2016).
  * (23) X. Zhang et. al. J. Phys.: Condens. Matter, 32 143001 (2020).
  * (24) C. H. Marrows and K. Zeissler Appl. Phys. Lett., 119, 250502 (2021)
  * (25) B. Dupé, et. al. Nat Commun., 7, 11779 (2016).
  * (26) Y. Zhou, Natl. Sci. Rev., 6 (2) 210–212 (2019).
  * (27) A. Fert, N Reyren and V. Cros, Nat. Rev. Mater, 2, 17031 (2017).
  * (28) Y. Ishida and K. Kondo, Jpn. J. Appl. Phys., 59 SGGI04 (2020).
  * (29) H. Zhang et. al., Mater. Futures, 2 032201 (2023).
  * (30) Y. Tokura and N. Kanazawa, Chem. Rev., 121 (5), 2857–2897 (2021).
  * (31) S. G. Je, et. al. Phys. Rev. B, 88, 214401 (2013).
  * (32) F. Ajejas, et. al. Appl. Phys. Lett., 111, 202402 (2017).
  * (33) J. Brandão, et. al. Sci. Rep., 7, 4569 (2017).
  * (34) W. Wang, et. al. Phys. Rev. Lett., 114, 087203 (2015).
  * (35) M. D. DeJong and K. L. Livesey Phys. Rev. B, 95, 054424 (2017).
  * (36) Y. Ren et. al., Rep. Prog. Phys. 79, 066501 (2016).
  * (37) A. Gao et al., Nature 595, 521–525 (2021).
  * (38) M. Huang et. al., Nano Lett. , 21, 10, 4280–4286 (2021).
  * (39) Z. F. Wang, Z. Liu, and F. Liu, Phys. Rev. Lett. 110, 196801 (2013)
  * (40) S. Acharjee et. al., J. Appl. Phys., 135(2), 023902 (2024).
  * (41) S. Acharjee et. al., Supercond. Sci. Technol. 36(12), 125014 (2023).
  * (42) N. Dutta et. al., Phys. B: Condens. Matter., 669, 415262 (2023).
  * (43) S. Acharjee et. al., J. Magn. Magn. Mater., 572, 170579 (2023).
  * (44) S. Acharjee et. al., Chaos, 33(1), 013136 (2023).
  * (45) H. F. Liu, Y. Z. Yang, Z. H. Dai and Z. H. Yu, Chaos 13, 839–844 (2003).
  * (46) A. Gottscholl et. al. Nat. Commun. 12, 4480 (2021).
  * (47) M. A. Khan et. al. Eng. Res. Express 3, 02200 (2021).
  * (48) C. Bonato et. al. Nature Nanotech 11, 247–252 (2016).
  * (49) E. L. Wong et. al. Nano Lett. 20(5) 3284–3290 (2020).
  * (50) C. J. Elkins and M. T. Alley, Exp. Fluids 43, 823–858 (2007).
  * (51) A. F. Bakuzis, P. C. Morais and F. A. Tourinho, J. Magn. Reson. Series A 122, 100–103 (1996).
  * (52) A. N. Anisimov et. al., Phys. Rev. Lett. 82, 2390 (1999).
  * (53) H. Y. Yuan, et al., EPL, 126 67006 (2019).
  * (54) R. Chen, et al., Nat. Commun. 12, 3113 (2021).
  * (55) H. Bauer et. al. Nat. Commun. 6, 8274 (2015).
  * (56) V. Flovik, B. H. Pettersen and Erik Wahlström, J. Appl. Phys. 119, 163903 (2016).
  * (57) I. Bertelli et. al. Sci. Adv. 6, eabd3556 (2020).
  * (58) S. Crooker et. al. Nature 431, 49–52 (2004).
  * (59) M. Sandri, Mathematica Journal. 6(3), 78-84 (1996).
  * (60) R. C. Wolff, J. R. Stat. Soc. Series B. 54(2), 353-71 (1992).


# Expressing turbulent kinetic energy as coarse-grained enstrophy or strain
deformations

Damiano Capocci<EMAIL_ADDRESS>Department of Physics and INFN,
University of Rome Tor Vergata

###### Abstract

In turbulent flows, the element of fluid gets deformed by chaotic motion due
to the formation of sharp velocity gradients. A direct connection between the
element of fluid stresses and the energy balance still remains elusive. Here,
an exact identity of incompressible turbulence is derived linking the velocity
gradient norm across the scales to the total kinetic energy. In the context of
three dimensional (3D) homogeneous turbulence, this relation can be
specialised obtaining the expression of total kinetic energy decomposed either
in terms of fluid element deformations due to strain motion or via the
resolved scale enstrophy. Applied to data from direct numerical simulations
(DNS), the decomposition reveals that, beyond the scales dominated by the
external forcing, extensional and contractile deformations account
approximately for $55\%$ and $40\%$ of the kinetic energy of the associated
scale saturating to these values at the small scales; less than the remaining
$5\%$ is carried by the indefinite-type ones. From both these two identities
one can derive an exact expression of the omni-directional energy spectrum
which is formulated via an ultraviolet-infrared mixing of scales. Data from
DNS indicated that this formulation of the energy spectrum reproduces both
analytically and numerically the power-law behaviour of the well known
Kolmogorov spectral scaling.

††preprint: APS/123-QED

Among the countless scientific contexts that deals with turbulent flows, two
characteristics are recurrent: the nonlinearity and the multiscale nature of
the velocity field leading to a spatio-temporal chaotic dynamics. Turbulence
is a common emerging property in water or air when the value of the viscosity
is _small_ compared to the related typical sizes and velocities. A turbulent
flow is characterised by the formation of small scale vortical motion from a
large scale injection of energy; this phenomenon is usually defined as the
_cascade_ of energy across the scales. This energy transfer is responsible of
the production of velocity gradient (VG) at every scale that has an impact on
the deformations that the element of fluid undergoes. To the author’s
knowledge, it is still unclear the functional dependence between the total
kinetic energy and both the intensity of strain and rotational motions.

In this work, two expressions for the total kinetic energy are proved in
homogeneous turbulence: one written in terms of the rate of strain deformation
of the element of fluid and the other based on the squared vorticity commonly
referred as _enstrophy_. As an extension of the latter, a corresponding exact
identity for the kinetic energy spectrum is derived through which it is
possible to connect the energy Fourier mode with the filtered enstrophy. The
derived identities are verified via data from direct numerical simulation
(DNS) and used to reveal the percentages of the purely contractile, the purely
extensional and the _indefinite-type_ stresses directions that constitute the
kinetic energy across the scales and the validity of the novel representation
of the kinetic energy spectrum.

The Navier-Stokes equation describing the evolution of an incompressible
velocity flow in an arbitrary dominion $\mathcal{D}$ is the following:

$\partial_{t}u_{i}=-u_{j}\dfrac{\partial u_{i}}{\partial
x_{j}}-\dfrac{1}{\rho}\dfrac{\partial p}{\partial
x_{i}}+\nu\nabla^{2}u_{i}+f_{i}$ (1)

where $\nu$ is the kinematic viscosity and $\rho$ the mass density. The linear
_viscous term_ $\nu\nabla^{2}u_{i}$ corresponds to the dissipation, being
related to the conversion of kinetic energy into heat. The pressure field
$p(\bm{x},t)$ enforces the incompressibility constraint $\partial_{i}u_{i}=0$.
The velocity gradient tensor $u_{i,j}$ describes the flow geometry in terms of
strain rate $S_{ij}=(u_{i,j}+u_{j,i})/2$ and rotation rate
$\Omega_{ij}=(u_{i,j}-u_{j,i})/2$, which can also be expressed as the
aforementioned vorticity vector $\omega_{i}=\varepsilon_{ijk}\Omega_{jk}$. Eq.
(1) admits one adimensional control parameter which is the Reynolds number
$Re=U\,L/\nu$, where $U$ and $L$ are the characteristic velocity and length-
scale of the flow, respectively.

By applying a low-pass filter, we can think of separate the large-scale
dynamics from the small scale ones [1]:

$\overline{u}^{\ell}_{i}=G^{\ell}*u_{i},\quad\quad\tilde{\overline{u}}^{\ell}_{i}={\tilde{G}^{\ell}}\cdot\tilde{u}_{i}$
(2)

where $\tilde{(\cdot)}$ indicates the Fourier transform and $*$ the
convolution product. As regards the properties of $G^{\ell}$, we require that
it is an even function with volume average equal to unity.

The application of a low-pass filter to eq. (1) provides:

$\partial_{t}\overline{u}^{\ell}_{i}=-\overline{u}^{\ell}_{j}\,\dfrac{\partial\overline{u}_{i}^{\ell}}{\partial
x_{j}}-\dfrac{1}{\rho}\dfrac{\partial\overline{p}^{\ell}}{\partial
x_{i}}+\nu\nabla^{2}\overline{u}_{i}^{\ell}+\overline{f}^{\ell}_{i}-\partial_{j}\tau^{\ell}_{ij}$
(3)

where
$\tau_{ij}^{\ell}=\overline{u_{i}u_{j}}^{\ell}-\overline{u}^{\ell}_{i}\overline{u}^{\ell}_{j}$
is the so-called Subgrid scale stress (SGS) tensor, representing the effective
stress due to the features smaller than $\ell$ on the resolved-scale velocity.
Therefore, we can define the kinetic energy related to scales larger than
$\ell$ i.e. the _resolved-scale_ energy
$E_{K}^{\ell}(\bm{x},t)=\frac{1}{2}\overline{u}^{\ell}_{i}(\bm{x},t)\,\overline{u}^{\ell}_{i}(\bm{x},t)$
and the kinetic energy at scales smaller than $\ell$ which is
$e^{\ell}(\bm{x},t)=\frac{1}{2}\tau^{\ell}_{ii}(\bm{x},t)$ known as the
_unresolved-scale_ energy. Note that the sum of these two quantities
determines the total kinetic energy:

$\dfrac{1}{2}\overline{u_{i}\,u_{i}}^{\ell}=\dfrac{1}{2}\overline{u}^{\ell}_{i}\overline{u}^{\ell}_{i}+\dfrac{1}{2}\tau^{\ell}_{ii}$
(4)

where we omit the spatio-temporal dependence. If the filtered kernel
$G^{\ell}(\bm{x})\geq 0$ $\forall\bm{x}\in\mathcal{D}$, then
$\tau_{ii}(\bm{x})\geq 0$ inside $\mathcal{D}$ [2], hence $\tau_{ii}$ can be
correctly intepreted as an energy. We also remind that
$\overline{u}^{\ell=0}(\bm{x},t)=u(\bm{x},t)$ implying that $\tau^{\ell=0}=0$.

At this stage, we note that for an arbitrary smooth function $f$, the relation
$\langle\overline{f}^{\ell}\rangle=\langle f\rangle$ is a consequence of the
property of our (generic) filter kernel. This observation becomes relevant for
the volume average of (4). In fact, as already noted by [3], we can write:

$\langle E_{K}\rangle=\dfrac{1}{2}\langle
E^{\ell}_{K}\rangle+\dfrac{1}{2}\langle\tau^{\ell}_{ii}\rangle$ (5)

where the filtering parameter belongs only to the RHS resembling the
_subtraction point_ of the renormalization theory.

At this point, for the reason expressed above, we center on positive-definite
filter kernels by specialising our discussion on the application of the
Gaussian kernel:

$G^{\ell}(\bm{r})=\left(2\pi\ell^{2}\right)^{-3/2}e^{-\bm{r}^{2}/2\ell^{2}},\quad\quad{\tilde{G}^{\ell}}(\bm{k})=e^{-\bm{k}^{2}\ell^{2}/2}$
(6)

As a specific feature of the Gaussian filter, there exists an exact
decomposition of the SGS stresses [4][5] that reads:

$\tau_{ij}^{\ell}=\int_{0}^{\ell^{2}}\,\mathrm{d}\theta\,\overline{\,\overline{u}^{\sqrt{\theta}}_{i,j}\,\overline{u}^{\sqrt{\theta}}_{i,j}}^{\sqrt{\ell^{2}-\theta}}$
(7)

expressing SGS stresses 111Note that the methodology of [4] requires
homogeneity or periodic boundaries in order to set to zero the boundary terms
and not for the derivation of the decomposition (7) itself.  as the sum of the
contributions of VG fields from the scales $\leq\ell$ that are in turn
projected into the larger scales by the complementary filter
$\sqrt{\ell^{2}-\theta}$. In this way we can recast eq. (5) into:

$\langle E_{K}\rangle=\dfrac{1}{2}\langle
E^{\ell}_{K}\rangle+\dfrac{1}{2}\int_{0}^{\ell^{2}}\,\mathrm{d}\theta\,\langle\overline{u}^{\sqrt{\theta}}_{i,j}\,\overline{u}^{\sqrt{\theta}}_{i,j}\rangle$
(8)

where, by virtue of the property expressed above, we can ignore the
complementary filtering at the scale $\sqrt{\ell^{2}-\theta}$ in the average
of (7) . In addition, the exchange of the volume averaging and the integral
over the scales is granted by the assumption that VG are bounded in the domain
$\mathcal{D}$ and that the system has finite energy. If we let our system to
be characterised by infinite scales of energy, this would immediately imply
that $\displaystyle\lim_{\ell\to\infty}\langle E^{\ell}_{K}\rangle=0$. This
prescription appears as a reasonable assumption for any classical physics
system which is well-defined in its range of scales 222The length contraction
in general relativity would already pose conceptual difficulties to the
present filtering approach.. A detailed and more generalised analysis on both
space and scale convergence of (8) will be provided elsewhere 333However, in
the applications, the interest obviously focuses on systems characterised by a
maximum allowed scale for which the integration upper bound of (7) does not
diverge ensuring a rigorous exchange of integrals to derive (8).. By applying
the previous considerations and calculating the limit for infinite scale in
the RHS of (8), we get:

$\langle
E_{K}\rangle=\dfrac{1}{2}\int_{0}^{\infty}\,\mathrm{d}\theta\,\langle\overline{u}^{\sqrt{\theta}}_{i,j}\,\overline{u}^{\sqrt{\theta}}_{i,j}\rangle$
(9)

that can be even used to reformulate (8) into:

$\langle
E_{K}^{\ell}\rangle=\dfrac{1}{2}\int_{\ell^{2}}^{\infty}\,\mathrm{d}\theta\,\langle\overline{u}^{\sqrt{\theta}}_{i,j}\,\overline{u}^{\sqrt{\theta}}_{i,j}\rangle.$
(10)

Eq. (9), and . (10) the resolved scale counterpart, can be conceived as an
exact gradient-tensor decomposition of kinetic energy.

At this point, in order to quantify eq. (10) we consider data from DNS. Our
database describes stationary homogeneous and isotropic turbulence [9] where
eq. (3) is evolved in a triply-periodic box with a stochastic forcing active
in the wavenumber band $k\in[0.5,2,4]$ where $k=\pi/\ell$ whose midpoint
determines the characteristic forcing scale $L_{f}$. In the following
quantification, we are relying on one instantaneous configuration of the
velocity field. In particular, fig. 1 compares the _usual_ expression of the
resolved-scale kinetic energy (i.e. the first addendum of the RHS of (4)) with
(10). Both these two expressions are normalised by the total kinetic energy
and displayed as a function of the adimensional parameter $L_{f}/\ell$. In
consequence of the validity of (10), we clearly observe that the two displayed
profiles are definitely indistinguishable.

Figure 1: Resolved-scale mean kinetic energy normalised by the unfiltered
value ${\langle E^{\ell}_{K}\rangle}/{\langle E_{K}\rangle}$ as a function of
the adimenesional parameter $L_{f}/\ell$ calculated from the standard
definition and from eq. (10). The two curves clearly overlap.

As a conclusion of this preliminary analysis, it is worth to underline that in
homogeneous turbulence, the first Betchov relation [10]:

$\langle\overline{S}^{\ell}_{ij}\overline{S}^{\ell}_{ij}\rangle=\dfrac{1}{2}\langle\,||\overline{\bm{\omega}}^{\ell}||^{2}\,\rangle,$
(11)

holding $\forall\ell$, implies that the mean squared VG can be either written
uniquely in terms of the (filtered) squared strain or via the (filtered)
enstrophy.

##  I Mean kinetic energy in terms of fluid element deformations

In this section, as introduced above, we are going to express the mean squared
VG appearing in the RHS of (9) only in terms of the resolved scale strain-rate
squared via the application of (11). Before proceeding to this direction, it
becomes important to remark that the strain-rate tensor, being symmetric, has
3 real eigenvalues $\lambda_{i}$ with $i\in[1,2,3]$ related to a set of
orthogonal eigenvectors. Moreover incompressibility imposes the pointwise
constraint
$\overline{\lambda}^{\ell}_{1}+\overline{\lambda}^{\ell}_{2}+\overline{\lambda}^{\ell}_{3}=0$.
As a consequence, the largest eigenvalue $\overline{\lambda}^{\ell}_{1}$ is
always positive hence extensional along the direction of its eigenvector while
the smallest one $\overline{\lambda}^{\ell}_{3}$ is always negative thus
related to a contractile direction and the intermediate
$\overline{\lambda}^{\ell}_{2}$ can be either positive or negative from which
the abovementioned definition _indefinite-type_.

Being the trace an invariant, we can express (10) via the (sorted) strain-rate
tensor eigenvalues obtaining:

$\displaystyle\langle E_{K}^{\ell}\rangle=$
$\displaystyle\int_{\ell^{2}}^{\infty}\,\mathrm{d}\theta\,\langle\,(\overline{\lambda}^{\sqrt{\theta}}_{1})^{2}+(\overline{\lambda}^{\sqrt{\theta}}_{2})^{2}+(\overline{\lambda}^{\sqrt{\theta}}_{3})^{2}\,\rangle$
(12) $\displaystyle=$
$\displaystyle\,P_{1}^{\ell}+P_{2}^{\ell}+P_{3}^{\ell}\,.$ (13)

where each $P_{i}^{\ell}$ represents the contribution of the $i$-th strain-
rate eigenvalue to the mean resolved-scale kinetic energy. In particular the
reference frame where $S_{ij}$ is diagonal is often called as _principal axis
frame_. The above equation goes beyond the usual representation of the kinetic
energy of the fluid element as the squared Euclidean norm velocity obtaining
the mean energy budget as the sum across the scales of the stresses intensity
along the strain-rate principal axes.

By applying eq. (12) to the already mentioned data, we can quantify the role
of each $P_{i}^{\ell}$ across the scales. in fact, fig. 2 displays the
adimensional $P_{i}^{\ell}$ normalised by $\langle E^{\ell}_{K}\rangle$ as a
function of the adimensional parameter $L_{f}/\ell$ where the former describes
the percentage of each $P_{i}^{\ell}$ in the kinetic energy balance at the
scale $\ell$. First of all we notice that, apart from the region
$L_{f}/\ell\leq 1$ where the flow geometry is governed by the external
forcing, each $P_{i}^{\ell}$ is definitely monotonic where
$\overline{\lambda}^{\ell}_{1}$ and $\overline{\lambda}^{\ell}_{2}$ increase
while $\overline{\lambda}^{\ell}_{3}$ decreases as we progress towards the
_small_ scales where each contribution saturates. In particular for these
scales the purely extensional motion carries about $54\%$ of the corresponding
resolved-scale energy while the purely contractile directions approximately
$41\%$ while the remaining part is due to $\overline{\lambda}^{\ell}_{2}$. The
impact of the intermediate eigenvalue in the present context is ostensibly
different from the one in the _energy flux_ where
$\overline{\lambda}^{\ell}_{2}$ determines the direction of the energy cascade
[5][11] while, from the above quantification, we conclude that in the energy
balance scenario its contribution accounts for less than $5\%$ 444This
quantification has been repeated on datasets with different $Re$ like HL1 from
[9]. We obtain qualitative the same saturation profile and the same asymptotes
as those displayed in fig. 2..

Figure 2: Strain-rate eigenvalues contributions $P_{i}^{\ell}$ normalised by
the resolved scale as a function of $L_{f}/\ell$. Both the axes show
adimensional quantities. The values of the curves sum to unity scale by scale.
The displayed percentages refer to the smallest resolved scale. The horizontal
thin lines, with their corresponding $y$-axis values, are added manually to
emphasize the (rounded) values at which $P_{i}^{\ell}$ saturates.

##  II Connecting resolved scale energy and enstrophy

In this section, as a complementary procedure, we are going to express the
mean squared VG appearing in the RHS of (9) only in terms of the resolved
scale enstrophy using (11). Hence, in homogeneous turbulence, (10) can be
simplified into:

$\langle
E_{K}^{\ell}\rangle=\dfrac{1}{2}\int_{\ell^{2}}^{\infty}\,\mathrm{d}\theta\,\langle\,||\overline{\bm{\omega}}^{\sqrt{\theta}}||^{2}\,\rangle.$
(14)

where if the mean resolved scale enstrophy is known in real space at each
scale, we can retrieve the mean kinetic energy. Clearly we recover the
expression of $\langle E_{K}\rangle$, like in (9), by setting to zero the
lower bound of the previous integral. This identity becomes even more
interesting in 2D turbulence where the enstrophy
$\zeta=\frac{1}{2}\langle||\overline{\bm{\omega}}||^{2}\rangle$ is an inviscid
conserved quantity, in this respect the above equation connects two conserved
quantities which is far from being trivial. It is even less trivial that in 3D
flows a conserved (inviscid) quantity can be expressed via the coarse-graining
of a non-conserved one.

However, it is possible to show that eq. (14), and in consequence (12),
setting $\ell=0$ are more general. For simplicity of notation, let us assume
that the velocity fields allows the continuous Fourier representation (see
[13] for details). In this regards, we consider the RHS of eq. (14) where we
consider that the filtering on the vorticity is implemented by a generic
filter kernel $G^{\ell}$ instead of Gaussian kernel on which the methodology
hereby presented is based. In this manner, the RHS of eq. (14) downgrades to:

$\displaystyle\dfrac{1}{2}\int_{\ell^{2}}^{\infty}\,\mathrm{d}\theta\,\langle\,||\overline{\bm{\omega}}^{\sqrt{\theta}}||^{2}\,\rangle=$
$\displaystyle\dfrac{1}{2}\int_{\ell^{2}}^{\infty}\mathrm{d}\theta\int_{\mathbb{R}}\,\,\,\mathrm{d}k\,k^{2}E(k)\,\tilde{G}^{\sqrt{\theta}}(k)^{2}$
$\displaystyle=$
$\displaystyle\dfrac{1}{2}\int_{\mathbb{R}}\,\,\,\mathrm{d}k\,E(k)\int_{\ell^{2}k^{2}}^{\infty}\mathrm{d}s\,\tilde{G}(\sqrt{s})^{2}$
(15)

where

$E(k)=\dfrac{1}{2}\int_{||\bm{k}||=k}\mathrm{d}\bm{k}\,\,\tilde{u}_{i}(\bm{k})\,\tilde{u}_{i}(-\bm{k})$
(16)

is the _omni-directional_ kinetic energy spectrum related to the velocity
Fourier modes $\tilde{u}_{i}(\bm{k})$ and its integration over
$k\in[0,\infty]$ provides the mean energy. Back to the RHS of (15), the kernel
dependency on the filtering-scale has been absorbed into its argument i.e.
$\tilde{G}^{\sqrt{\theta}}(k)=\tilde{G}(k\sqrt{\theta})$. This effective
rescaling holds for the standard filter functions, like the ones listed in
table 13.2 of [14] but in principle this would remove from our discussion
_exotic_ and substantially _intractable_ filter kernel types for which the
above property does not hold. Once again, the inversion of the integrals in
the derivation of (15) is a consequence of the assumed convergence of the
fields involved.

It is worth to highlight an observation [15, 16] about which if one sets
$\ell=0$ in (15), the second integral can be factorised from the first one; in
this way the RHS of (15) collapses to the product of $C_{G}\,\langle
E_{K}\rangle$ where $C_{G}$ is a number depending on the filter choice. In
such a way, we have generalised eq. (14) to a generic class of filters without
introducing the SGS stresses decomposition of (7) and without requiring the
associated positive definiteness of its trace. Note that a similar
factorization can be derived for the $nth$-derivative of the velocity gradient
[16].

However the generalisation originating from the previous observation does not
hold when (15) is evaluated at $\ell\neq 0$. In order to show this, we
consider the sharp filter whose Fourier transform is
$G^{\ell}(k)=H(\pi-k\,\ell)$ for which (15) becomes
$\int_{\mathbb{R}}\mathrm{d}k\,E(k)(\pi^{2}-\ell^{2}k^{2})/2$ which does not
equal the resolved-scale kinetic energy $\langle E^{\ell}_{K}\rangle$.

As a conclusion of this section, it is relevant to underline that
approximating the integral of eq. (14) at a given scale $\ell$, meaning
$\langle
E_{K}\rangle\approx\frac{1}{2}\ell^{2}||\overline{\bm{\omega}}^{\ell}||^{2}$,
resembles in solid mechanics the expression of the kinetic energy in terms of
moment of inertia of a ring of radius $\ell$ i.e.
$E_{K}=\frac{1}{2}\ell^{2}||\hat{\bm{\omega}}||^{2}$ where $\hat{\bm{\omega}}$
corresponds, in this case, to the angular velocity.

##  III Addressing the non-stationarity

Having obtained an exact expansion of the mean kinetic energy allows for the
dynamical characterization of non-stationary flows $\partial_{t}\langle
E_{K}(t)\rangle\gtrless 0$, respectively describing the flow energy increase
and the turbulence decay where the time dependence of the energy has been
purposefully reinstated. For this reason, applying a time derivative to eq.
(9) one can obtain an expression involving VG quantities. For simplicity we
focus on the homogeneous turbulence case, where the calculation of the time
derivative in (14) reads:

$\displaystyle\partial_{t}\langle
E_{K}(t)\rangle=\frac{1}{2}\int_{0}^{\infty}\mathrm{d}\theta\,\partial_{t}\,\langle\,||\overline{\bm{\omega}}^{\sqrt{\theta}}||^{2}\,\rangle$
(17)
$\displaystyle=\int_{0}^{\infty}\mathrm{d}\theta\left(-\dfrac{4}{3}\langle\overline{S}_{ij}^{\sqrt{\theta}}\overline{S}_{jk}^{\sqrt{\theta}}\overline{S}_{ki}^{\sqrt{\theta}}\rangle-\nu\langle\overline{u}^{\sqrt{\theta}}_{i,jk}\,\overline{u}^{\sqrt{\theta}}_{i,jk}\rangle-\langle\overline{S}^{\sqrt{\theta}}_{ij,\,kk}\tau_{ij}^{\sqrt{\theta}}\rangle+\langle
F^{\sqrt{\theta}}\rangle\right)$ (18)

where
$F^{\sqrt{\theta}}=\overline{u}^{\sqrt{\theta}}_{i,j}\overline{f}^{\sqrt{\theta}}_{i,j}$
depends on the external forcing and it is set to zero in the case of
turbulence decay. As concerns the other terms, the one containing the
contraction of three strain-rate tensors is linked to a contribution of
_energy transfer_ across the scales [4]; being also positive $\forall\,\ell$,
at least in the stationary case, we expect that it maintains the same sign in
the present non-stationary case as well. It is also clear that the term
proportional to the viscosity is negative. Especially in the context of
turbulence decay, the third term is formed by the contraction between the
strain-rate tensor Laplacian and the SGS stresses; being a priori sign
indefinite, it may play a relevant role in the above energy evolution
equation. Such a term also shows an intrinsic _multi-scale_ nature from the
presence of the SGS stress tensor which, by virtue of (7), depends on the
scales smaller than the integration scale $\sqrt{\theta}$. Applying (18) to
data from DNS, will asses the importance of such terms in unsteady
configurations.

Independently from the quantification of the terms, in case of stationarity
i.e. $\partial_{t}\langle E_{k}^{\ell}\rangle$, the above expression indicates
that the stationarity can be reached in two ways: a scale-by-scale
stationarity where the contribution of each integrated scale is identically
zero or the _scale-time ergodicity_ [17] where the RHS of (17) is zero because
of a compensation of the non-stationarity of the corresponding scales
contributions.

##  IV An exact identity for the kinetic energy spectrum

From the knowledge of the resolved-scale kinetic energy it is possible to
obtain the kinetic energy spectrum. Although very clear in the presence of the
aforementioned sharp filter, the following formula

$E_{K}(k)=-\dfrac{\pi}{k^{2}}\dfrac{\partial}{\partial k}\langle
E^{\ell}_{K}\rangle$ (19)

holds even for the Gaussian filter, where $E_{K}(k)$ is the already defined
kinetic energy spectrum. In order to get an expression in terms of velocity
gradient, we can plug eq. (14) into the previous equation obtaining:

$E_{K}(k)=\pi^{2}k^{-3}\,\langle||\bm{\overline{\omega}}^{\pi/k}||^{2}\rangle$
(20)

which is an exact expression for the kinetic energy flux where the wavenumber
dependency is both kept by the power of the wavenumber and by the filter scale
that the vorticity is coarse-grained at. Note that the previous expression is
an example of Infrared-Ultraviolet mixing (_UV-IR_) where the determination of
an observable at at the wavenumber $k$ corresponds to the another another one
at $1/k$. On top of that, the previous equation makes manifest that the energy
spectrum can be determined just by real space quantities like the resolved-
scale enstrophy. Moreover it is straightforward to verify that the integration
over all the wavenumber equals the mean kinetic energy $\langle E_{K}\rangle$.
The above expression must be compared with the common definition of the energy
spectrum of (16). Thus, we can consider the abovementioned data to leverage
this result. In homogeneous turbulence theory, when the study of energy
spectra is of interest, it is customary to discuss the corresponding
Kolmogorov’s scaling [18, 19, 13] which is
$E^{KOL}(k)=C\,\varepsilon^{2/3}\,k^{-5/3}$ where $C$ is the Kolmogorov
constant, according to the employed dataset, its value is approximately $1.6$,
see [9] while $\varepsilon=2\nu\langle\,||\bm{\omega}||^{2}\,\rangle$ is the
energy dissipation rate.

Figure 3: Omni-directional kinetic energy spectrum related to one
instantaneous configuration and calculated respectively via two different
formulations. The dotted line refers to the Kolmogorov power-law spectral
scaling.

Fig. 3 shows the comparison between eqs (16), (20) and the Kolmogorov scaling
in the evaluation of the kinetic energy spectrum. In order to discriminate the
former from the one of (20), we indicate as $E^{U}(kl)$ where the superscript
stands for _usual_. As regards $E(k)$ from (20), the first feature we notice
is the good agreement with Kolmogorov scaling 555Note that by [21],
$||\overline{\bm{\omega}}^{\ell}||\sim\ell^{-2/3}$ hence eq. (20) scales as
$\ell^{5/3}\sim k^{-5/3}$ like Kolmogorov formula. in the wavenumber band
$1.5\lesssim k/k_{f}\lesssim 30$. Moreover we also appreciate that it
reproduces $E^{U}_{K}(k)$ even in the forcing range; the overall smoothed
profile is an _inheritance_ of the Gaussian filter acting to the averaged
resolved-scale enstrophy. For those scales where $E^{KOL}(k)=E(k)$, one can
tentatively write an expression for the Kolmogorov constant:

$C=\pi^{2}\,\dfrac{\langle||\bm{\overline{\omega}}^{\pi/k}||^{2}\rangle}{\varepsilon^{2/3}\,k^{4/3}}$
(21)

that requires a further investigation.

##  V Conclusions

In conclusion, from a more general identity relating the total kinetic energy
and the coarse-grained velocity gradient norm based on the Gaussian filter
kernel, two (sub)identities are derived for homogeneous turbulence. One
expresses the mean energy balance in terms of the intensity of the strain
motion, outside the range of scale dominated by the forcing, the purely
extensional and contractile deformations provides approximately $54\%$ and
$41\%$ of the resolved scale energy. The remainder is linked to the
intermediate eigenvalues which can be either related to contractions and
extensions. Using kinematic constraints due to incompressibility, one can
obtain a cognate identity which is reshaped in terms of the enstrophy. This
new relation shows that the sum of enstrophy across the scales determines the
total kinetic energy. Both the identities can be used to characterise the mean
velocity gradient dynamics of unsteady flows where the interpretable terms
governing the unsteadiness have been identified. In consequence, one can
obtain an exact identity for the kinetic energy spectrum providing a deeper
fundamental understanding on the energy spectrum and an expansion in real
space interpretable observables together with a guidance for the analytical
determination of the Kolmogorov constant.

As future directions, the filtered enstrophy based identity can be potentially
used in the study of vortex flow where the dynamics of isolated vortices is of
interest. On top of that it can used to improve the generation of synthetic
data in data-driven approaches. The analysis dynamical equation regarding
unsteady flows can be deeply expanded and quantified to unveil the governing
physical mechanisms that lead the decay or the increase of energy.

##  VI Acknowledgments and fundings

The author thanks Prof. Luca Biferale and Dr. Massimo Cencini for the
financial support provided by both European Research Council (ERC) under the
European Union’s Horizon 2020 research and innovation programme (grant
agreement No 882340) and FieldTurb experiment from the Istituto Nazionale di
Fisica Nucleare (INFN). This work used the ARCHER2 UK National Supercomputing
Service (www.archer2.ac.uk) with resources provided by the UK Turbulence
Consortium (EPSRC grants EP/R029326/1). The author thanks Dr. Moritz Linkmann
for the interesting discussions, encouragement and the support. The author
wants to thank Prof. Michael Wilczek and Lukas Bentkamp for their precise
observations.

## References

  * Germano [1992] M. Germano, Turbulence: the filtering approach, Journal of Fluid Mechanics 238, 325 (1992).
  * Vreman  _et al._ [1994] B. Vreman, B. Geurts, and H. Kuerten, Realizability conditions for the turbulent stress tensor in large-eddy simulation, Journal of Fluid Mechanics 278, 351 (1994).
  * Aluie [2017] H. Aluie, Coarse-grained incompressible magnetohydrodynamics: analyzing the turbulent cascades, New Journal of Physics 19, 025008 (2017).
  * Johnson [2020] P. L. Johnson, Energy transfer from large to small scales in turbulence by multiscale nonlinear strain and vorticity interactions, Phys. Rev. Lett. 124, 104501 (2020).
  * Johnson [2021] P. L. Johnson, On the role of vorticity stretching and strain self-amplification in the turbulence energy cascade, Journal of Fluid Mechanics 922, A3 (2021).
  * Note [1] Note that the methodology of [4] requires homogeneity or periodic boundaries in order to set to zero the boundary terms and not for the derivation of the decomposition (7) itself.
  * Note [2] The length contraction in general relativity would already pose conceptual difficulties to the present filtering approach.
  * Note [3] However, in the applications, the interest obviously focuses on systems characterised by a maximum allowed scale for which the integration upper bound of (7) does not diverge ensuring a rigorous exchange of integrals to derive (8).
  * Biferale  _et al._ [2024] L. Biferale, F. Bonaccorso, M. Linkmann, and D. Capocci, Turb-hel: an open-access database of helically forced homogeneous and isotropic turbulence, arXiv preprint arXiv:2404.07653 (2024).
  * Betchov [1956] R. Betchov, An inequality concerning the production of vorticity in isotropic turbulence, Journal of Fluid Mechanics 1, 497 (1956).
  * Tsinober [2001] A. Tsinober,  _An informal introduction to turbulence_ , Vol. 63 (Springer Science & Business Media, 2001).
  * Note [4] This quantification has been repeated on datasets with different $Re$ like HL1 from [9]. We obtain qualitative the same saturation profile and the same asymptotes as those displayed in fig. 2.
  * Frisch [1995] U. Frisch,  _Turbulence: the legacy of AN Kolmogorov_ (Cambridge University Press, 1995).
  * Pope [2000] S. B. Pope,  _Turbulent Flows_ (Cambridge University Press, 2000).
  * [15] L. Bentkamp, Private communication.
  * [16] M. Wilczek, Private communication.
  * [17] M. Linkmann, Private communication.
  * Kolmogorov [1941] A. N. Kolmogorov, Dissipation of Energy in Locally Isotropic Turbulence, Akademiia Nauk SSSR Doklady 32, 16 (1941).
  * Kolmogorov [1991] A. N. Kolmogorov, The Local Structure of Turbulence in Incompressible Viscous Fluid for Very Large Reynolds Numbers, Proceedings of the Royal Society of London Series A 434, 9 (1991).
  * Note [5] Note that by [21], $||\overline{\bm{\omega}}^{\ell}||\sim\ell^{-2/3}$ hence eq. (20) scales as $\ell^{5/3}\sim k^{-5/3}$ like Kolmogorov formula.
  * Eyink [1995] G. L. Eyink, Local energy flux and the refined similarity hypothesis, Journal of Statistical Physics 78, 335 (1995).