Protocol for generating an arbitrary quantum state of the magnetization in cavity magnonics

We propose and numerically evaluate a protocol to generate an arbitrary quantum state of the magnetization in a magnet. The protocol involves repeatedly exciting a frequency-tunable superconducting transmon and transferring the excitations to the magnet via a microwave cavity. To avoid decay, the protocol must be much shorter than magnon lifetime. Speeding up the protocol by simply shortening the pulses leads to non-resonant leakage of excitations to higher levels of the transmon accompanied by higher decoherence. We discuss how to correct for such leakages by applying counter pulses to de-excite these higher levels. In our protocol, states with a maximum magnon occupation of up to $\sim9$ and average magnon number up to $\sim4$ can be generated with fidelity $>0.75$.

Magnets have found commercial applications in magnetic field sensing and storing of information, and are a promising building block for long-range information transfer [1,2] and low-power logic devices [3]. Recently, there has been an interest in bringing such applications into the quantum domain, known as 'quantum magnonics' [4,5,6,7]. This is partly fuelled by the extremely low magnetic dissipation found in the ferrimagnetic insulator Yttrium Iron Garnet (YIG) [8], along with evidence for macroscopic (mm-long) coherence lengths [9,10]. In the limit of small deviations of spins from their equilibrium value, the collective excitations of a magnet can be modelled as a set of harmonic oscillators called magnons with a lifetime characterized by the Gilbert damping constant of the material that is particularly low for YIG [8]. Magnons can be probed via their interaction with microwaves (MWs) [11,12,13] or optical light [14,15,16,17]. This interaction can be enhanced by employing an electromagnetic cavity to confine the photons, whose decay rate is given by the quality factor of the cavity. The strong coupling regime between cavity MW photons and magnons, where information can be coherently exchanged between the two excitations at rates > 10 MHz much faster than each of their decay rates < 1 MHz, can be routinely achieved in experiments. The strong coupling regime in the optical domain is notoriously more difficult to achieve, however there are promising developments in this direction [18,19,20,21,22,23].
A quantum platform involving magnets requires generation, manipulation, and detection of non-classical states of the magnetization. Possible applications in the quantum regime include magnon-based quantum transducers between light and MWs, as well as transduction involving other degrees of freedom such as phonons [24,25] and electrons [26,27,28]. Furthermore, non-classical magnetization states can be useful in magnetic field sensing reaching the Heisenberg limit [29]. While single-magnon detection has been demonstrated [30,31,32] by exploiting the coupling of magnets to superconducting circuits mediated by MW cavities [4,33,34,35], experiments so far have probed coherent or thermal magnon states whereas "true" non-classical states (defined e.g. by a Wigner function with negative regions) remain to be realized. Theoretical proposals in this direction include the generation of cat states of the magnetization using MW cavity photons, [36] and the generation of single magnon Fock states [37] and cat states [38] via optical means. In thin films, the generation of entangled pairs of traveling magnons was proposed [39]. All of these proposals are probabilistic for the state generation (except for the last) and are specialized to their respective target states, while there is no known method to deterministically generate an arbitrary magnon state.
In this manuscript we fill in this gap by proposing a protocol to this end, tailored for magnets coupled to a superconducting transmon via a MW cavity, as depicted in Fig. 1. Our method is based on the protocol proposed in Ref. [40], which shows that the deterministic generation of an arbitrary quantum state of a harmonic oscillator can be accomplished by coupling the oscillator to a qubit. A superconducting transmon is a weakly anharmonic oscillator whose first two levels can be treated as a qubit. We show how the protocol of Ref. [40] can be modified to take into account the anharmonicity of the superconducting transmon, in order to generate arbitrary quantum magnon states with high fidelity. The paper is organized as follows. We discuss the system and effective Hamiltonian in Sec. 1. In Sec. 2 we approximate the superconducting transmon by a two-level system and review the protocol from Ref. [40] for the cavity mangnonic system. In Sec. 3 we discuss how to correct for errors due to deviations from the two-level approximation, in particular due to interferences from higher energy levels present in the transmon. In Sec. 4, we discuss the fidelity of state generation for paradigmatic examples of non-classical states. We discuss possible improvements to our protocol that can increase the fidelity of generation further in Sec. 5. We review the main findings qualitatively in Sec. 6.

