Tunable Non-Markovianity for Bosonic Quantum Memristors

We studied the tunable control of the non-Markovianity of a bosonic mode due to its coupling to a set of auxiliary qubits, both embedded in a thermal reservoir. Specifically, we considered a single cavity mode coupled to auxiliary qubits described by the Tavis–Cummings model. As a figure of merit, we define the dynamical non-Markovianity as the tendency of a system to return to its initial state, instead of evolving monotonically to its steady state. We studied how this dynamical non-Markovianity can be manipulated in terms of the qubit frequency. We found that the control of the auxiliary systems affects the cavity dynamics as an effective time-dependent decay rate. Finally, we show how this tunable time-dependent decay rate can be tuned to engineer bosonic quantum memristors, involving memory effects that are fundamental for developing neuromorphic quantum technologies.


Introduction
In open quantum systems, the Markovian approximation is widely used due to its mathematical simplicity and because it provides a good description of the phenomenology observed in the lab. The Markovian approximation implies that the state of the reservoir is not correlated at different times, which can be interpreted as a memoryless bath [1,2]. Nonetheless, in several real-world systems and applications, memory effects play an important role and have to be explicitly accounted for, such as in biological systems [3][4][5], quantum metrology [6,7], quantum simulation [8,9], and quantum memdevices [10][11][12][13][14]. Therefore, the manipulation of non-Markovianity is an important step for the development of new technologies, where quantum memristive devices could particularly benefit the realization of neuromorphic quantum computing [15][16][17].
On the other hand, the definition and quantification of non-Markovianity in quantum systems is still an open question [18][19][20][21][22]. Nevertheless, there are two widely accepted cases in the scientific community. One is based on the distinguishability of the states of a quantum system [21] undergoing dissipative evolution. This definition takes into account that when a quantum system interacts with a Markovian environment, the system's information will flow unidirectionally to the environment. States evolving under these conditions can only lose distinguishability during the dynamics. Therefore, we can distinguish non-Markovian behavior if, for a dissipative evolution, the states recover some of their distinguishability during certain time intervals. The second definition of non-Markovianity is based on the behavior of entanglement between the system and an auxiliary system [22]. If at the initial moment of the evolution the system is correlated to some auxiliary system, then these quantum correlations can only monotonically decrease in time if the system is coupled to a Markovian environment. Therefore, we identify non-Markovianity with cases where the quantum correlations do not decrease monotonically with time. Both of these definitions of non-Markovianity provide a method to quantify its degree and involve an optimization process over all possible initial states for the evolution.
In experiments using an all-optical setup, the transitions between Markovian and non-Markovian regimes can be reached, controlling the information backflow of the system [25] as well as the observation of the so-called weak non-markovianity regime [26].
In this article, we focus on dynamical non-Markovianity (DnM), which means the degree of non-Markovianity presented in a given dynamic. Specifically, we will focus on a system composed of a cavity mode (main system) coupled to a set of qubits (auxiliary systems) described by the Tavis-Cummings (TC) model [30][31][32] embedded in a Markovian bath. We are interested in studying the DnM that arises in the main system dynamics by tracing the auxiliary qubits, creating a tunable bosonic quantum memristor. We will also explore how the DnM can be manipulated by external control over the auxiliary systems. We found that by tuning the energy gap of the set of qubits, we can simulate a time-dependent decay rate in the cavity going from a regime with maximal DnM to another with minimal DnM and Markovian evolution. This tunable dynamical non-Markovianity allows us to define variables that follow a memristive behavior, obtaining an experimentally feasible, scalable and general framework to implement switchable memory devices that are useful for neuromorphic quantum computing.

