# Switching of Magnetization in Quantum Antiferromagnets with Time-Dependent
Control Fields

Asliddin Khudoyberdiev<EMAIL_ADDRESS>Condensed Matter
Theory, TU Dortmund University, Otto-Hahn-Straße 4, 44221 Dortmund, Germany
Götz S. Uhrig<EMAIL_ADDRESS>Condensed Matter Theory, TU Dortmund
University, Otto-Hahn-Straße 4, 44221 Dortmund, Germany

###### Abstract

Ultrafast manipulation of magnetic states is one of the necessities in modern
data storage technology. Quantum antiferromagnets are promising candidates in
this respect. The orientation of the order parameter, the sublattice
magnetization, can be used to encode “0” and “1” of a bit. Then, the switching
of magnetization and the full control of its reorientation are crucial for
writing data. We show that the magnetization can be switched efficiently by an
external short THz pulses with relatively low amplitude. The coupling to the
spin degrees can be direct via the magnetic field or indirect via the electric
field inducing spin-polarized charge currents. Our description is based on
time-dependent Schwinger boson mean-field theory which includes the intrinsic
dephasing mechanisms beyond a macrospin description. The findings help to
introduce the use of antiferromagnets in data storage technology.

##  I Introduction

Miniaturization and acceleration are two key goals in modern data technology.
Since long, ferromagnetism has been exploited for data storage and
manipulation [1]. Recently, antiferromagnetism has attracted much interest for
two reasons essentially [2, 3, 4]. First, the absence of stray fields for
dense packing of the magnetic domains [5] representing single bits because
they do not perturb one another beyond the length scale set by the domain wall
thickness. Second, their generic energy scales are larger by about a factor
1000 so that ultrafast manipulation comes into reach [6, 7].

If one wants to use the orientation of the sublattice magnetization to encode
information it is crucial to be able to control it. On the one hand,
completely spin isotropic systems allow for ideal, effortless rotation of the
orientation. On the other hand, however, they possess no robustness against
any perturbation, for instance due to noise fields, so that the information
encoded in the orientation is likely to be lost. Thus, one focuses on systems
with a certain degree of spin anisotropy. This, in turn, sets a certain
activation energy which needs to be overcome to rotate the magnetization. How
this can be done is topic of a plethora of current experimental and
theoretical studies, see Ref. 8 for a review.

In the present theoretical study, we consider the simple model of an easy-axis
anisotropic spin-$1/2$ Heisenberg model on the square or cubic lattice, i.e.,
in two or three dimensions. For simplicity, we incorporate the external
control field as an effective, uniform magnetic field coupling to the
localized spins via a Zeeman term. It can be thought of being the magnetic
field component of a THz pulse as in Refs. [7, 9, 10] or being generated by
the electric field inducing currents which in turn exert torques on the spins,
spin-transfer torques [11] or spin-orbit torques [12, 13, 14]. We stress that
the time-dependence of the control field is crucial because it will allow us
to adjust the pulses to be resonant with the spin anisotropy gap. Being on
resonance is a key aspect in order to manipulate quantum states efficiently
and to overcome activation energies [15, 16].

Describing the time evolution of the reorientation of magnetization reliably
is far from an easy task. Brute force numerics quickly comes to its limits
because only about 20 to 30 quantum spins can be treated which is not
sufficient for the treatment of long-range order. In many studies, the
sublattice magnetization is described by a macrospin so that
antiferromagnetism requires two macrospins, see, e.g., Refs. [17, 18, 19, 6,
20, 21]. Any dephasing and/or relaxation is incorporated in the standard form
of Landau-Lifshitz-Gilbert (LLG) equation with a phenomenological Gilbert
damping $\alpha$. Then, one can deal with the precessional time evolution
described by six variables (twice three components) or suitable further
approximate reductions.

An extension consists in the use of LLG equations for each spin [22, 23].
Furthermore, finite temperature can be accounted for by the Landau-Lifshitz-
Bloch equation [22, 24] which comprises a stochastic Langevin term to account
for thermal fluctuations. Again, these are essentially classical equations,
but with a number of variables scaling linearly with system size. Indeed,
classical spin variable are well suited to capture the dynamics of quantum
spins as well, at least in disordered states [25] which can be justified by
the truncated Wigner approximation [26, 27].

Since the candidate antiferromagnets for data storage are used deep in their
long-range ordered phases, the application of classical equations of motion
can be questioned. Quantum fluctuations also play a role. The usual way to
capture the leading quantum and low-temperature thermal fluctuations is spin
wave theory about the long-range ordered ground state. The two most common
bosonic representations of spins are the Holstein-Primakoff representation
[28] and the Dyson-Maleev representation [29, 30] which both describe small
deviation around a specific ordered state. Thus they are not suitable to
describe large deviations as they necessarily occur in switching processes
where the orientation of the magnetization is rotated by 90∘ or 180∘ so that
the state after the switching is as far away from the initial state as
possible. Here a third representation comes in handy which is due to Schwinger
[31] allowing to describe disordered spin states or orders with arbitrary
orientation [32, 33]. Very recently, it was shown that a time-dependent
version of Schwinger boson mean-field theory can capture switching processes
[34]. This approach has the advantage that the leading quantum and thermal
fluctuations are captured. Not only the modes at zero momentum appear in the
formalism, but the modes at all momenta so that the important intrinsic
dephasing mechanism is included in the theoretical treatment without adhoc
assumptions and phenomenological parameters such as the Gilbert damping.

In the preceding study [34] the switching was investigated as it is induced by
a constant, homogeneous magnetic field $h=g\mu_{\text{B}}B$ applied to a two
dimensional square lattice. This switching field needs to be sufficiently
strong to overcome the anisotropy. At zero temperature, we found that the
threshold field is fairly precisely given by the spin gap
$\Delta=h_{\mathrm{thresh}}$. But at temperatures close to the Néel
temperature, the threshold field becomes small.

Experimentally the samples are dealt with at temperatures well below the
ordering Néel temperature $T_{\text{N}}$. Heating up close to $T_{\text{N}}$
necessitates heating pulses which require time and will slow down any
switching. Thus, it is attractive to look for other means to lower the minimum
switching fields. Generally, overcoming some energy barrier in quantum
mechanics can be achieved by exciting in resonance, see e.g. Ref. [16]. Hence,
it is interesting to extend the previous treatment to time-dependent control
fields which can be tuned in resonance to the energy gap of the magnetic
system to be switched. With realizations in mind, such control fields can be
generated by THz pulses because their frequencies lie well in the range of
magnetic gaps.

In this work, we show that the switching of antiferromagnetic sublattice
magnetization can be achieved by applying experimentally accessible THz pulses
with moderate amplitudes even at zero temperature. The optimum angular carrier
frequency of the pulses lies slightly below the spin gap $\Delta/\hbar$. The
theoretical approach and equations stay the same as in Ref. [34] except that
the control magnetic field acquires a time dependence. Furthermore, we extend
previous results from two to three dimensions, i.e., we complement results for
the square lattice by results for the cubic lattice. This important extension
becomes possible by methodological progress: generalizing the equations from
momentum space to a sort of frequency space weighted by the density-of-states.