System
The hybrid system we considered is depicted in Fig. 1: inside a MW cavity, magnons can couple to a superconducting transmon [4,33]. The effective coupling Hamiltonian can be obtained by tracing out the MW field, as we review in this section.
When the magnet is much smaller than the wavelength of the MWs (∼ cm), only the total magnetization M couples to the MWs, so we need to consider only the uniform magnon mode known as the Kittel mode. It can be quantized via a Holstein-Primakoff transformation [41,42] (see Appendix A for details) that, in the limit of M x,y M z , reduces to whereM is the magnon annihilation operator satisfying the bosonic commutation relations [M ,M † ] = 1, M s is the saturation magnetization, and the zero-point fluctuations of the magnetization are given by with γ > 0 being the absolute gyromagnetic ratio and V m being the volume of the magnet. Assuming a spherical magnet and an external magnetic field H app , the Hamiltonian becomesĤ m =hΩ mM †M with Ω m = γµ 0 H app (see Appendix A) which can be tuned to be in the range of ∼ 2π × 5-10 GHz. For a YIG sample with volume V m > 1 µm 3 , M s = 140 kA m −1 , and γ = 2π×28 GHz T −1 , we have M ZPF /M s < 8×10 −6 , implying that the deviations from the ground state are small.
A flux-tunable transmon consists of two Josephson junctions in parallel with Josephson energy E J1,J2 respectively, forming a SQUID loop along with a capacitor with capacitance C also in parallel [43,44,45,46]. When the charging energy E C = e 2 /2C is such that E C E J1,J2 , the energy levels of the circuit resemble that of a weakly anharmonic oscillator, where A = E C /4 and Ω q = E C E J,eff with the effective Josephson energy E J,eff being Here Φ is the flux inside the SQUID loop and Φ 0 = h/2e is the superconducting flux quantum. Such a set of energy levels can be labelled using an annihilation operatorQ and its corresponding number operatorQ †Q . The number of quanta exchanged between the two Josephson junctions is given bŷ The Hamiltonian becomeŝ The frequency Ω q can be tuned dynamically via the external flux in SQUID Φ, to be in a typical range of 2π ×5-20 GHz. Tuning via Φ can be achieved at a fast timescale 10 ns [47,48], but it comes at the cost of increased decoherence proportional to the allowed tunability |E J1 − E J2 |, as a higher tunability implies larger susceptibility to flux noises [49,45]. For our case, we assume that Ω q is tunable within a range of ∼ 2π × 1 GHz.
We assume that the magnet and the transmon are both placed inside a MW cavity which can be modelled as a set of harmonic modes with frequencies Ω a,r and annihilation operatorÂ r . The details of the modes depend on the shape and size of the cavity. The Hamiltonian for the unloaded cavity isĤ mw = rh Ω a,rÂ † rÂ r . Its coupling to both the magnet and the transmon corresponds to exchange of quanta and can be modelled as [see Appendix A] where the coupling G m,r ∝ |B(r magnet )| with B(r magnet ) the magnetic field at the location of the magnet, and G q,r ∝ |E(r transmon )| where E(r transmon ) is the electric field at the location of the transmon. If for all r modes, Ω a,r −Ω m > G m,r and Ω a,r −Ω q > G q,r , we can define cavity-dressed magnon and transmon modes with annihilation operatorŝ m andq respectively. Up to quadratic order in the couplings [4], Details about this Hamiltonian are provided in Appendix A, while here we list down the features relevant for the rest of the article. The renormalized magnon's frequency ω m and transmon's frequency ω q stay close to their bare values Ω m and Ω q respectively. A schematic representation of each term of the Hamiltonian is shown in Fig. 1.
The functionα gives the anharmonicity of the transmon andĝ gives the magnontransmon coupling. Below, we use the same characters α(n) = n|α|n and g(n) = n|ĝ|n for the values of anharmonicity and coupling. Typically α(2) ∼ 2π × 150-300 MHz increasing with n and the coupling g can achieve ∼ 2π × 20-30 MHz depending on the system parameters as discussed in Appendix A.
Restricting the analysis to only the first two levels of the transmon, we get the Jaynes-Cummings Hamiltonian where g ≡ g(0), {|g , |e } are the first two excited states of the transmon,σ + = |g e| andσ − = |e g|. As mentioned before, the frequency ω q is tunable. When |ω q − ω m | g, the magnons and transmon are decoupled besides a negligible dispersive coupling, while at resonance, ω q = ω m , the two subsystems can exchange quanta.