Model and Methods
We consider a system consisting of a single bosonic mode (resonator) coupled to a set of n qubits in contact with a thermal reservoir at zero temperature, as shown in Figure 1. The interaction between the qubits and the resonator is described by the TC model are the Hamiltonians for the bosonic mode, the qubits, and the resonator-qubits interaction, respectively. Here, ω R , ω Q , g, andh denote the resonator frequency, the qubit frequency, the qubit-resonator coupling strength, and the Planck constant, respectively. The operator a(â † ) is the annihilation (creation) operator for the bosonic mode,σ z,j is the Pauli z matrix for the jth qubit, andσ is the lowering (raising) operator for the jth qubits. In order to ensure the validity of our model, we consider ω Q /ω R ∼ 1 and g/ω R < 0.1. From now on we will considerh = 1.
We consider that the total system undergoes Markovian dynamics described by the following master equation,ρ whereÔ 0 = a andÔ j = σ − j for j > 0 and Γ j is the decay rate of the jth channel. We are interested only in the dynamics of the resonator; thus, we focus on its reduced state by tracing out the qubits, ρ R (t) = Tr Q (ρ(t)). In this way, the set of qubits act as an auxiliary system that introduces non-Markovian properties to the dissipative evolution of the resonator. We remark that the TC model has been experimentally realized in several platforms, such as quantum dots [33], trapped ions [34], and superconducting circuits [35]. We highlight that in superconducting circuits, it has been achieved up to 10 qubits coupled to a single resonator, which is even higher than the number of qubits we consider in this work to induce non-Markovian dynamics in the resonator. Figure 1. Diagram of the model: a cavity (bosonic mode) coupled to a set of qubits embedded in a Markovian reservoir. Each auxiliary qubit can be dynamically tuned, and the cavity can be classically driven.

Environment
We want to characterize the degree of non-Markovianity of a particular evolution of our system (the resonator) determined by its initial state. We look for a figure of merit that can be understood as a degree of non-Markovianity of the particular dynamics of the system that result from a given initial condition. To this end, we notice that when the dynamics of the system are Markovian and purely dissipative, then its quantum state will monotonically approach the corresponding steady state of the environment. We can characterize this behavior by calculating the trace distance between the instantaneous state of the system and the steady state of the evolution, where the subindex S denotes that the trace distance is taken with respect to the steady state. For a Markovian evolution, this quantity will decrease monotonically to zero [21], where |ρ| = Tr[ ρ † ρ] and ρ SS is the steady state of the system. In our case, the temperature of the environment is zero, and therefore, ρ SS = |0 0|. Now, the quantity D(t) allows us to detect when the evolution deviates from Markovian behavior whenever it is no longer monotonically decreasing. Therefore, we can characterize the non-Markovianity of a particular system evolution by considering all the time intervals with non-monotonic behavior. In this way, we define the DnM as where ζ(t) = (d/dt)D S (ρ R (t)), and the integration time goes from 0 to infinite, which means until a sufficiently long time as the system reaches the steady state. Notice that for a given interval t ∈ [τ i , τ j ], in which the trace distance detects information backflow (derivative is positive), the DnM over that interval is simply D S (ρ R (τ j )) − D S (ρ R (τ i )), which provides a simple way of calculating N D . In addition, this definition is closely related to the non-Markovianity measure for dissipative channels based on distinguishability [21]. However, our definition considers only the dynamics under study and not an optimization over all the initial conditions. We emphasize that the evolution of the total system is Markovian, and thus, it ultimately reaches the thermal steady state corresponding to the temperature considered. Nonetheless, the reduced dynamics of the resonator will be non-Markovian since the qubit-resonator coupling will introduce an information backflow from the set of qubits. In the next section, we will characterize how a given configuration of the set of qubits affects the behavior of the DnM for such a bosonic quantum memristor.