This paper is organized as follows. In Section II we briefly recall the main
equations. The magnetization dynamics of antiferromagnets under short pulses
will be analyzed in Sect. III. Additionally, the set of optimal parameters
will be established. Sect. IV summarizes our results.

##  II Model and theoretical tools

The approach focuses on the two dimensional square and the three dimensional
cubic lattices with spin $S=1/2$. The first step in the theoretical treatment
is to determine the initial conditions, i.e., the thermodynamic equilibrium.
This has been done in literature many times [32, 33, 35, 34]. We recall the
equilibrium below to the extent that it fixes the initial values as in Ref.
[34].

###  II.1 Model and Schwinger boson mean-field theory

We consider a Heisenberg model with easy-axis anisotropy defined by the
anisotropy parameter $\chi=J_{xy}/J_{z}\in[0,1]$ so that the Hamiltonian takes
the form

$\mathcal{H}_{0}=J\sum_{\langle
i,j\rangle}\left\\{\frac{\chi}{2}\left(S^{+}_{i}S^{-}_{j}+S^{-}_{i}S^{+}_{j}\right)+S^{z}_{i}S^{z}_{j}\right\\}$
(1)

in terms of the spin ladder operators. The switching is induced by the time
dependent magnetic field included in the formalism by adding the term

$\mathcal{H}(t)=-\mathbf{h}(t)\cdot\sum_{i}\mathbf{S}_{i}.$ (2)

In the simulations, we choose the magnetic field to point along the $y$ axis
because we rotate sublattice B about the spin $y$ axis by 180∘ to restore
translational invariance of the original lattice even in the long-range
ordered state, see below. In principle, the $x$ and the $y$ axis are
equivalent in this problem, but due to the sublattice rotation the application
of the control field along the $y$ axis is convenient because this preserves
the translational invariance of the resulting Hamiltonian.

Figure 1: Illustration of the switching after a THz magnetic field pulse. The
delay $\Delta t=t-t_{0}$ is of the order of picoseconds. The magnetization of
the $A$ sublattice is plotted as redish arrow while the magnetization of the
$B$ sublattice is plotted as blue arrow.

Figure 1 sketches the process of switching in our case by the help of
ultrashort magnetic field pulses. One of our goals is to obtain switching with
$\Delta t$ being of the order of picoseconds. If this succeeds one can claim
benefits of antiferromagnets over ferromagnets.

We begin with a representation of quantum spin operators in terms of Schwinger
bosons [31, 36, 37]. Two flavors of bosons $a_{i}$ and $b_{i}$ can be present
at one site of the lattice and express the spin operators by

$\displaystyle S_{i}^{+}$ $\displaystyle=a_{i}^{\dagger}b_{i}$ (3a)
$\displaystyle S_{i}^{-}$ $\displaystyle=b_{i}^{\dagger}a_{i}$ (3b)
$\displaystyle S_{i}^{z}$
$\displaystyle=\frac{1}{2}\left(a_{i}^{\dagger}{a}_{i}-b_{i}^{\dagger}{b}_{i}\right).$
(3c)

The bosonic Hilbert space exceeds the physical relevant subspace so that a
constraint needs to be included which reads

$2S=a_{i}^{\dagger}{a}_{i}+b_{i}^{\dagger}{b}_{i}.$ (4)

The sublattice magnetization reads

$m=\langle S_{i}^{z}\rangle=\frac{1}{2}\left(\langle
a_{i}^{\dagger}{a}_{i}\rangle-\langle b_{i}^{\dagger}{b}_{i}\rangle\right),$
(5)

so that the mean occupation of the bosons determines the magnetization. For
further analysis, the bipartite lattice is divided into two interpenetrating
sublattices, A and B as is common for unfrustrated antiferromagnets so that
the nearest neighbors of the A spins are all on sublattice B and vice versa.
Additionally, we introduce bond operators $A_{ij}=a_{i}b_{j}-b_{i}a_{j}$ and
$B_{ij}=a_{i}b_{j}+b_{i}a_{j}$. Then, it is convenient to apply a canonical
transformation to the Hamiltonian in which the spins on one sublattice, say B,
are rotated by $180^{\circ}$ about the $y$ axis with the effect that
$a_{j}\rightarrow-b_{j}$ and $b_{j}\rightarrow a_{j}$. Consequently, the bond
operators transform into the simpler, more symmetric form
$A_{ij}=a_{i}a_{j}+b_{i}b_{j}$ and $B_{ij}=a_{i}a_{j}-b_{i}b_{j}$. Remarkably,
the Hamiltonian keeps translational invariance in this representation even in
the antiferromagnetically ordered phase. Finally, we simplify the resulting
quartic Hamiltonian by applying a mean-field theory using the expectation
values [33, 34]

$\displaystyle A$ $\displaystyle\coloneqq\langle A_{ij}\rangle$ (6a)
$\displaystyle B$ $\displaystyle\coloneqq\langle B_{ij}\rangle.$ (6b)

The total mean-field Hamiltonian $\mathcal{H}=\mathcal{H}_{0}+\mathcal{H}(t)$
from (1) and (2) takes the following form after a standard Fourier
transformation

$\displaystyle\mathcal{H}$
$\displaystyle=E_{0}-\frac{z}{8}\sum_{\mathbf{k}}\gamma_{\mathbf{k}}\Big{(}C_{-}a_{\mathbf{k}}^{\dagger}a_{\mathbf{-k}}^{\dagger}+C_{+}b_{\mathbf{k}}^{\dagger}b_{\mathbf{-k}}^{\dagger}$
$\displaystyle\qquad+C_{-}^{*}a_{\mathbf{k}}a_{\mathbf{-k}}+C_{+}^{*}b_{\mathbf{k}}b_{\mathbf{-k}}\Big{)}+\lambda\sum_{\mathbf{k}}\left(a_{\mathbf{k}}^{\dagger}{a}_{\mathbf{k}}+b_{\mathbf{k}}^{\dagger}{b}_{\mathbf{k}}\right)$
$\displaystyle\qquad-\frac{h_{y}(t)}{2i}\sum_{\mathbf{k}}\left(a_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}-b_{\mathbf{k}}^{\dagger}a_{\mathbf{k}}\right),$
(7)

where the superscript $*$ stands for the complex conjugate. The momentum
dependence enters through the parameter

$\gamma_{\mathbf{k}}\coloneqq\frac{1}{d}\sum_{i=1}^{d}\cos(k_{i})$ (8)

where we set the lattice constant to unity and $k_{i}$ is the component of the
wave vector in $i$ direction. The coordination number, i.e., the number of
nearest neighbors of each site, is denoted by $z$ given by twice the dimension
$2d$ for hypercubic lattices and we use

$C_{\pm}\coloneqq A(1+\chi)\mp B(1-\chi).$ (9)

The energy constant reads

$E_{0}\coloneqq\frac{zN}{4}\left\\{|A|^{2}+|B|^{2}+\chi(AB^{*}+BA^{*})+2S^{2}\right\\}.$
(10)

The term with the site independent Lagrange multiplier $\lambda$ has been
added by hand in order to ensure that the constraint (4) is fulfilled, at
least on average over the whole lattice.