Ideal State Generation
In this section we review the protocol to generate an arbitrary state of the magnetization via its coupling to a two-level system. This is discussed in [40] theoretically for a general harmonic oscillator. The protocol consists of: (1) exciting the qubit; (2) partially transfering the excitation to the oscillator and (3) repeating the process multiple times to construct the Fock-state superposition that describes the target state. In this setup, the transmon can be excited via a time-dependent drive, while the transmon and the magnons can be put in and out of resonance by tuning the magnon frequency (via an external magnetic field) or the transmon frequency (via an applied flux), effectively changing the coupling rate (see Fig. 2). The procedure was experimentally implemented to generate quantum states of MW photons coupled to a superconducting transmon [48].

Magnon Transmon
Magnon Transmon If we start from the qubit's state c 0 |e + c 1 |g while the magnon is in the ground state |0 , switching the coupling on for a period t = π/(2g) in the absence of dissipation evolves the state as where |n is n-magnon state. Thus, we can create an arbitrary superposition of |0 and |1 . If we repeat this process N times, we can, in principle, generate an arbitrary magnon state containing up to a maximum of N magnons. The protocol then consists of a series of operations in which first the transmon is excited, and then the qubit's frequency is changed to bring it in resonance with the magnons for a time of the order π/g. To make the process concrete, consider the Hamiltonian restricted to the first two levels of the transmonĤ (2) (t) =Ĥ (2) s +Ĥ = ω mm †m + ω q |e e| + g mσ + +m †σ − , and the dynamic part iŝ Here and in what follows, the superscript (2) indicates the approximation of the transmon as a two-level system. In this expression, ∆(t) is the externally induced change in the transmon's frequency, andε(t) is the transmon drive amplitude. We assume that the qubit and the magnons are far detuned, |ω m − ω q | g, implying that magnons and qubit are decoupled at ∆ = 0.
A non-zeroε can excite a qubit. We define the time evolution operatorQ (2) for ∆ = 0 andε(t) = ε(t)e −iωqt+iφ , where t g is a pre-determined gate time, ε(t) is chosen satisfying tg 0 dτ ε(τ ) = θ and Texp is the time-ordered exponential. As derived in Appendix C, which applies a Bloch rotation of magnitude 2θ around the axis (cos φ, − sin φ, 0). After the qubit is excited, we want to (partially) transfer the quanta to magnons. We then define the time evolution operator for a detuning ∆ 0 = ω m − ω q and no qubit drivinĝ An additional 'off' time, t off , adjusts the relative phases in the wave-function ‡. As derived in Appendix B,T (2) =T (2) phT (2) ex where the phases evolve aŝ T and the quanta are exchanged viâ with the total number operatorN =m †m + |e e|. A given target state |ψ T containing a maximum of N m magnons, can be systematically achieved by adding magnons one by one as described above. For that, we need a set of time evolution operatorsQ (2) j ≡Q (2) and each |Ψ j has a maximum of j excitations. Explicitly, for |Ψ j = ψ e j , e + ψ g j , g there are a maximum of (j − 1) and j magnons in ψ e j and ψ g j respectively. The coefficients {θ j , φ j , t off,j , t on,j } are found by reversing the problem in order to systematically remove magnons one by one from |Ψ Nm . Inductively, we assume that |Ψ j+1 is found. Explicitly, we find {t on,j , t off,j } to ensure j j + 1, g|Ψ j+1 cos gt on,j j + 1 + i j, e|Ψ j+1 sin gt on,j j + 1 e −i∆ 0 t off,j = 0. (19) Next, we find {θ j , φ j } to ensure j, e Q (2), † jT (20) ‡ In [40], the coupling was assumed to be a freely varying complex parameter, so the phases could be adjusted via the phase of g.  Then, we find |Ψ j = T (2)

jQ
(2) j † |Ψ j+1 and the induction continues until all quanta are removed. The procedure to inductively obtain the protocol's parameters is schematically shown in Fig. 3.

Realistic State Generation
The protocol discussed in the previous section is suitable to generate states with only a small number of excitations. States with higher number occupation require longer protocols and hence are spoiled by dissipation. The protocol duration can be decreased via two ways. Firstly, the gate time to excite the qubit can be decreased. Nevertheless, for a transmon qubit, a shorter gate time, t g ∼ 1/α(n) where α is the level-dependent anharmonicity [see below Eq. (8)], implies a broader frequency spectrum allowing for significant non-resonant transitions to higher transmon levels. Secondly, a stronger coupling between magnons and the transmon speeds up the excitation transfer step. However, a very strong coupling, g(n) √ n ∼ α(n) where g is the level-dependent magnontransmon coupling [see below Eq. (8)], implies a leakage of magnon excitations into higher transmon states. Such interferences limit the protocol time and eventually the largest size of the harmonic oscillator state that can be created. In this section, we show that transitions to higher levels can be incorporated by adding extra pulses for canceling out such transition amplitudes. The extra pulses are analogous to the ones used for correcting qubit gates [50,51] and can significantly decrease the protocol time allowing for states containing a large magnon number.
For clarity in presentation, we consider the case of g(n) ≡ g and while the calculations for the general case are presented in Appendix B and Appendix C. The Hamiltonian, Eq.
and the dynamic part iŝ For a target magnon state |ψ , we want to find a set of operationsQ j andT j satisfying (see Fig. 2) where |ψ 1 , ψ 2 refers to the magnon state |ψ 1 and transmon state |ψ 2 . As before, starting from |Ψ Nm , we inductively find {T j ,Q j , |Ψ j } where the transfer operator ensures j + 1, 0 T † j Ψ j+1 = 0 (removing a magnon) and the gate operation ensures for s ≥ 1 (removing transmon excitations).