Dynamical Non-Markovianity (DnM)
For our first case, we will focus on the resonator interacting with one qubit (Jaynes-Cummings model) and interacting with n = 5 qubits. We will analyze how the DnM depends on the qubit frequency and coupling strength. It is important to mention that the set of auxiliary qubits is always initialized in the ground state in order not to introduce energy into the resonator, since it would undermine the interpretation of the DnM. For our calculations, we numerically solve the master equation of Equation (3) in Python using the mesolve method in QuTiP [36]. Then, the DnM is calculated from the reduced state of the resonator. First, we consider the initial state |ψ 0 = |1 R 0 Q . In Figure 2a, we show the DnM of the resonator when varying the coupling strength g/ω R and the frequency ratio ω Q /ω R . We can see that the DnM is largest when the qubit and resonator are in resonance and when g increases. Notice that for larger values of g/ω R , the qubit-resonator detuning can yield significant values of DnM. Figure 2b shows the DnM for the case of five auxiliary qubits. We can observe that the effect of enlarging the set of auxiliary qubits is relaxing the resonance condition and increasing the value for the DnM. In both cases, the decay rate of qubit and resonator is Γ Q = Γ R = 0.005. We consider resonator frequency ω R = 1, qubit frequency ω Q /ω R ∈ [0.5, 1.5], the coupling strength g/ω R ∈ [0, 0.1] and the initial state |ψ 0 = |1 R 0 Q . This behavior is to be expected since the resonance condition allows for maximal information transfer and information backflow due to the complete Rabi oscillations (in the case of n = 1). In addition, the coupling strength g/ω R is related to the speed of the information transfer and the information backflow. Then, for a small g (slow information transfer), a stronger resonance condition is needed to have information backflow before the system reaches the stationary state. If g increases, the communication between the auxiliary qubits and the bosonic mode is faster, and a more relaxed resonance condition will still have information backflow. Increasing the number of qubits increases the channels for information backflow, which leads to larger values of DnM even at higher detuning.
Next, we study the scaling of the DnM with the number of qubits under fixed conditions. In Figure 3, we study N D as a function of the number of qubits (until n = 8) for different coupling strengths when the resonator is initialized with one excitation. Figure 3a shows how N D scales with the number of particles (n), with a monotonically increasing behavior reminiscent of a power law. In Figure 3b, we do a log-log plot of the quantities of Figure 3a, which confirms the power-law dependence. We perform a linear fit of the data in Figure 3b and obtain in all instances an R 2 coefficient larger than 0.995, which allows us to approximate the scaling of the DnM by where the exponent k depends on the coupling strength g/ω R and the qubit frequency ω Q /ω R . In Figure 3c,d, we show the dependence of the exponent k on the coupling strength and the qubit frequency, respectively, which has been numerically calculated by the same fitting procedure as with Figure 3b. We can see that when the DnM is maximal, that is, for a large g and ω Q = ω R , the value of k is at a minimum. The monotonically increasing behavior of the DnM obtained in Equation (7) is due to the frequency of the Rabi oscillations between the cavity and the set of qubits, which increases with the number of qubits and is caused by the collective enhancement of the coupling strength of the Tavis-Cummings model with n qubits [30]. This phenomenon behaves effectively as a cavity interacting with a single qubit with coupling strength √ ng. Therefore, with faster Rabi oscillations, there are more time intervals where the trace distance is positive from the information backflow and yields a higher value of DnM. Increasing the number of qubits also reduces the relaxation time, and the initial excitation is dissipated faster, but with the parameters we have considered, increasing the number of qubits increases the Rabi oscillations faster than the relaxation. Thus, we obtain the behavior of Equation (7). We need to mention that as the Hilbert space of our system grows exponentially with the number of qubits, we can only calculate Figure 3a until eight qubits. Nevertheless, it is interesting in that it provides more numerical evidence of the power-law dependence given in Equation (7). In addition, an analysis of the thermodynamic limit of the proof of Equation (7) for a large qubit limit would be interesting but is out of the scope of the present study We have seen that the DnM of the resonator strongly depends on the parameters of the set of auxiliary qubits. It is then interesting to consider whether we can have dynamic control over the DnM by manipulating the set of auxiliary qubits. In what follows, we apply a driving term in the z-direction to the set of auxiliary qubits in order to dynamically modulate the qubit gap and control the degree of DnM in the evolution. The driving is chosen so that it does not introduce energy into the qubits, which could excite the resonator and be interpreted as information backflow by the DnM. This situation is described by the following Hamiltonian:Ĥ where H TC is the Hamiltonian of Equation (1), and Ω Q and µ Q are the amplitude and frequency of the driving over the qubits, respectively. Notice that we consider that each qubit is driven by the same signal. We numerically calculate N D for different values of the driving frequency and amplitude, which we show in Figure 4. In Figure 4a, we show the case of one auxiliary qubit. Here, there is a non-zero DnM over the whole range of parameters. However, it is interesting to notice the dark lines that are spanned from near the origin, where the DnM is almost completely suppressed. A similar behavior occurs when we increase the number of qubits, as is shown in Figure 4b for the five qubits case, where the DnM is suppressed over thin lines in the frequency/amplitude plane. Although the suppression is not as strong as in the one-qubit case, these lines show a significant decrease in the DnM. This indicates that by modulating either the frequency or the amplitude of the driving, we can enhance or suppress the DnM of the resonator.
Such suppression of the DnM could be explained by the phenomena of coherent destruction of tunneling [37]. It can be understood as the suppression of coherent evolution between two states of a system due to a coherent driving. We note that we are considering a TC model with at most one excitation evolving only under an amplitude damping channel. In such a situation, our model can be described as a three-level system, where we only populate the states |1 R |0 Q , |0 R |1 Q and |0 R |0 Q , where the state |1 Q is the uniform superposition of one excitation in the set of qubits, and |0 Q is the ground state for all the qubits. It has been shown that coherent destruction of tunneling is presented in a three-level system [38], as well as for the Jaynes-Cummings model [39], which suggests that for a particular value of frequency µ q and amplitude Ω q , the coherent information transfer between the set of qubits and the cavity is suppressed, and therefore, the DnM goes to zero.
Similarly, we study the DnM in terms of the coupling strength for different numbers of qubits in the auxiliary system. Here, we will consider the driving parameters that yield the maximum and minimum DnM obtained from Figure 4 and the analogous calculations for n = 2, 3, 4, 5. In Figure 5a, we plot the minimal DnM for different coupling strengths g. We can see that up to g = 0.02, we can essentially completely suppress the non-Markovian behavior by a suitable choice of driving parameters. Increasing the number of qubits decreases the necessary value of the coupling strength that allows for a completely suppressed DnM. On the other hand, in Figure 4b, we plot the maximal DnM for different coupling strengths g. Here, we can see that N D has a linear dependence on the coupling strength, except for a small range around zero. From these results, we have that, provided we choose a suitable value of the coupling strength, we can switch between Markovian and non-Markovian dynamics for the resonator just by controlling the auxiliary set of qubits. In both cases, the decaying rate of qubit and resonator is Γ Q = Γ R = 0.005ω R . The driving frequency of qubit µ Q /ω R ∈ [0, 1], the driving amplitude of the qubit Ω Q /ω R ∈ [0, 1], the qubit frequency ω Q /ω R = 1, and the coupling strength g/ω R = 0.1. The initial state is |ψ 0 = |1 R 0 Q . Figure 5. The minimum (a) and maximum (b) of DnM for the resonator with a different number of auxiliary qubits. In both cases, the decaying rate is Γ Q = Γ R = 0.005ω R . The frequency of qubit ω Q /ω R = 1, and the initial state is of resonator |ψ 0 = |1 R . The qubits are all initialized in the ground state.
In Figure 6, we show how we can dynamically switch the non-Markovian behavior on and off just by changing the driving frequency of the qubits. Here, we plot the trace distance as a function of time. At the start of the evolution, we choose a driving frequency that yields maximum non-Markovianity (µ Q = ω R ); later, at t = 350ω −1 R , we switch the driving frequency to µ Q = 0.75ω R , which yields the minimum DnM. As can be seen in the figure, at t = 350ω −1 R , the trace distance switches from non-monotonic to monotonically decreasing behavior, which characterizes Markovian evolution. Figure 6. Transition from non-Markovian to Markovian dynamics by changing the driving frequency over the auxiliary qubits. Parameters: the driving amplitude of qubit is Ω q = 0.5, the number qubit n = 1, decaying rate is Γ Q = Γ R = 0.005, frequency ω Q = ω R , and coupling strength g = 0.05.
Finally, we analyze the effect of the decay rates on the dynamical non-Markovianity. We consider varying the decay rates in three cases: varying it only for the qubits, only for the resonator, or varying both simultaneously. We plot the three cases in Figure 7; in general, we can see that smaller decay rates will yield larger values of dynamical non-Markovianity since the original excitation will survive longer in the system before being dissipated. Fixing all the other parameters, in all three cases considered, the DnM appears to exponentially decay as the decay rate increases, with the strongest effect being the case where the qubits and the cavity decay rate are changed simultaneously. Figure 7. Dynamical non-Markovainity of the resonator varies with decaying rate. In the blue line, the decay rate of resonator Γ R = 0.005, and qubit decay rate Γ Q = Γ. In the green line, the decay rate of qubit Γ Q = 0.005, and resonator decay rate Γ R = Γ. In the orange line, the decay rate of qubit and resonator are the same, Γ Q = Γ R = Γ. We consider that qubit and resonator are in resonance (ω Q = ω R = 1), the coupling strength g/ω R = 0.05, and the initial state |ψ 0 = |1 R 0 Q .