In the case of equilibrium without time-dependent Zeeman term, the gapped
dispersion relations can be obtained from standard Bogoliubov transformations
to diagonal $\alpha$\- and $\beta$-bosons with dispersion
$\omega^{-}_{\mathbf{k}}$ and $\omega^{+}_{\mathbf{k}}$, respectively, which
are given by

$\omega^{\pm}_{\mathbf{k}}=\sqrt{\lambda^{2}-\left({z|C_{\pm}|\gamma_{\mathbf{k}}}/{4}\right)^{2}}.$
(11)

Here, $C_{\pm}$ is still given by Eq. (9), but the expectation values $A$ and
$B$ and thus $C_{\pm}$ can be chosen real. In the above dispersion relations,
both types of bosons acquire a energy gap as

$\Delta^{\pm}\coloneqq\omega^{\pm}_{\textbf{k}=0}.$ (12)

If we iterate the above equations, see Appendix A, to reach convergence in the
expectation values starting with more $a$ bosons than $b$ bosons we obtain a
positive sublattice magnetization and the $b$ bosons have the larger energy
gap $\Delta^{+}$. This is the initial scenario from which we start. Then, the
physical spin gap is given by the difference

$\Delta\coloneqq\Delta^{+}-\Delta^{-}$ (13)

because the boson number in the physical subspace can only change by two, not
by one. Of course, we could also start from majority $b$-bosons; but this
would only make a difference in the sign of the sublattice magnetization.

For a static control field we found that the minimum field
$h_{\mathrm{thresh}}$ is almost equal to the spin gap $\Delta$. In particular
for $\chi\rightarrow 1$, $h_{\mathrm{thresh}}$ and $\Delta$ both vanish
proportional to $\sqrt{1-\chi^{2}}$ [34]. In Schwinger boson representation
the switching can be understood as converting the majority $a$ bosons to $b$
bosons, i.e., to make the $b$ bosons the majority bosons. In parallel, the
size of the spin gap determines the robustness of the system against any kind
of perturbation. In this work, we consider oscillating fields of short pulses
with angular frequencies close to $\Delta/\hbar$. Indeed, such coupling in
resonance facilitates the spin reversal.

For the time dependent mean-field theory we need the time dependence of
expectation values which can be computed by Heisenberg’s equation of motion in
momentum space. They can be efficiently expressed as functions of the
parameter $\gamma$ because the momentum dependence only enters via this
parameter, see Sect. II.2. The resulting differential equations (20) are given
in in Appendix A. Eventually, the dynamics of sublattice magnetization follows
from

$m(t)=\frac{1}{2N}\sum_{\mathbf{k}}\big{(}\langle
a_{\mathbf{k}}^{\dagger}{a}_{\mathbf{k}}\rangle-\langle
b_{\mathbf{k}}^{\dagger}{b}_{\mathbf{k}}\rangle\big{)}.$ (14)

###  II.2 Equations of motion based on the density-of-states of
$\gamma_{\mathbf{k}}$

In order to calculate the quantities in the systems under study and to solve
the differential equations (20) fulfilling Eqs. (19), one needs to carry out
sums of the type $\sum_{\mathbf{k}}F_{\mathbf{k}}$ where the wave vectors
$\mathbf{k}$ run over the Brillouin zone. In the thermodynamic limit, these
sums become integrals. These integrals can be simplified further because the
dependence on the wave vectors only enters through $\gamma_{\mathbf{k}}$
defined in (8), see for instance Eq. (7), Eq. (19), or Eqs. (20). Hence, one
can replace the integrals over all wave vectors by a single integral over
$\gamma\in[-1,1]$. But it is important to weight the contributions according
to how many wave vectors contribute to a certain value of $\gamma$ [33, 36].
This is achieved as usual by a density-of-states $\rho_{d}$ in dimension $d$
as function of $\gamma$. We define it as

$\rho_{d}(\gamma)\coloneqq\frac{1}{(2\pi)^{d}}\int_{\text{BZ}}\delta(\gamma-\gamma_{\mathbf{k}})dk^{d}.$
(15)

With the help of $\rho(\gamma)$ sums over the two or three dimensional
Brillouin zone can be reduced to one dimensional integrals in the
thermodynamic limit

$\lim_{N\to\infty}\frac{1}{N}\sum_{\mathbf{k}}F(\gamma_{\mathbf{k}})=\int_{-1}^{1}d\gamma\rho(\gamma)F(\gamma).$
(16)

In the numerical calculations, we discretize these integration again. But the
calculation remains very efficient because only a one dimensional integration
needs to be dealt with. In practice, we use the summed Newton-Cotes midpoint
interval rule. We checked the results obtained via (16) on a sample basis by
sums over the first Brillouin zone of the square lattice for linear system
size $L=500$ and obtained the same results within numerical accuracy.

##  III Effect of short pulses

We aim at switching the antiferromagnetic magnetization by pulses of low
magnetic fields in the THz range. Thus, we want to model and to optimize
realistic pulses so that our results can guide future experiments. Short
magnetic pulses are technologically well possible in contrast to continuous
excitation with continuous radiation. Thus, we focus on pulses of the form

$h_{y}(t)=h_{0}\cos(\alpha\Delta(t-3\tau)+\phi_{0})\cdot
e^{-\frac{(t-3\tau)^{2}}{2\tau^{2}}}$ (17)

where $h_{0}$ is the maximum amplitude, $\Delta$ is the spin gap (13) which
essentially sets the excitation frequency (setting $\hbar=1$) although it is
slightly modified by the prefactor $0<\alpha\leq 1$; a phase shift $\phi_{0}$
is also accounted for. It helps to tune the field induced torque in phase with
the spin precessions as discussed in Ref. [7]. The envelope is chosen Gaussian
to treat a generic pulse of essentially finite duration $\tau$. We start our
simulations at $t=0$; in order to capture the significant part of the pulse it
is shifted by $3\tau$. We strive for the optimum values of the parameters in
(17) for switching. In Appendix B we provide results also for two other pulse
shapes.

###  III.1 Dynamics of the magnetization for the square lattice at zero
temperature

Fig. 2 demonstrates how the magnetization evolves in time as a result of the
pulse in (17). The magnetic amplitudes $h_{0}=0.158\,J$ and $h_{0}=0.107\,J$
are slightly higher than the threshold fields for $\chi=0.9$ given by
$h_{0}^{\mathrm{thresh}}=0.1578\,J$ and for $\chi=0.98$ given by
$h_{0}^{\mathrm{thresh}}=0.1065\,J$. For comparison, the threshold fields
required in the static case are $h_{\mathrm{thresh}}=0.863\,J$ and
$h_{\mathrm{thresh}}=0.364\,J$, correspondingly. One sees that the control at
resonance frequencies truly pays. The blue curves depicting $h_{y}(t)$ clearly
show quicker oscillations for a larger spin gap comparing Fig. 2(a) and in
Fig. 2(b). Another notable point is that the switching takes a bit longer for
lower anisotropy with lower spin gap. In parallel, the oscillations of the
transversal component $\langle S^{x}\rangle$ is also slower closer to isotropy
because the spin gap is smaller and thus also the chosen pulse frequency.
Based on these results, we emphasize that the switching is obtained with
relatively low maximum amplitudes of the THz pulse. The other parameters are
chosen as $\alpha=0.85$, $\tau=10$ $J^{-1}$, and $\phi_{0}=\pi/3$ because our
analyses indicate that these values are best suited for switching, for further
discussion see below.