Transfer of excitation
In Appendix B, we deriveT j satisfying the condition j + 1, 0 T † j Ψ j+1 = 0, while we discuss the salient features of the derivation here. We start with the same ansatz as before [see Eq. (15)] except for the fact that the Hamiltonian includes all higher transmon levels. This ansatz can be simplified by noticing that the number of excitationsN =m †m +q †q are conserved during the evolution, i.e. N ,T ans =0. An expression forT ans restricted to the space of {|j + 1 − s, s } s≥0 can be found by numerically diagonalizing the two Hamiltonians,Ĥ s andĤ s + ∆ 0q †q . Finally, we get a condition on {t on , t off } of the form where {λ (p) on , λ The above ansatz yields a solution in most cases, nevertheless, it can fail for some states, for example, when |Ψ j+1 = |j + 1, 0 . Intuitively, when the coupling is on, the excitation from |j + 1, 0 transfers to all |j + 1 − s, s but with different oscillation frequencies ∼ g 2 + α 2 s 2 (s−1) 2

4
. Unless these oscillation frequencies are commensurate, the amplitude of the wave-function in |j + 1, 0 will stay non-zero.
In cases where this ansatz fails, we can apply a second operation, where we choose {t off,j } as the values found from a two-level approximation [see Eq. (20)], and {t (2) on,j , t (2) off,j } are found by applying the above optimization tô T ans, † t off,j |Ψ j+1 . In principle, this process can be continued until a solution is found. However, in our simulations, we never had to go beyond a second step.

Transmon gates
As noted before, the application ofQ j ensures thatQ † j Ψ j has no excitations above the level j. The results are derived in Appendix C and here, we discuss the salient features of the derivation. First, we want to apply a gate between transmon states |0 and |1 ,Q ideal 01 (θ 01 j , φ 01 j ), to "transfer" the amplitude of |j, 1 to |j, 0 , i.e. ensure j, which can be obtained using the same arguments as in the two-level approximation, see Eq. (20). To implement this with realistic signal, we use the ansatẑ where {ε(t), φ app } are to be determined. We can separate the phases bŷ where we ignore the interaction between magnons and photons (they are far off-resonant) andÛ In the two-level approximation, we can choose any For concreteness, we choose where t g is a pre-determined gate time. In the two-level approximation, the timedependence is given bŷ Protocol for generating an arbitrary quantum state of the magnetization in cavity magnonics11 As the term in the exponential commutes with itself at all times, we can remove the time-ordering operator and simplifyÛ (0) 01 to, where Θ(t) = t 0 ε (0) (τ )dτ . However, a drive amplitude with a general form ε(t) can induce transitions to higher transmon levels, e.g. transitions between |j, 1 and |j, 2 . This is problematic as a significant transition amplitude j, 2 Q † j Ψ j implies excitation number j + 2, while we wantQ † j Ψ j to have no excitations above j. Thus, we consider are higher frequency components that cancel out higher order transitions, The error in the transmon gate is given bŷ where U err cannot be found analytically. To linear order in the perturbation, we ignore the corrections due to non-commutativity, and impose tg 0 dτ δĥ(τ ) =0. (38) This in turn implies the following set of equations that the drive amplitude ε(t) must satisfy The first condition is a no-bias condition ensuring that tg 0 ε(τ )dτ = θ 01 j . The second condition, with Θ(t) = t 0 dτ ε (0) (τ ) ensures that the transitions |1 ↔ |2 cancel out at the end of the gate. Similarly, the transitions |n ↔ |n + 1 are cancelled by the third condition. In our simulations, we find that for typical parameters, the population beyond level 2 is negligible, so we can ignore the conditions for n ≥ 3.
The above conditions leave sufficient freedom in choosing ε(t). In what follows, we choose t g = 3π/α and where the parameters ε µ are governed by an underdetermined set of linear equations that can be found from the above conditions [shown explicitly in Appendix C]. Given ε(t), we can find φ app by minimizing j, 1 Q † 01 Ψ j as discussed in Appendix C.3.
We point out that although we choose a linear approximation, the corrections are not always small. However, in our simulations we find that the error matrixÛ err is dominated by extraneous phases for all transmon levels. One might think that these can be removed by adding higher order corrections, however, we observed that such higher order corrections typically shifts the errors from |2 to |n ≥ 3 instead of removing them §. Thus, we account for these extra phases by modifying the next steps of the induction.
If j − 1, 2 Q † 01 Ψ j is significant (we take the threshold as 0.1 in our simulations), a second gate between |1 and |2 , sayQ ideal 12 (θ 12 j , φ 12 j ), can be applied to ensure that Essentially, the same arguments presented above hold for the gateQ 12 defined bŷ where we can choose and ε(t) satisfying the same conditions as before. The levels |n ≥ 3 are found to be approximately unoccupied, so we do not need to apply more gates.

Results
To evaluate the above algorithm, we numerically solve the time evolution that implements the protocol for the density matrix, using QuTiP [52,53], satisfying For dissipation, we take the Lindblad operators pertaining to magnon decayĈ 1 = √ κ mm , transmon decayĈ 2 = √ κ qq and transmon dephasingĈ 3 = √ γ qq †q . We did not include magnon dephasing as it has not yet been experimentally observed, although there are theoretical reasons for pure dephasing of magnons [54]. We choose typical parameters ω m = 2π×6 GHz, ω q = 2π×5 GHz, κ m = 10 −4 ω m , κ q = γ q = 10 −4 ω q , A = 2π×300 MHz, and coupling g = 2π × 25 MHz. For a given target state |ψ T , the algorithm described in the previous sections gives the signals ∆ Finally, we calculate the fidelity of state generation F [ψ T ] = ψ T |ρ mag [ψ T ]|ψ T which can be interpreted as the probability of finding the magnetization in the target state |ψ T .
We compare the Wigner function of the target state |ψ T and the achieved statê ρ mag [ψ T ] for some special cases. The Wigner function is defined as The Wigner function can be interpreted as follows. Define a symmetric product m †,pmq S as an average over all permutations of pm † 's and qm's. For example, In a given stateρ, the average of the above operator is given by [55,56] Tr ρ m †,pmq We can convert this into expectation values of the magnetization operatorsM x andM y via the substitution 2M ZPFm →M x − iM y (note that we are using the dressed operator m instead of the bare magnetization operatorM as the distinction between them is small) getting whereW (M x , M y .ρ) = W Mx−iMy 2M ZPF ,ρ . The above equation suggests an interpretation ofW as a quasi-classical probability distribution. However, note thatW is not always positive and a negativeW implies a lack of a probabilistic interpretation considered to be a signature of non-classicality. The unit M ZPF , defined in Eq. (2), is given by where V m is the volume of the magnet. This is a very small angle but it can be measured via the magnet's coupling to the transmon [31,32]. In Fig. 4, we consider the case of Fock states |ψ Fock (n) = |n and their superposition with the ground state |ψ Sup (n) ∝ |0 + |n where n ∈ {1, 2, . . . , 14}. The fidelity decreases with increasing n. This is caused by two factors. The dominant factor is an increase in the protocol time with the number of magnons leading to higher decay. Secondly, each gate introduces a small amount of error, including residual transfer to higher levels, that compounds with the number of gates, linearly increasing with n. It appears that Fock states are more difficult to generate than the superposition states, which could be due to a larger average number of magnons in |ψ Fock , compared to |ψ Sup , leading to a faster decay. For n ∈ {6, 10} and T ∈ {Fock, Sup}, we also show Wigner functions of the target state,W (M x , M y , |ψ T (n) ψ T (n)|) and the achieved statẽ W (M x , M y ,ρ mag [ψ T (n)]). The essential feature of |ψ Sup (n) are the azimuthal fringes  of order n, which is well captured for n = 6 and weakly reproduced for n = 10. The essential feature of |ψ Fock (n) are the n circular fringes, that are destroyed by dephasing.
A superposition of two semi-classical states, known as a cat state, is useful for error correction protocols because of their insensitivity to particle loss noise [57], particularly as carriers of information [58,59,60] and in metrology [61,62,63,64]. In case of a magnet, we can define semi-classical states as the minimum-uncertainty states with isotropic quantum fluctuations, given by coherent states [65] for an arbitrary complex number ξ. A coherent state has average m = ξ and can be interpreted as the quantum state 'closest' to the classical state with magnetization M x = 2M ZPF Re (ξ) and M y = −2M ZPF Im (ξ). An 'odd' cat state given by |ψ odd (ξ) ∝ |ξ/2 − |−ξ/2 has two semi-classical components separated by the magnetization 2M ZPF |ξ|; the same holds for an 'even' cat state |ψ even (ξ) ∝ |ξ/2 + |−ξ/2 . Fig.  5 shows the fidelity of generation as a function of ξ. The Wigner function for an ideal cat state has two blobs at the expected peaks along with fringes in the center. The presence of the fringes signify the quantum coherence between the two semiclassical components. In contrast, the Wigner function for a classical distribution, ρ cl (ξ) ∝ |ξ/2 ξ/2| + |−ξ/2 −ξ/2| would still have the blobs but no fringes.
In the achieved states, we find such fringes for cat sizes up to 4 until which the 1−e −ξ 2 /2 , than even cats, ξ 2 4 1−e −ξ 2 /2 1+e −ξ 2 /2 , an effect that is negligible for ξ > 2. Fidelity for generation of any other state can be found using the available codes [66].

Improvements
Here, we focussed on adapting and characterizing the protocol for a system with a magnet and a transmon. There are several methods by which the protocol can be improved, as we now elaborate.
Pre-processing: We considered only a specific type of operations, where the transmon is excited and the excitation is transferred to the magnons. With the help of the MW cavity, one can also apply other operations, in particular, a displacement given byD(α) = e αm † −α * m caused by a coherent excitation of the magnet by MWs. This can be used to decrease the time of the protocol as one can choose an appropriate α, such thatD † (α) |ψ T has fewer magnons than |ψ T . For example, if ψ T |m|ψ T = C, then implying thatD † (C) |ψ T has fewer average magnons than |ψ T . A better α can be found numerically that decreases the highest amplitude in |ψ T to within some given tolerance. Then, we run the protocol to createD † (α) |ψ T and apply a MW pulse to displace keeping a high fidelity. If there are other operations that can be easily performed, they can be used for pre-processing to decrease the size. High fidelity gates: In our analysis, we find that non-ideal transmon gates are responsible for a significant part of fidelity loss which increases linearly with the number of gates necessary to generate the target state. We took a simple case of correcting at first order perturbation similar to [50,51]. Simply adding higher order terms in a Dyson or a Magnus expansion does not necessarily improve gate fidelities. When we correct for second order term in Magnus expansion, we find that we can decrease the transitions |0 ↔ |2 and |1 ↔ |2 significantly, but at the expense of a high |0 ↔ |3 transition probability . This gets corrected only at the third order in Magnus expansion, at which point the |0 ↔ |4 transition might play a role. However, there is a rich literature on speeding up gates to ∼ 2π/α with low loss of fidelity, such as DRAG [67,68] or optimal control [69,70,71,72], that can not only speed up but also improve the accuracy of gates.
Combining steps: We followed a semi-analytical approach where we systematically add an excitation one by one. In principle, some steps can be combined, e.g. exciting the transmon and transferring the excitation to magnon can occur simultaneously. This becomes infeasible to handle analytically but relevant pulses can, in principle, be found using optimal control theory. A similar case has been discussed in the twolevel approximation [73] for a general harmonic oscillator states with 2 − 3 quanta. Although this approach should result in the least possible errors, it is computationally very expensive and scaling such an approach to high magnon numbers [74,70] would still require a good initial guess for time-dependent functions describing the transmon driving and the transmon frequency.

Conclusions
We proposed and numerically demonstrated a protocol to generate an arbitrary quantum state of the magnetization using its coupling to a superconducting transmon, which is a weakly anharmonic oscillator. The protocol involves repeatedly applying transmon gates and transferring an excitation to the magnons. As longer protocol times lead to higher dissipation, we can improve fidelity of generation by reducing the gate time, t g , and the transfer time, t on , in the notation of Secs. 2 and 3. Fast gate pulses with duration ∼ 2π/α (α is defined in Eq. (21)) would imply leakage of excitation into higher transmon levels that can be reduced by modifying the input pulses [47,69,70,67]. We use correction pulses by removing errors at first order in a Magnus expansion, similar to Refs. [50,51], bringing the gate time down to 3π/α. Correcting for second order term in the expansion does not seem to decrease errors, but rather transfer the population to higher levels. The transfer time, t on , can be reduced with a higher magnon-transmon coupling. As the magnon-MW coupling can reach very high values, it is likely to be In Refs. [50,51], only a three level system was considered which does not have this problem experimentally feasible to have high magnon-transmon couplings. However, a coupling g ∼ α causes leakage of magnons into higher transmon levels leading to loss of fidelities, so we choose g < 0.1α. For such couplings, the leakage is small and can be incorporated by renormalizing the duration of the coupling, t on,j .
The main source of errors is the magnon dissipation putting a limit on the size of the states. We find that the protocol faithfully produces states up to a largest occupation of ∼ 9 magnons and an average of ∼ 4 magnons. We note that our protocol creates states in dressed magnons, given bym instead of bare magnonsM . As discussed in Appendix A, for a sufficiently large detuning between MWs and magnons, the distinction is small. We demonstrated our protocol for cases of Fock states |n , Fock superpositions |0 + |n , and cat states |α ± |−α . While our protocol gives good fidelities for states of moderate sizes, we expect that advanced numerical techniques such as optimal quantum control [74,70,71] can improve the state generation (likely by reducing the protocol time) by using our protocol as an initial guiding guess.
With a setup to generate an arbitrary magnetization, there still lies a question of how to experimentally measure it, which we leave to future work. The figures can be generated via the publicly available codes [66].
In this appendix, we derive the effective Hamiltonian for a magnet coupled to a transmon indirectly via a MW cavity.

Appendix A.1. Microwaves
The MW cavity is modeled by the classical Hamiltonian density It contains several modes labelled by r giving the quantization and similar holds for E. Here, H r andÂ r are respectively the mode profile and the annihilation operator of mode r. With appropriate normalization, this gives the where Ω a,r is the resonance frequency of mode r.

Appendix A.2. Magnet
If the magnet is much smaller than MW wavelengths (∼ cm), then only the total magnetization M couples to the MWs. Under a DC field along z with magnitude H app , its classical Hamiltonian density is where H is the magnetic field due to the cavity. The magnetization is quantized by the Holstein-Primakoff transformation [41,42], given by a form where M ZPF and β are to be found. We know that a spin in an applied field undergoes a Larmor precession with frequency Ω m = γµ 0 H app where γ is the absolute value of the gyromagnetic ratio. This implies that the Kittel mode should have the resonance frequency, i.e. after quantization Typically M sing /M s < 10 −10 , so we can resort to the leading order corrections, For M †M 1/β ∼ 10 10 , we can simply substitute where Ω m = γµ 0 H app and We assumed rotating wave approximation which is valid when G Ω m + Ω a .

Appendix A.3. Transmon
We consider a flux-tunable transmon that consists of two Josephson junctions in parallel forming a SQUID loop. Its Hamiltonian is [45,44] H tr,mw = E C (n +n mw Here,n is the number of quanta exchanged between the two junctions,n mw is the number of quanta excited by MWs, andφ 1,2 are the fluxes through each junction. Under weak coupling, we haven mw ∝ |E|. Flux quantization constraintsφ 1 −φ 2 = 2πm + 2πΦ/Φ 0 where m is an arbitrary integer, Φ is the flux through the SQUID loop, and Φ 0 = h/2e is the superconducting flux quantum. Defining 2φ =φ 1 +φ 2 , where the phase is and the effective junction energy is In the transmon limit E J,eff E C , the fluctuations inφ + φ 0 are small, so we can expand the cosine and get an approximately quadratic Hamiltonian. This suggests the quantization,n that reduces the Hamiltonian tô ignoring higher order terms. HerehΩ q = 2E C E J,eff − E C /4 andhA = E C /4.

Appendix A.4. Dressed magnons and transmons
The dressed states are found by diagonalizing the total HamiltonianĤ =Ĥ 0 +Ĥ int whereĤ 19) and the interaction However, exact diagonalization is not possible analytically, so we resort to a perturbative treatment, in particular a Schrieffer-Wolff (SW) transformation [75]. The SW transformation perturbatively gives the rotation operatorŜ = −Ŝ † such that e −ŜĤ eŜ is diagonal to a given order. We findŜ to linear order in coupling constants G that gives the rotated Hamiltonian, To remove first order corrections, we choosê where the transmon susceptibility is , where δΩ a,r [x] = |G m,r | 2 Ω a,r − Ω m + |G q,r | 2 ((x + 1) χ r (x + 1) − xχ r (x)) . Here,Ĥ mw,eff is the effective MW Hamiltonian that is not relevant for our discussion. The renormalized magnon's frequency is Stark shifted by the MWs: Similarly, the transmon's frequency is The Stark shift of higher transmon levels renormalizes the anharmonicity tô Finally, the coupling of each transmon transition to the magnons is given bŷ

Appendix B. Transfer of Quanta
In this appendix, we show how to obtain the periods required for transferring quanta from transmon to magnons. In what follows, |m, q denotes a state with m magnons and transmon at qth level. Given a state |Ψ having a maximum of j excitations, i.e. m, q|Ψ = 0 for m + q > j and a tunable Hamiltonian g nm † |n n + 1|+g * nm |n + 1 n| ,(B.1) we want to find {∆(t), t tot } s.t. j, 0 T † Ψ = 0, where the time-evolution operator iŝ with Texp being the time-ordered exponential operator. Here, α n ≡ n|α|n and g n ≡ n|ĝ|n as derived in the previous section. First, consider the ansatz ∆(t) = ∆ 0 [u(t − t on ) − u(t)] with ∆ 0 = ω m − ω q , s.t. the magnons and transmons are resonant for a time-period t on , effectively switching on the coupling. Adding an off-period of t off afterwards, we get t tot = t on + t off and we want to find {t on , t off } to satisfy the above condition. We can derive the time-evolution unitary explicitly using [Ĥ(t),N ] = 0 where the total number of excitations are given bŷ Thus, we divide the wave-function as During the off-time, we have the equations of motion, [using the Schrödinger equation This is a set of linear differential equations which can be solved by diagonalizing the dynamical matrix. For the case N = j, we can write where the vector ξ j = (ξ j0 ξ j1 . . . ξ jj ) T , and the dynamical matrix is diagonalized by {−iΛ off , V off }. Similar analysis with ω q → ω m works for the on-time, We want to ensure Ψ T j, 0 = 0. We findT |j, 0 by putting ξ j (0) = (1 0 . . . 0) T ≡ u j and then, the condition Ψ T j, 0 = 0 becomes where the vector Ψ j = ( j, 0|Ψ j − 1, 1|Ψ . . . 0, j|Ψ ) T . This can be numerically solved for {t on , t off }.
In the special case of a two-level approximation (ignoring higher levels) and for sufficiently large detuning g 0 ω m − ω q , the above condition reduces to In some special cases, Eq. (B.8) will have no solution. Then, we try the ansatz with a pulse of width t (0) on , followed by an off-time t off and a second pulse of an undetermined time t on , off . After another off-period of time t off , the optimization condition becomes (following the same analysis as before) which can be solved for {t on , t off }. This process can be repeated indefinitely, however, in practice, we did not observe the need to go beyond a second iteration.

Appendix C. Transmon Gates
We now describe how the transmon gates can be implemented in our scheme. In what follows, |m, q denotes a state with m magnons and transmon at qth level. Given an initial state |Ψ with a maximum of j magnons and (j + 1) excitations, i.e. m, q|Ψ = 0 for m + q > j + 1 and j + 1, 0|Ψ = 0, and a tunable Hamiltonian α n |n n| +ε(t)q † +ε * (t)q, (C.1) we want to findε(t) s.t. j, 1 Q † Ψ = 0, where the time-evolution operator with Texp being the time ordered exponential operator and t g being a pre-determined gate time. Note that we have ignored the coupling between magnons and transmons assuming large enough detuning ω m − ω q . Then, the evolution is given bŷ Here, α n |n n| +ε(t)q † +ε * (t)q Û g . barring an unimportant global phase Φ global . The parameters for this gate can be found by ensuring j, 1 Q † ideal Ψ = 0, whereQ ideal = e −iωmtgm †m ⊗Û ideal giving and To implement this, we try the ansatzε(t) = −ε(t)e −iωqt+iφapp with real {ε(t), φ app } to be found. We separate out the phases by the transformation, where by definition α 0 = α 1 = 0. First, consider the two-level approximation, where we want to solve d dtÛ where ε (0) (t) is the 'bare' solution to be chosen. The unitary becomeŝ For such a case, we havê This resembles the ideal gateÛ ideal except for the extra phases due to the free evolution of the states. For this two-level approximation, we can ensure (C.14) Notice that we can choose any ε (0) (t) satisfying tg 0 dτ ε (0) (τ ) = θ ideal . For concreteness, we choose a sinusoidal
We try for a pulse shapẽ Again, we have the 'bare' pulse ε (0) which we take as half-sines. To remove the errors at first order in the Magnus expansion, we find that ε(t) here satisfies the same set of conditions Eqs. (C.20)-(C.22). Thus, the same calculations as above apply.