Time-Dependent Decay Rate
The observed memory effect can be understood as the system effectively interacting with an environment with a time-dependent decay rate, which becomes negative during some time intervals, favoring the information backflow [40]. To understand this statement, consider a resonator with stateρ undergoing dissipative dynamics with a time-dependent decay rate and without any interaction with an auxiliary system; the system is then described by the following master equation:ρ where H =hω R a † a, and Γ(t) is a time-dependent decay rate that can be negative. Here, ρ represents the state of the resonator undergoing dynamics as described above and is different from ρ R , which is the reduced state of the resonator as described by Hamiltonian (1). For Γ(t) > 0, the energy of the resonator dissipates to the environment, meaning that the information in the resonator is continuously lost. Meanwhile, for Γ(t) < 0, there is energy entering the resonator, giving way to information backflow and therefore to a non-Markovian process. We consider the time-dependent decay rate parametrized as Γ 1 (t) = A(sin(Bt) + C). Notice that the master equation in Equation (9) has the same steady state as that of our original system in Equation (3). Therefore, for a given dynamic induced by the set of auxiliary qubits, we can find the closest non-Markovian dynamic corresponding to a negative decay rate by finding A, B, and C that optimize the cost function |D(ρ R (t)) − D(ρ(t))| 2 , whereρ(t) is the density operator calculated using the time-dependent decay rate. In Figure 8, we plot D S (ρ R (t)) and D S (ρ opt (t)), where ρ R (t) is for one qubit case and ρ opt (t) is the resonator evolved with the optimal parameters for the decay rate. Bottom (d), the DnM of resonator in different driving frequencies µ Q ∈ (0, 1). Parameters: the number qubit n = 1, decaying rate is Γ Q = Γ R = 0.005, frequency ω Q = ω R , coupling strength g = 0.05, and driving amplitude Ω Q = 0.5ω R .
We consider three cases:in Figure 8a, the effective decay rate is time-independent Γ = 0.005, corresponding to Markovian behavior, whereas for the time-dependent decay rate we have in Figure 8b is Γ(t) = 0.05[sin(0.023t) + 0.09] and in Figure 8c it is Γ(t) = 0.25[sin(0.079t)+0.021]. Finally, Figure 8d shows the DnM as a function of the qubit-driving frequency, where it displays the qubit-driving frequency corresponding to Figure 8a-c. We can see that for time-dependent decay, the behavior of both trace distance is almost the same, which means that the set of auxiliary qubits can switch on/off highly non-Markovian dynamics.
Finally, it is interesting to study how this simulated time-dependent decay rate can affect the response of the cavity over external driving, in order to control the memristive properties of the dynamics.