Figure 2: Temporal evolution of the sublattice magnetization $m(t)$ (red
curve) for $\chi=0.9$, $h_{0}=0.158\,J$ in panel (a) and $\chi=0.98$,
$h_{0}=0.107\,J$ in panel (b). The evolution of the pulse $h_{y}(t)$ is
included as blue curves. The dynamics of the transversal spin component
$\langle S^{x}\rangle$ is shown in light blue. The spin gap takes the values
$\Delta(\chi=0.9)=0.837\,J$ and $\Delta(\chi=0.98)=0.357\,J$. Other parameters
are chosen as $\alpha=0.85$, $\phi_{0}=\pi/3$, and $\tau=10\,J^{-1}$.

To show the dynamics of the magnetization for various amplitudes, we plot it
in Fig. 3. This plot illustrates nicely the threshold behavior which occurs
also for the time-dependent pulses. Only for maximal amplitudes above a
certain value the magnetization switches sign; this value is about
$h_{0}^{\mathrm{thresh}}\approx 0.1065$ $J$. For lower amplitudes no sign
reversal takes place.

Figure 3: Dynamics of the magnetization $m(t)$ at $\chi=0.98$ for
$h_{0}\in[0.100\,J,\,0.112\,J]$. For amplitudes above
$h_{0}^{\mathrm{thresh}}=0.1065\,J$ switching takes place. All other
parameters are the same as in Fig. 2. Figure 4: Threshold value
$h_{0}^{\mathrm{thresh}}$ as function of the renormalization $\alpha$ of the
excitation frequency in (17) for two anisotropies. The minima are found by
fitting parabolae. Clearly, the optimum excitation is obtained for
$\alpha\approx 0.85$. The other parameters of the pulse are $\tau=10\,J^{-1}$,
$\phi_{0}=\pi/3$.

In Fig. 4 one can see that the excitation frequency plays an important role in
switching. Generically, it helps if the pulse is in resonance with the spin
gap $\Delta$. But the activation energy induced by the anisotropy acts like a
hill which needs to be overcome. Thus, in the course of switching the slope of
this hill diminishes and the required frequency diminishes as well. This has
been shown explicitly for a large spin [16]. To account for this effect, at
least partially, we introduced the frequency renormalization $\alpha$ in the
pulse (17). Indeed, the data depicted in Fig. 4 shows that the lowest
threshold field is not obtained for $\alpha=1$, but for values around $0.8$.
In order to determine the minimum we fitted the data by parabolae and obtained
for both anisotropies an optimum value of $\alpha\approx 0.85$.

Another parameter in the ansatz (17) is $\tau$. It appears as delay and as
pulse duration. We refrain from considering the delay separately because for a
Gaussian pulse one can be quite sure that it has essentially all its effect if
it is applied over a full-width interval of six times its standard deviation
which is given by $\tau$ in (17). Thus, delaying the pulse peak to $3\tau$
captures almost all the effect of the pulse.

Figure 5: Threshold fields (symbols) vs. pulse duration $\tau$ including the
indicated fits. The other parameters are $\alpha=0.85$ and $\phi_{0}=\pi/3$.
Dashed horizontal lines correspond to the threshold values for uniform
magnetic field with $h_{\mathrm{thresh}}=0.863\,J$ and
$h_{\mathrm{thresh}}=0.364\,J$ for $\chi=0.9$ and $\chi=0.98$, respectively.
The fits (solid lines) are power laws $h_{0}^{\mathrm{thresh}}=a\tau^{b}+c$
with parameters $a=1.007\,J^{b+1}$, $b=-1.158$, $c=0.093\,J$ for $\chi=0.9$
and $a=1.195\,J^{b+1}$, $b=-1.167$, $c=0.026\,J$ for $\chi=0.98$.

Still, the duration of the pulse matters so that we analyze the dependence on
$\tau$ in Fig. 5. Obviously, the threshold field is very high for small values
of $\tau$. For increasing $\tau$ the threshold field decreases as one would
naively expect because a longer lasting pulse has more effect or it can have
the same effect at lower amplitude. We fit the threshold fields by a power law
$a\tau^{b}+c$ and find the values indicated in the caption. The exponent is of
the order of $-1$, but the finite values $c>0$ indicates that even a very long
pulse requires some finite minimum amplitude to achieve switching. This is in
line with what we found for static fields [34]. The dashed horizontal lines
indicate the threshold fields required for static control fields. From the
data in Fig. 5 we choose $\tau=10/J$ as a reasonable tradeoff between a small
threshold amplitude and a not too long switching time of about $40/J$.

Figure 6: Threshold fields (symbols) vs. phase $\phi_{0}$ in (17) for
$\alpha=0.85$ and $\tau=10$ $J^{-1}$.

Last, but not least, we study the effect of the phase on switching. Fig. 6
shows the threshold amplitude of the pulse as function of the phase
$\phi_{0}$. The overall dependence is rather weak, i.e., the phase does not
have a large effect on the threshold field. Quite unexpectedly, however, there
is one jump occuring where the threshold field suddenly drops upon changing
the phase. This jump is related to a sudden shift in time at which the
switching actually takes place. Just below a certain phase $\phi_{0}$ the
switching occurs at a zero of the control field while just above that value
the switching occurs almost a control cycle later. For an illustration of this
subtlety we refer the reader to Fig. 12 in Appendix B. We have chosen the near
optimum value of $\phi_{0}=\pi/3$; it is applicable for any value of a weak
anisotropy parameter in our analysis. But generally, the influence of the
phase is rather weak.

###  III.2 Dynamics of the magnetization for the cubic lattice at zero
temperature

Since we solve the equations of motion not in reciprocal space, but in
frequency space, more precisely in $\gamma$ space, there is not major
difference in the numerical effort between treating the two dimensional and
the three dimensional cases. The only necessary modification is to choose the
appropriate density-of-states $\rho_{d}$. Explicit results can be taken from
literature [38]. As a first step, we compute and compare the spin gap at zero
temperature and the static threshold field in Fig. 7. As for the square
lattice, both quantities agree very well. This supports the view that the spin
gap is a measure for the activation energy to be overcome for switching in
case of anisotropies.

Figure 7: Threshold field $h_{\mathrm{thresh}}$ and spin gap $\Delta$ both at
zero temperature as function of the anisotropy $\chi$. They almost coincide as
was found previously in two dimensions as well [34].

Next, we proceed with the analysis of the effects of short Gaussian pulses
(17). First, we determine the optimum value of the frequency renormalization
$\alpha$ in Fig. 8. The results are very similar to the results in two
dimensions. A slight lowering of the spin gap frequency by $\alpha\approx
0.85$ appears to allow for the lowest amplitudes for switching.

Figure 8:  The threshold value $h_{0}^{\mathrm{thresh}}$ as function of the
renormalization $\alpha$ of the excitation frequency in (17) for two
anisotropies. The minima are found by fitting parabolae. The other parameters
of the pulse are $\tau=10\,J^{-1}$, $\phi_{0}=\pi/3$.

Also the dependence of the threshold amplitude $h^{\mathrm{thresh}}_{0}$ on
the pulse duration $\tau$ is qualitatively the same as for the square lattice.
The only quantitative difference is the value of the threshold amplitude
$h_{0}^{\mathrm{thresh}}$ which is higher for the cubic lattice in line with
the spin gap. This is a direct consequence of the higher coordination number
stabilizing the ordered phase and enhancing the effect of the exchange
anisotropy. We conclude that the theoretical findings will carry over to many
other comparable types of lattices.

Figure 9: Threshold fields (symbols) vs. pulse duration $\tau$ including the
indicated fits. The other parameters are $\alpha=0.85$ and $\phi_{0}=\pi/3$.
Dashed horizontal lines correspond to threshold values at uniform magnetic
field with $h_{\mathrm{thresh}}=1.281\,J$ and $h_{\mathrm{thresh}}=0.578\,J$
for $\chi=0.9$ and $\chi=0.98$, respectively. Fitting is done with the same
power law $h_{0}^{\mathrm{thresh}}=a\tau^{b}+c$ as in 2D. The fitting
parameters are $a=0.937\,J^{b+1}$, $b=-1.068$, $c=0.084\,J$ for $\chi=0.9$ and
$a=1.127\,J^{b+1}$, $b=-1.248$, $c=0.051\,J$ for $\chi=0.98$.  Figure 10:
Threshold fields (symbols) vs. phase $\phi_{0}$ in (17) for $\alpha=0.85$ and
$\tau=10$ $J^{-1}$.

Fig. 10 shows the threshold field dependence on the phase in three dimensions.
Again, we find the phase of the control field has only a weak effect on the
threshold value as for the square lattice. We again found a jump in the phase
dependence resulting from a shift in the switching time as shown and discussed
in Appendix B for the square lattice. Due to the small influence of the
precise phase value, we fixed the phase to the value $\phi_{0}=\pi/3$ for
cubic lattice although this does not represent the best possible value, cf.
Fig. 10 .

###  III.3 Dynamics at finite temperature

In the case of switching with static fields it was found that elevated
temperatures in the vicinity of the Néel temperature reduced the threshold
field $h_{\mathrm{thresh}}$ significantly [34]. Thus, it suggests itself to
study the effect of finite temperatures also for time-dependent short pulses.
This is achieved in Fig. 11 for the square lattice where the threshold
amplitudes are depicted as function of $\alpha$ for $\chi=0.98$ for four
different temperatures. Interestingly enough, the overall shape of this
dependence becomes flatter and flatter as the temperature is increased. This
implies that for control pulses far off resonance the increased temperature
facilitates switching as found before in the static case [34]. But in the
optimum region around $\alpha=0.85$ the changes are very small. This is an
interesting finding carrying the message to experiment that it pays off better
to choose the optimum resonance conditions, perhaps accounting for some
renormalization by 10 to 20%, than to work at elevated temperatures. If the
Néel temperature lies well above room temperature, all experimental studies
can well be conducted at room temperature.

Figure 11: The threshold value $h_{0}^{\mathrm{thresh}}$ as function of the
renormalization $\alpha$ of the excitation frequency in (17) as in Fig. 4 for
$\chi=0.98$ for various temperatures. Recall that the Néel temperature takes
the value $T_{\text{N}}=0.578\,J$ [34]. All other parameters are chosen as
before.

##  IV Conclusions

We studied the time-dependent control of the sublattice magnetization in two
and three dimensional anisotropic quantum antiferromagnets. The motivation is
to use such systems for data storage; the anisotropy being crucial to make the
long-range order robust. Our study uses spin wave theory based on the
representation of the spins by Schwinger bosons and the corresponding mean-
field theory. The representation by Schwinger bosons has the fundamental
advantage that it can capture any orientation of the magnetization and vary it
in time without losing its justification. In this way, magnon modes at _all_
wave vectors ${\bf k}$ are included in the description, not only the
macrospins corresponding to ${\bf k}={\bf 0}$ and ${\bf k}=(\pi,\pi)^{\top}$
or ${\bf k}=(\pi,\pi,\pi)^{\top}$, respectively. Thereby, the intrinsic
dephasing effects are included and no Gilbert damping needs to be introduced
to account for it.

The main focus was to study short pulses as they can be generated nowadays in
the THz regime. As expected such pulses are most effective if their carrier
angular frequency is chosen to be in the range of the spin gap $\Delta/\hbar$.
But we found that the optimum effect is obtained for a frequency lowered by
about 15%. We attribute this finding to the reduction of the spin gap in the
course of switching similar to the activation energy of a large anisotropic
spin [16]. We found that a short Gaussian pulse is a very suitable choice to
achieve switching of the magnetization fairly quickly. Unexpectedly, finite
temperature does not help much to improve the switching efficiency if the
pulses have been optimized.

Other choices of pulse shapes, for instance a pulse shaped like a sinus
cardinalis [16], require stronger fields. A continuous oscillating control
field can switch the magnetization with low amplitudes, but it takes fairly
long for the switching to take place. Also the Gaussian pulses require less
maximum amplitude if they are applied longer, i.e., if the pulse duration is
larger. On the downside, this makes the switching process slower.

A quantitative estimate strongly depends on the numbers. If we assume
$J\approx 10$ meV and $\chi=0.98$ the required Gaussian amplitude for
$\tau=10/J\approx 1$ ps is as large as about 10 T. If the pulse is applied
about four times longer and the anisotropy is smaller, let us say
$\chi=0.995$, then the estimate yields about 1 T. Still, these are fairly
large values so that further reductions are desirable.

A promising route consists in the use of alternating (Néel) effective magnetic
fields. Theoretical treatments based on macrospin theory [19, 39] and recent
experimental evidence [40] indicate that in this case the effective control is
not longer proportional to $h$, but to $\sqrt{hJ}$. This exchange enhancement
is similar to the enhancement of the spin gap in antiferromagnets which is
proportional to $\sqrt{\delta JJ}$ where $\delta J$ is the deviation of the
relevant couplings from the spin isotropic case, for instance $\delta
J=(1-\chi)J$. The extension of the present calculation to alternating Néel
fields is ongoing research. Furthermore, it is easy to conceive further
important extensions such as the treatment of other lattices in two and three
dimensions and/or higher spins $S>1/2$. Such studies will decisively help to
introduce antiferromagnets as systems for data storage and manipulation.

###### Acknowledgements.

We are grateful to Timo Gräßer, Dag-Björn Hering, and Joshua Althüser for very
useful discussions. This work has been financially supported by the Deutsche
Forschungsgemeinschaft (German Research Foundation) in project UH 90/14-1 and
by the Stiftung Mercator in project Ko-2021-0027.