Bosonic Quantum Memristor
One interesting application of our results is to induce memristive behavior into the bosonic mode, which can be tuned by the set of auxiliary qubits. In Ref. [41], it was shown that a certain kind of time-dependent decay rate produces a quantum memristor, which could be reached in a superconducting circuits platform. Later, in Ref. [42], a memristive dynamic was obtained in a quantum computer by the simulation of a non-Markovian bath.
In this line, we analyze the response of the cavity under an external driving, obtaining a Hamiltonian of the form Now, we define the variables I = − i(a − a † ) and O = Ṅ + α N with α a constant. If we consider α = Γ c as the natural decay rate of the cavity, we have that O = F(t)I + G(t), (11) and for more details see Appendix A. The function G(t) depends on the DnM of the system, which means that we can control the memristive relation by controlling the value of G(t). Now, if we choose F(t) = Ω c [1 − sin(cos µ c t)], it is possible to obtain the typical pinched hysteresis loop that characterizes a quantum memristor. This situation is shown in Figure 9, where we obtain the pinched hysteresis loop (green curve) in a similar way as it has been obtained for previous proposals of quantum memristos [12,14,41] as a signature of memristive behavior. It is interesting to note that such bosonic quantum memristor dynamics appear when the auxiliary qubits are not driven and off-resonant with the cavity, which means that in an effective way, the qubits are decoupled from the cavity. In contrast, we can observe that when we drive the qubits, the memristive behavior can be destroyed for different cases, obtaining a way to go from memristive dynamics to non-memristive dynamics. It means that we also can control the memory properties induced by the decay rate in the cavity, which can be helpful in neuromorphic computing, considering that the proposed system can be implemented in many platforms such as optical devices, trapped ions, and superconducting circuits, among others. We also need to remark that our proposal can work as a switchable bosonic quantum memristor. This suggests that our formalism allows for implementing devices with controllable and switchable memory properties only by tuning the energy gap of auxiliary qubits. This proposal opens the door for the experimental implementation of memristive devices, providing a general, platform-free, and scalable model for the next generation of neuromorphic quantum computing technology. In the present study, we have shown results for the homogeneous TC model, where all auxiliary qubits are identical in frequency, dissipation, and coupling with the cavity. The extension of our study to an inhomogeneous system is an interesting possibility for future research, which may focus on controlling the response of the cavity, opening the door to designing memristive devices with a given response. Nevertheless, the extension of our study to such disordered systems is computationally demanding and is left for future studies.

Conclusions
We have considered a cavity coupled to a set of auxiliary qubits, which induce a controllable dynamical non-Markovianity (DnM). We show that by the dynamical tuning of the energy gap of the auxiliary qubits, we can go from high to low values of DnM. The dynamics induced by the auxiliary qubits can be considered as an effective timedependent and tunable decay rate. We also showed that the induced DnM in the cavity mode follows a power-law dependence with the number of auxiliary qubits, at least for a low number of qubits. Finally, we showed as an application that we can define memristive response in the cavity mode, which can be switched off by controlling the energy gap of the auxiliary qubits. This means that we can control the dynamical response of the cavity by external control of the auxiliary system, obtaining a switchable bosonic quantum memristor. These results provide a general protocol to obtain controllable bosonic quantum memristors, which can be useful in neuromorphic quantum computing. As our proposal only considers a Tavis-Cummings model where qubits are considered with a tunable energy gap, our work is experimentally feasible in different platforms, such as trapped ions and a superconducting circuit. Finally, the present work opens the door for the implementation of quantum memristors with minimal experimental requirements, paving the way for implementable neuromorphic quantum computing in the near future. Data Availability Statement: The data in this article are available upon reasonable request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: DnM Dynamical non-Markovianity