## References

  * Chappert  _et al._ [2007] C. Chappert, A. Fert, and F. Van Dau, The emergence of spin electronics in data storage, Nature Materials 6, 813 (2007).
  * Marrows [2016] C. Marrows, Addressing an antiferromagnetic memory, Science 351, 558 (2016).
  * Jungwirth  _et al._ [2016] T. Jungwirth, X. Marti, P. Wadley, and J. Wunderlich, Antiferromagnetic spintronics, Nature Nanotechnology 11, 231 (2016).
  * Gomonay  _et al._ [2017] O. Gomonay, T. Jungwirth, and J. Sinova, Concepts of antiferromagnetic spintronics, Physica Status Solidi - Rapid Research Letters 11, 1700022 (2017).
  * Loth  _et al._ [2012] S. Loth, S. Baumann, C. P. Lutz, D. M. Eigler, and A. J. Heinrich, Bistability in atomic-scale antiferromagnets, Science 335, 196 (2012).
  * Gomonay and Loktev [2014] E. V. Gomonay and V. M. Loktev, Spintronics of antiferromagnetic systems (review article), Low Temperature Physics 40, 17 (2014).
  * Kampfrath  _et al._ [2011] T. Kampfrath, A. Sell, G. Klatt, A. Pashkin, S. Mährlein, T. Dekorsy, M. Wolf, M. Fiebig, A. Leitenstorfer, and R. Huber, Coherent terahertz control of antiferromagnetic spin waves, Nature Photonics 5, 31 (2011).
  * Song  _et al._ [2018] C. Song, Y. You, X. Chen, X. Zhou, Y. Wang, and F. Pan, How to manipulate magnetic states of antiferromagnets, Nanotechnology 29, 112001 (2018).
  * Mukai  _et al._ [2016] Y. Mukai, H. Hirori, T. Yamamoto, H. Kageyama, and K. Tanaka, Nonlinear magnetization dynamics of antiferromagnetic spin resonance induced by intense terahertz magnetic field, New Journal of Physics 18, 013045 (2016).
  * Jin  _et al._ [2023] Z. Jin, Z. Mics, G. Ma, Z. Cheng, M. Bonn, and D. Turchinovich, Single-pulse terahertz coherent control of spin resonance in the canted antiferromagnet YFeO3, mediated by dielectric anisotropy, Physical Review B 87, 094422 (2023).
  * Katine  _et al._ [2000] J. A. Katine, F. J. Albert, R. A. Buhrman, E. B. Myers, and D. C. Ralph, Current-driven magnetization reversal and spin-wave excitations in co/cu/co pillars, Physical Review Letters 66, 3149 (2000).
  * Liu  _et al._ [2012] L. Liu, C.-F. Pai, Y. Li, H. W. Tseng, D. C. Ralph, and R. A. Buhrman, Spin-Torque Switching with the Giant Spin Hall Effect of Tantalum, Science 336, 555 (2012).
  * Železný  _et al._ [2014] J. Železný, H. Gao, K. Výborný, J. Zemen, J. Mašek, A. Manchon, J. Wunderlich, J. Sinova, and T. Jungwirth, Relativistic Néel-Order Fields Induced by Electrical Current in Antiferromagnets, Physical Review Letters 113, 157201 (2014).
  * Wadley  _et al._ [2016] P. Wadley, B. Howells, J. Železný, C. Andrews, V. Hills, R. P. Campion, V. Novák, K. Olejník, F. Maccherozzi, S. S. Dhesi, S. Y. Martin, T. Wagner, J. Wunderlich, F. Freimuth, Y. Mokrousov, J. Kuneš, J. S. Chauhan, M. J. Grzybowski, A. W. Rushforth, K. W. Edmonds, B. L. Gallagher, and T. Jungwirth, Electrical switching of an antiferromagnet, Science 351, 587 (2016).
  * Kampfrath  _et al._ [2013] T. Kampfrath, K. Tanaka, and K. A. Nelson, Resonant and nonresonant control over matter and light by intense terahertz transients, Nat. Phot. 7, 680 (2013).
  * Miyashita and Barbara [2023] S. Miyashita and B. Barbara, How to cross an energy barrier at zero kelvin without tunneling effect, Physical Review Letters 131, 066701 (2023).
  * Bar’yakhtar and Ivanov [1980] I. V. Bar’yakhtar and B. A. Ivanov, Nonlinear waves in antiferromagnets, Solid State Commun. 34, 545 (1980).
  * Gomonay and Loktev [2008] H. Gomonay and V. Loktev, Spin-plarized Current-induced instability in the Spin-Valve Antiferromagnetic Layer, J. Magn. Soc. Jpn. 32, 535 (2008).
  * Gomonay and Loktev [2010] H. V. Gomonay and V. M. Loktev, Spin transfer and current-induced switching in antiferromagnets, Physical Review B 81, 144427 (2010).
  * Kim  _et al._ [2014] T. H. Kim, S. Y. Hamh, J. W. Han, C. Kang, C.-S. Kee, S. Jung, J. Park, Y. Tokunaga, Y. Tokura, and J. S. Lee, Coherently controlled spin precession in canted antiferromagnetic YFeO3 using terahertz magnetic field, Appl. Phys. Express 7, 093007 (2014).
  * Kim  _et al._ [2017] T. H. Kim, P. Grüenberg, S. H. Han, and B. K. Cho, Field-driven dynamics and time-resolved measurement of Dzyaloshinskii-Moriya torque in canted antiferromagnet YFeO3, Science Reports 7, 4515 (2017).
  * Nowak [2007] U. Nowak, Classical spin models, in  _Magnetism_, Handbook of Magnetism and Advanced Magnetic Materials, edited by H. Kronmüller and S. Parkin (J. Wiley & Sons, New York, 2007) p. 858.
  * Weißenhofer  _et al._ [2023] M. Weißenhofer, F. Foggetti, U. Nowak, and P. M. Oppeneer, Néel vector switching and terahertz spin-wave excitation in mn2au due to femtosecond spin-transfer torques, Physical Review B 107, 174424 (2023).
  * Kazantseva  _et al._ [2008] N. Kazantseva, D. Hinzke, U. Nowak, R. W. Chantrell, U. Atxitia, and O. Chubykalo-Fesenko, Towards multiscale modeling of magnetic materials: Simulations of fept, Physical Review B 77, 184428 (2008).
  * Stanek  _et al._ [2014] D. Stanek, C. Raas, and G. S. Uhrig, From quantum-mechanical to classical dynamics in the central-spin model, Physical Review B 90, 064301 (2014).
  * Polkovnikov [2010] A. Polkovnikov, Phase space representation of quantum dynamics, Ann. of Phys. 325, 1790 (2010).
  * Davidson and Polkovnikov [2015] S. M. Davidson and A. Polkovnikov, $su(3)$ semiclassical representation of quantum dynamics of interacting spins, Physical Review Letters 114, 045701 (2015).
  * Holstein and Primakoff [1940] T. Holstein and H. Primakoff, Field dependence of the intrinsic domain magnetization of a ferromagnet, Physical Review 58, 1098 (1940).
  * Dyson [1956] F. J. Dyson, General theory of spin-wave interactions, Physical Review 102, 1217 (1956).
  * Maleev [1958] S. V. Maleev, Scattering of slow neutrons in ferromagnets, Soviet Physics JETP 6, 766 (1958).
  * Schwinger [1952] J. Schwinger, On angular momentum, US Atomic Energy Commission 10.2172/4389568 (1952).
  * Auerbach and Arovas [1988] A. Auerbach and D. P. Arovas, Spin dynamics in the square-lattice antiferromagnet, Physical Review Letters 61, 617 (1988).
  * Auerbach [1994] A. Auerbach,  _Interacting Electrons and Quantum Magnetism_, Graduate Texts in Contemporary Physics (Springer, New York, 1994).
  * Bolsmann  _et al._ [2023] K. Bolsmann, A. Khudoyberdiev, and G. S. Uhrig, Switching the magnetization in quantum antiferromagnets, PRX Quantum 4, 030332 (2023).
  * De Silva  _et al._ [2002] T. N. De Silva, M. Ma, and F. C. Zhang, Pathology of Schwinger boson mean-field theory for Heisenberg spin models, Physical Review B 66, 104417 (2002).
  * Marder [2010] M. P. Marder,  _Condensed Matter Physics_ (John Wiley & Sons, New York, 2010).
  * Pires [2021] A. S. T. Pires,  _Theoretical Tools for Spin Models in Magnetic Systems_, 2053-2563 (IOP Publishing, 2021).
  * Hanisch  _et al._ [1997] T. Hanisch, G. S. Uhrig, and E. Müller-Hartmann, Lattice dependence of saturated ferromagnetism in the hubbard model, Physical Review B 56, 13960 (1997).
  * Roy  _et al._ [2016] P. E. Roy, R. M. Otxoa, and J. Wunderlich, Robust picosecond writing of a layered antiferromagnet by staggered spin-orbit fields, Physical Review B 94, 014439 (2016).
  * Behovits  _et al._ [2023] Y. Behovits, A. L. Chekhov, S. Y. Bodnar, O. Gueckstock, S. Reimers, Y. Lytvynenko, Y. Skourski, M. Wolf, T. S. Seifert, O. Gomonay, M. Kläui, M. Jourdan, and T.Kampfrath, Terahertz Néel spin-orbit torques drive nonlinear magnon dynamics in antiferromagnetic Mn2Au, Nature Communications 14, 6038 (2023).

##  Appendix A Dynamics of expectation values

First, we give the expectation values in equilibrium at temperature $T$. They
provide the initial values for the subsequent dynamic calculation. The thermal
equilibrium results from standard bosonic statistics after diagonalization of
the mean-field Hamiltonian (7) for $h_{y}=0$ by the appropriate Bogoliubov
transformation [34].

$\displaystyle\langle a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}$
$\displaystyle=\frac{z\gamma
C_{-}}{8\omega_{k}^{-}(\gamma)}\coth(\frac{\omega_{k}^{-}(\gamma)}{2T}),$
(18a) $\displaystyle\langle b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}$
$\displaystyle=\frac{z\gamma
C_{+}}{8\omega_{k}^{+}(\gamma)}\coth(\frac{\omega_{k}^{+}(\gamma)}{2T}),$
(18b) $\displaystyle\langle
a_{\mathbf{k}}^{\dagger}{a}_{\mathbf{k}}\rangle_{\gamma}$
$\displaystyle=\frac{\lambda}{2\omega_{k}^{-}(\gamma)}\coth(\frac{\omega_{k}^{-}(\gamma)}{2T})-\frac{1}{2},$
(18c) $\displaystyle\langle
b_{\mathbf{k}}^{\dagger}{b}_{\mathbf{k}}\rangle_{\gamma}$
$\displaystyle=\frac{\lambda}{2\omega_{k}^{+}(\gamma)}\coth(\frac{\omega_{k}^{+}(\gamma)}{2T})-\frac{1}{2},.$
(18d)

The constants $C_{\pm}$ are defined in Es. (9) and the dispersions in Eq.
(11). The variables $A$ and $B$ required to compute $C_{\pm}$ are defined by

$\displaystyle A$ $\displaystyle=\langle a_{i}a_{j}\rangle+\langle
b_{i}b_{j}\rangle$ (19a)
$\displaystyle=\int_{-1}^{1}\gamma\rho_{d}(\gamma)\left(\langle
a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}+\langle
b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}\right)d\gamma,$ (19b)
$\displaystyle B$ $\displaystyle=\langle a_{i}a_{j}\rangle-\langle
b_{i}b_{j}\rangle$ (19c)
$\displaystyle=\int_{-1}^{1}\gamma\rho_{d}(\gamma)\left(\langle
a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}-\langle
b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}\right)d\gamma,$ (19d)
$\displaystyle 2S$ $\displaystyle=\langle
a_{i}^{\dagger}{a}_{i}\rangle+\langle b_{i}^{\dagger}{b}_{i}\rangle$ (19e)
$\displaystyle=\int_{-1}^{1}\rho_{d}(\gamma)\big{(}\langle
a_{\mathbf{k}}^{\dagger}{a}_{\mathbf{k}}\rangle_{\gamma}+\langle
b_{\mathbf{k}}^{\dagger}{b}_{\mathbf{k}}\rangle_{\gamma}\big{)}d\gamma,$ (19f)

where the last equation is the average constraint on the number of Schwinger
bosons which needs to be fulfilled as well by adjusting the Lagrange
multiplier $\lambda$. The density $\rho_{d}(\gamma)$ is the density-of-states
in $d$ dimensions for $\gamma_{\mathbf{k}}$ in Eq. (15), see Ref. [38] for
concrete expressions in two and three dimensions.

Finally, the temporal evoluation is determined from the equations of motion
for the introduced expectation values. As stated before, this dynamics only
depends on the value $\gamma_{\mathbf{k}}=\gamma$

$\displaystyle\partial_{t}\langle
a_{\mathbf{k}}^{\dagger}a_{\mathbf{k}}\rangle_{\gamma}$
$\displaystyle=-i\frac{z}{4}\gamma\big{(}C_{-}^{*}\langle
a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}-C_{-}\langle
a_{\mathbf{k}}^{\dagger}a_{\mathbf{-k}}^{\dagger}\rangle_{\gamma}\big{)}$
$\displaystyle\quad+\frac{h_{y}(t)}{2}\big{(}\langle
a_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}\rangle_{\gamma}+\langle
b_{\mathbf{k}}^{\dagger}a_{\mathbf{k}}\rangle_{\gamma}\big{)},$ (20a)
$\displaystyle\partial_{t}\langle
b_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}\rangle_{\gamma}$
$\displaystyle=-i\frac{z}{4}\gamma\big{(}C_{+}^{*}\langle
b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}-C_{+}\langle
b_{\mathbf{k}}^{\dagger}b_{\mathbf{-k}}^{\dagger}\rangle_{\gamma}\big{)}$
$\displaystyle\quad-\frac{h_{y}(t)}{2}\big{(}\langle
a_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}\rangle_{\gamma}+\langle
b_{\mathbf{k}}^{\dagger}a_{\mathbf{k}}\rangle_{\gamma}\big{)},$ (20b)
$\displaystyle\partial_{t}\langle
a_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}\rangle_{\gamma}$
$\displaystyle=-i\frac{z}{4}\gamma\big{(}C_{-}^{*}\langle
a_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}-C_{+}\langle
a_{\mathbf{k}}^{\dagger}b_{\mathbf{-k}}^{\dagger}\rangle_{\gamma}\big{)}$
$\displaystyle\quad+\frac{h_{y}(t)}{2}\big{(}\langle
b_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}\rangle_{\gamma}-\langle
a_{\mathbf{k}}^{\dagger}a_{\mathbf{k}}\rangle_{\gamma}\big{)},$ (20c)
$\displaystyle\partial_{t}\langle
a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}$
$\displaystyle=i\frac{z}{4}\gamma\big{[}C_{-}(2\langle
a_{\mathbf{k}}^{\dagger}{a}_{\mathbf{k}}\rangle_{\gamma}+1)\big{]}$
$\displaystyle\quad-2\lambda i\langle
a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}+h_{y}(t)\langle
a_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma},$ (20d)
$\displaystyle\partial_{t}\langle
b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}$
$\displaystyle=i\frac{z}{4}\gamma\big{[}C_{+}(2\langle
b_{\mathbf{k}}^{\dagger}{b}_{\mathbf{k}}\rangle_{\gamma}+1)\big{]}$
$\displaystyle\quad-2\lambda i\langle
b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}-h_{y}(t)\langle
a_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma},$ (20e)
$\displaystyle\partial_{t}\langle
a_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}$
$\displaystyle=i\frac{z}{4}\gamma\big{(}C_{-}\langle
a_{\mathbf{k}}^{\dagger}b_{\mathbf{k}}\rangle_{\gamma}+C_{+}\langle
b_{\mathbf{k}}^{\dagger}a_{\mathbf{k}}\rangle_{\gamma}\big{)}$
$\displaystyle\quad-2\lambda i\langle
a_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}$
$\displaystyle\quad-\frac{h_{y}(t)}{2}\big{(}\langle
a_{\mathbf{k}}a_{\mathbf{-k}}\rangle_{\gamma}-\langle
b_{\mathbf{k}}b_{\mathbf{-k}}\rangle_{\gamma}\big{)}.$ (20f)

This concludes the set of required equations.

##  Appendix B Sudden shift in switching time due to a change in phase

Figure 12: Evolution of the magnetization for phases just before the jump with
$h_{0}=0.187\,J$, $\phi_{0}=0.22\pi$ (a) and just after the jump with
$h_{0}=0.158\,J$, $\phi_{0}=0.23\pi$ (b) in Fig. 6 for $\chi=0.9$. The time
shift where switching occurs due to the minute change in phase is denoted
$\delta t$.

We emphasize that the phase of the Gaussian pulse only has a weak effect on
the required field strength. Still, we found a small jump in Fig. 6. To
understand its origin, we provide two panels in Fig. 12 to show the dynamics
of the magnetization with phases just before and just after the jump. One can
see that the slightly higher phase with a lower value of the control field
allows for a switching at a later instant in time, shifted by $\delta t$. The
instant of switching occurs about 3/4 of a pulse cycle later. Hence, the pulse
takes a little longer to induce the sign flip of $m$ which in return allows
for the lower amplitude.

##  Appendix C Effect of alternative control fields in 2D

###  C.1 Equally weighted superposition of frequencies

A typical energy landscape with a parabolic activitation maximum between the
energy minima has been studied recently by Miyashita and Barbara [16]. They
showed that it is advantageous to design the pulse such that in its Fourier
transform the frequencies of all necessary level transitions are comprised
with equal weight. If one generalizes this to an infinite spin the Fourier sum
becomes an integral and the resulting pulse shape is a sinus cardinalis well-
known from the diffraction of a single slit; it reads

$h_{y}(t)=h_{0}\frac{\sin(\alpha\Delta(t-\tau))}{\alpha\Delta(t-\tau)}.$ (21)

Thus, we studied this pulse shape as well. Fig. 13 displays a representative
result for the evolution of the magnetization for the chosen anisotropy.

Figure 13: Evolution of the magnetization $m(t)$ and the expectation value
$\langle S^{x}\rangle$ induced by the pulse Eq. (21) for $h_{0}=0.36\,J$,
$\alpha=1.1$ and $\tau=30\,J^{-1}$ at zero temperature for the square lattice.

Indeed, this pulse as well can achieve switching of the magnetization. On the
one hand, the required amplitude $h_{0}$ is lower than the threshold field for
a static control field $h_{\mathrm{thresh}}=0.863\,J$ from Ref. 34. On the
other hand, the amplitude needed for switching is significantly larger than
the one for the optimized Gaussian pulse (17)
$h_{0}^{\mathrm{thresh}}=0.1578\,J$, cf. Fig. 2. Thus, we conclude that the
pulse shape (21) is not optimum. In addition, it is experimentally not easily
realizable where short pulses consist of very few cycles only.

For completeness, Fig. 14 shows the threshold field dependence on the
frequency renormalization $\alpha$ in (21). One can see that there is almost
no shift away from the resonance condition $\alpha=1$ resonance in contrast to
the Gaussian pulse in (17). It appears that the optimum $\alpha$ lies a bit
above unity; for instance, we found $\alpha\approx 1.1$ with
$h_{0}^{\mathrm{thresh}}=0.357\,J$ for $\chi=0.9$. We attribute the finding
that $\alpha<1$ is not optimum to the fact that the pulse shape (21) already
comprises all frequencies below $\alpha\Delta$ in its Fourier transform. Thus,
no additional lowering is needed.

Figure 14: The threshold field dependence on $\alpha$ for pulse (21). Here,
the delay is chosen to be $\tau=30\,J^{-1}$.

###  C.2 Continuous oscillating control field

Another conceivable pulse applies a continuous electromagnetic wave to the
sample. Within the framework of our model this corresponds to the application
of

$h_{y}(t)=h_{0}\cos(\alpha\Delta t).$ (22)

Figure 15: Switching in square lattice by long lasting control field with
amplitude $h_{0}=0.055\,J$ in (22), $\Delta=\Delta(\chi=0.9)$ and
$\alpha=0.95$ . Figure 16: The threshold value $h_{0}^{\mathrm{thresh}}$
dependence on renormalization $\alpha$ of the continuous oscillating field in
(22).

Fig. 15 displays a generic result for the continuous driving (22) which
illustrates successful switching for a fairly small amplitude, even lower than
the threshold values given in Fig. 5. But the switching takes very long so
that this choice of control does not appear appealing in practice, besides its
difficulty in realization.

For completeness, we also analyzed the effect of a renormalization $\alpha$ in
Fig. 16. The effect is the similar to the effect for the Gaussian pulse (17).
A value $\alpha\lessapprox 1$ allows us to reduce the amplitude slightly. But
the overall dependence on $\alpha$ is not very regular and deviates from a
parabolic shape even around its minimum. We refrain from a more detailed
analysis since the continuous driving is not easily accessible experimentally
and our findings indicate unnecessary long switching times.