Exploring synchrony and chaos of parahydrogen-pumped two-compartment radio-frequency ampliﬁcation by stimulated emission of radiation

A nuclear-spin-based RASER (radio-frequency ampliﬁcation by stimulated emission of radiation) is an ideal experimental system to explore nonlinear interaction phenomena of nuclear spins coupled via virtual photons to a resonator. This is due to the RASER being stable for several hours, allowing for extended observation of these phenomena. Nonlinear phenomena in multimode RASERs range from mode oscillations in synchrony, frequency shifts, frequency combs, period doublings, and even chaos. These phenomena are observed in a parahydrogen-pumped two-compartment proton RASER. In two independently pumped compartments, the separation in frequency space between the two RASER modes is precisely controlled with a magnetic ﬁeld gradient. By controlling the mode separation, we can select the type of nonlinear phenomena observed. A key ﬁnding is that the ranges of mode separation where chaos and synchrony occur are very close together. The experimental results are supported by numerical simulations, based on two-mode RASER equations. DOI: 10.1103/PhysRevA.108.022806


I. INTRODUCTION
RASER (radio-frequency amplification by stimulated emission of radiation) emission in nuclear spin 1/2 systems requires a population inversion between two Zeeman levels, d 0 = N 2 -N 1 .This population inversion can be generated via hyperpolarization techniques, where a large negative nuclear spin polarization can be created, corresponding to a positive population inversion d 0 .
SABRE is an efficient way to continuously hyperpolarize (or pump) a liquid state 1 H RASER at room temperature and at low and high magnetic field [15,31].In SABRE, p-H 2 and target substrates reversibly coordinate to a metalorganic catalyst with a system-specific contact time.During this time, the singlet state of p-H 2 is transferred and converted via J couplings into a highly polarized state of the target nuclei [15,16].Using this approach, the first parahydrogen-pumped 1 H RASER demonstrated a substantial increase in the precision of J-coupled nuclear magnetic resonance (NMR) spectroscopy [32].Subsequent work on SABRE-and PHIP-pumped 1 H RASERs at high and low magnetic fields have shown strong connections to nonlinear science and have extended the range of possible applications.These include high-resolution NMR spectroscopy, gyroscopes and magnetometry, membrane-pumped RASERs, RASERinduced multinuclear signal enhancement, and RASER MRI [31,[33][34][35][36][37][38][39].
In the present work, we analyze nonlinear phenomena in these systems including synchrony and chaos.These two are fundamental phenomena which play an essential role in many disciplines of science and technology [40][41][42][43][44].In the context of a two-mode RASER, synchrony means that both modes oscillate at a single frequency with a fixed phase relationship between them [44].This results in a single line in the RASER spectrum.For a two-mode RASER, chaotic motion is associated with a continuum of spectral lines, and the temporal characteristics of the RASER signal are very sensitive to the initial conditions [35,44].
Chaos arises in many nonlinearly interacting systems, including nuclear spin 1/2 systems.[23,[45][46][47] Recently, 1 H RASER experiments have shown that chaos and intermittence can occur in a p-H 2 -pumped NMR RASER at 1.45 tesla [35].However, for these experiments, the population inversion d 0 was generated outside of the magnet through a hydrogenation reaction of vinyl acetate and hydroxyethyl propionate.This initial population inversion was very high at d 0 ∼ 10 19 , and consequently the system "passed" through several nonlinear regimes while d 0 was decaying over time.Spectra corresponding to two RASER modes of J-coupled lines were observed for about 1 min.Period doubling and chaos were also observed over short time windows (≈2 s).This motivates the precise study of nonlinear phenomena in RASER systems.Limitations of previous work include the transient (non-steady-state) nature of the acquired signals as well as the fixed frequency difference between the modes given by the J coupling [35].
Here we describe an experimental system designed for precise measurements of chaos, multiple-period doubling, and synchrony in RASER systems.Stable continuous RASER emission is generated using SABRE pumping and proton inductive detection at 7.8 mT at a proton resonance frequency of 333.3 kHz.Under these conditions, crucial experimental parameters for RASER emission remain constant [33,39], as will be discussed in Secs.II and III.A frequency separation between two RASER modes is achieved and controlled by a magnetic field gradient G z over two individually pumped compartments.Both compartments are pumped independently by a continuous parahydrogen (p-H 2 ) gas flow, generating two stable 1 H RASER modes which oscillate at two distinct frequencies.The separation between the two modes can be adjusted to arbitrary values by controlling the strength of a magnetic field gradient G z .Here the nonlinear phenomena are dominated by radiation damping.Distant dipolar field effects are averaged out through motion induced by the p-H 2 flow and primarily lead to an equal frequency shift of both modes (see Appendix D).The resulting stable continuous RASER signals can be observed over periods of minutes and allow for a precise Fourier analysis of the signal to identify all mentioned nonlinear phenomena.

II. EXPERIMENTAL
The experimental realization of the SABRE-pumped twocompartment RASER is illustrated in Fig. 1.All experiments were conducted in an electromagnet with a static magnetic B 0 field ranging from 0 to 20 mT.Four shimming coils (G x , G y , G z , and G 2 z ) were used to attain a magnetic field homogeneity of about 1 ppm over a volume of 0.5 cm 3 .The gradient strength was varied throughout the experiments by changing the current I z .For sensitive NMR detection, a cylindrical coil (inner diameter = 10 mm, height = 10 mm) was connected to an external ferrite LC-resonator with a quality factor of Q ext ≈ 300.The combined input coil with external resonator had a total Q ≈ 100 (EHQE-NMR [48]) at a 1 H detection frequency of 333.3 kHz (B 0 = 7.8 mT).The 1 H RASER signal was fed into a low-noise preamplifier and a custom-made lock-in amplifier.Typical signal acquisition times were 100 s with a data acquisition rate of a few kHz.In these RASER experiments, off-resonance frequencies varied between 40 and 160 Hz.
The liquid sample was contained in a cylindrical dualchamber glass sample cell (inner radius r = 4 mm) with a volume of 0.21 cm 3 per chamber.The two chambers were separated by a glass slide with s = 1 mm thickness (see Fig. 3).The total sample volume was V s = 0.42 cm 3 , and the sensitive volume of the detection coil was ≈0.64 cm 3 , resulting in a filling factor (the fraction of the coil detection volume filled with sample) estimated as η = 0.66.
The required p-H 2 delivery was realized by a custommade glass setup (see Fig. 1) in combination with one glass capillary injecting gas into each chamber filled with the SABRE solution (∼100 µm outer capillary diameter, ∼30 µm inner diameter).The p-H 2 flow rate and pressure were controlled with a needle valve in combination with a flow meter and a pressure gauge.The p-H 2 gas pressure was varied between 2 and 7 bar with nearly equal p-H 2 flow rate in both chambers ranging from 10 to 100 sccm (standard cubic centimeters per minute).The H 2 gas (grade 5; purity >99.999%) was enriched to 92%−93% p-H 2 using a Bruker BPHG 90 p-H 2 generator.The sample solution was based on a composition known to yield high proton hyperpolarization using SABRE [15,17] and prepared and handled under inert gas conditions.A methanol-d 4 solution contained the catalyst precursor ([IrCl(COD)(IMes)]) at 5 mmol/L as well as the target molecule pyrazine at c pyr = 100 mmol/L [49].
Pyrazine was chosen as the target molecule because it only has a single-proton NMR resonance with n pyr = 4 chemically and magnetically equivalent protons.The catalyst was activated in situ with a continuous p-H 2 flow of 40-70 sccm at 4 to 5 bar overpressure for several minutes prior to each experiment.p-H 2 overpressure also helps minimizing solvent evaporation, although the p-H 2 bubbling does cause methanol-d 4 loss over time, limiting the measurement times to about 1.5 h at a p-H 2 flow rate of 70 sccm.The build-up of population inversion d 0 during the activation was monitored by applying multiple 90°excitation pulses until the RASER threshold is reached, which is marked by continuous RASER oscillations.The pumping rate of the SABRE process was measured during the initial startup of the RASER signal.Typical values for in the SABRE experiments were between 0.08 and 0.16 s −1 .At flow rates <70 sccm, the polarization is not high enough to observe period doubling or chaos.T 1 , T * 2 for the pyrazine sample were measured, ranging from 5 to 12 s and from 0.4 to 0.7 s, respectively, depending on the actual concentration of pyrazine and catalyst which both change slowly over time due to methanol evaporation.
Each chamber generates its own mode in the presence of a gradient, which is elaborated in detail in Sec.III.To generate a mode separation ν = ν 2 − ν 1 , the current I z was varied in the range from −20 to +16 mA.For the given C z , this corresponds to a range of gradients −3.84 mG/cm < G z < 3.07 mG/cm.Continuous RASER signals at different mode separations were acquired for 50-100 s at each given gradient G z = C z I z .

III. THEORETICAL MODEL
A two-mode RASER model based on the RASER theory [33,35] and with the additional inclusion of a magnetic field gradient is described here to analyze and understand the experimental results.In Sec.II A, a frequency and phase analysis is presented for two nonlinearly coupled oscillators moving in the one-dimensional angular space (1D model).Next, in Sec.II B, the gradient-controlled two-compartment RASER equations moving in the three-dimensional (d, A, ) space are presented (3D model).An intermediate case in the (A, ) space (2D model) is derived in Appendix B. All three models are compared in Sec.IV.

A. Two nonlinear coupled oscillators in 1D space
In this subsection we present a Fourier and phase analysis of two nonlinear coupled oscillators in one dimension.The 022806-3 model is used to explain several phenomena, such as nonlinear frequency shifts, frequency comb, and synchronous motion (i.e., spectral line collapse).Furthermore, the results are used as a reference to compare with the exact two-compartment RASER model (Sec.III B).The model is represented by two rotating beads [Fig.2(a)] at angular positions φ 1 and φ 2 , which oscillate with angular frequencies ω 1 and ω 2 .The two oscillators are assumed to be frictionless and are coupled by a nonlinear element described by R cos(φ 2 − φ 1 ), where R is the strength of the coupling in units (s −1 ).The governing equations of motion [44] are given in the bottom of Fig. 2(a).
Here we focus on a Fourier and phase analysis of the two nonlinearly coupled oscillators and on the expected dependence of the predicted mode separation ν as a function of the frequency difference ν 2 − ν 1 .It is convenient to introduce the phase difference = φ 2 − φ 1 , which represents the angular variable of motion in one dimension.The equation of motion with respect to is obtained by subtraction of the two equations in the bottom of Fig. 2(a), which results in ˙ = (ω 2 − ω 1 ) − 2R sin . ( Figures 2(c)-2(e) show the numerical evaluation of the projected total signal amplitudes cosφ 2 (t ) + cosφ 1 (t ) (top), with the corresponding Fourier transformed spectrum (middle trace) and the evolution of the two transformed phases φ * 2 (t ) and φ * 1 (t ) (bottom).The star indicates the transformation of the phases φ 2 and φ 1 into a frame which rotates with the angular frequency ω 0 = (ω 2 + ω 1 )/2, which means explicitly φ * 2 = φ 2 − ω 0 t and φ * 1 = φ 1 − ω 0 t.In the first case in Fig. 2(c), where holds, the two modes oscillate nearly independently, and the corresponding spectrum (middle) is characterized by two lines at frequencies ν 2 and ν 1 centered at ν 0 = 0.5(ν 2 + ν 1 ) and separated by ν ≈ ν 2 − ν 1 .The two lines are accompanied on both sides by smaller lines, and two consecutive lines are all separated by ν.The transformed phases (bottom) evolve as a linear function of time, R/π , all side lines of the frequency comb vanish, and the spectrum is characterized by exactly two lines at frequencies ν 2 and ν 1 and separated by ν = ν 2 − ν 1 .In the second case in Fig. 2 the spectrum is a dense frequency comb with several side lines and with a center frequency at ν 0 .The separation between two consecutive lines is ν |ν 2 − ν 1 |.The two phases φ * 2 and φ * 1 are linear functions of time on average, and the slope is smaller than |ω 1,2 − ω 0 |.Each phase is superimposed by periodic wiggles which reflect phase modulations caused by the nonlinear interaction.The observed frequency comb in the spectrum is a result of this phase modulation.In the third case in Fig. 2(e), |ω 2 − ω 1 | < 2R, the signal (top) oscillates like a pure cosine, and the corresponding spectrum has collapsed into one line at the center frequency ν 0 .After a short transient evolution of about 100 ms (see bottom) the two phases become constant but not equal.The two modes oscillate in synchrony at one frequency ν 0 with a fixed phase difference, s = (φ FIG. 3. Cross section of a cylindrical sample cell with two separated compartments in the shape of equal circular segments.Both compartments are exposed to a total magnetic field B z = B 0 + zG z , where B 0 is a homogeneous magnetic field and G z the magnetic field gradient.Both compartments are RASER active since they contain sufficiently negative polarized 1 H spins (d 0 > 10 17 ).Parahydrogen (p-H 2 ) bubbles (red) move the liquid fast enough such that all 1 H spins in each compartment experience an average magnetic field strength given by B av = B 0 ± z c G z .The center of gravity z c is given by Eq. ( 3).The averaging process produces two 1 H Larmor frequencies, ν 2 and ν 1 (black dotted lines on the right) in each compartment.The two red lines on the right correspond to the observed frequencies in the RASER spectrum, which are separated by ν.
For values of ν 2 − ν 1 between the linear and the collapse regime, ν vs ν 2 − ν 1 is a nonlinear function with a corresponding differential slope in absolute units larger than one.
In summary, for the two nonlinear coupled oscillators moving in 1D space, either a frequency comb with an associated nonlinear frequency shift or a synchronized (collapsed) regime is observable.Multiple frequency doubling and chaos could not be observed in the 1D space.This is expected by the Poincaré-Bendixson theorem [44], which predicts the existence of chaos to require at least three dimensions.

B. Gradient-controlled two-compartment RASER model in the (d, A, ) space
We derive here the simplest possible model for a RASER with two independently pumped compartments, which in the presence of a magnetic field gradient G z oscillate at two average angular frequencies ω 1 and ω 2 .The corresponding scheme is shown in Fig. 3.A cross section of a cylindrical sample is shown (bottom left), which is composed of two equal circular segments with radius r and separated by a slide of thickness s.Both compartments contain a liquid with chemically equivalent RASER active 1 H spins, evolving in a magnetic field B(z) = B 0 + z G z (top left).The 1 H angular frequency along the central position at z = 0 is given by ω 0 = γ H B 0 , where the gyromagnetic ratio of the 1 H spins is given by γ H = 2π × 4.257 Hz/mG.The gradient coil produces a field z G z = C z I z z.Without motion of the liquid, the 1 H Larmor frequency would be a linear function of z.We assume that by the pumping mechanism, all voxels of the sample are moving fast enough (a few centimeters per second; see Appendix D) such that the 1 H spins of each compartment experience an average magnetic field strength.This average field, indicated by red arrows in the top of Fig. 3, is given by the magnetic field at the two centers of gravity at ±z c in each compartment, i.e., by Assuming equal sample geometries and pumping conditions, both compartments experience the same average dipolar field B dip = 6.4 nT (see Appendix D).The value for B dip adds to the static field B 0 and generates only a small offset.Although distant dipolar fields play a significant role in several nonlinear NMR and MRI experiments [18,23,29,30,45,46,50,51] they are insignificant for the phenomena as presented here due to motional averaging.The average angular 1 H frequencies ω μ , μ = 1, 2 in both compartments can be expressed as ( The center of gravity z c can be calculated for the known geometry of the sample.Assuming a homogeneous spin density, the number of spins at a slice of thickness dz at position z is S d (z)dz, where the function S d (z) in Fig. 3 describes a circular arc with a maximum value at angle θ 0 at position z = s/2.For z > 0 the center of gravity for the right compartment is given by z c = (1/N ) ∫ r 0 S d (z) z dz, where N = ∫ r 0 S d (z) dz is the total number of spins in the right compartment.By symmetry arguments, the center of gravity in the left compartment is -z c .After a coordinate transformation z = r cos θ the integral over dz can be rewritten as an integral over dθ , and the center of gravity becomes z c = r ∫ θ 0 0 sin 2 θ cosθ dθ/ ∫ θ 0 0 sin 2 θ dθ .The maximum angle θ 0 , as indicated in the bottom of Fig. 3, is connected to the radius r and thickness s by θ 0 = arccos(s/2r).An evaluation of the integrals leads to an analytical expression for the center of gravity As a result of the averaging process, the two compartments of the RASER will oscillate at two frequencies which are indicated by the dashed lines at the bottom right of Fig. 3.According to Eq. (3) and with s = 1 mm, r = 4 mm (θ 0 = 1.4454) the center of gravity is calculated as z c = 0.197 cm.For magnetic field gradients used in our experiments ranging from −3.84 mG/cm < G z < 3.07 mG/cm the range of frequency separations between the two compartments is given by −6.4 Hz < ν 2 − ν 1 < 5.12 Hz.
The simplest possible equation of motion for the twocompartment RASER, which can be compared to experiments of two interacting RASER modes with equal amplitudes, is given by ḋ In Eqs. ( 4)-( 6), d, A, and are the population inversion, transverse spin component, and phase difference = φ 2− φ 1 .Both modes have equal amplitudes A 1 = A 2 = A and population inversions d 1 = d 2 = d (see Appendix A for the derivation).In the following, the model represented by Eqs. ( 4)-( 6) will be termed the 3D model, since the motion is in the three-dimensional (d, A, ) space.1/T * 2 and 1/T 1 are the effective transverse and the longitudinal relaxation rates, respectively.d 0 and are the equilibrium population inversion and the pumping rate, which are assumed to be equal for both modes.The coupling constant between the two modes is given by β = μ 0 hγ 2  H ηQ/(4V s ), where μ 0 is the vacuum permeability, h is Planck's constant, Q is the quality factor of the LC resonator, V s denotes the sample volume, and η = [0, 1] denotes the filling factor of the detection coil.
The parameters used for the simulations were measured or calculated based on the actual experimental conditions.The coupling parameter β is calculated from experimental parameters (Q = 100, V s = 0.42 cm 3 , η = 0.66), i.e., β = 3.7 × 10 −16 s −1 .For T * 2 = 0.4 s, T 1 = 6 s, and = 0.12 s −1 , the RASER threshold is d th = 1.6 × 10 16 (see Appendix A).For each compartment with V s /2 = 2.1 × 10 −4 L the equilibrium population inversion d 0 can be estimated from the 1 H polarization P H ≈ −3 × 10 −3 , the concentration of the pyrazine molecules c pyr ≈ 0.1 mol/L, the number of protons n pyr = 4, and N A = 6.023 × 10 23 /mol.The estimated value is d 0 = c pyr (V s /2) (−P H )n pyr N A = 1.52 × 10 17 .A minus sign in front of the negative-valued polarization P H ensures that the condition d 0 > 0 for RASER action holds.In a RASER experiment with 100 s duration typical fluctuations of and d 0 are on the order of 10%-20%.For the two compartments, and on the timescale of our measurements (∼30 s), this fluctuation in d 0 corresponds to a dipolar-field-induced frequency difference of about ±27 mHz, which is superimposed by frequency drifts of ∼100 mHz induced by B 0 field fluctuations in the field of the electromagnet.
The 3D model Eqs.( 4)-( 6) opens the possibility for chaotic motion, since according to the Poincaré-Bendixson theorem the minimum number of dimensions for chaos is three [44].The chaotic regime might eventually be accompanied by multiple period doubling and by mode collapse.In order to explore experimentally at which frequency differences all these phenomena occur, it is crucial that the mode separation ω 2 − ω 1 = γ H 2z c G z can be controlled precisely while the other parameters T 1 , T * 2 , d 0 , β, and are known and kept constant.4)-( 6), all nonlinear phenomena occur at the much lower values |ν 2 − ν 1 | < 2.5 Hz, and for |ν 2 − ν 1 | > 2.5 Hz a linear dependence ν = ν 2 − ν 1 = (γ H /π )2z c G z is expected.This is experimentally proven by Fig. 4(b), where the measured values for ν (blue circles) are plotted versus the applied G z .An exact proportionality between ν and G z is observed outside the range −1 mG/cm < G z < 0.67 mG/cm.The two slopes with positive and negative sign are determined by linear regressions (solid blue lines), which have the same absolute values |(γ H /π )2z c | = 1.771Hz cm/mG.From this slope, the center of gravity is given as z c = 0.208 cm.This value z c is slightly higher than the theoretically expected value of z c = 0.197 cm, which follows from Eq. ( 3) using θ 0 = arccos(s/2r) with r = 4 mm, s = 1 mm.The difference of about 5% between the measured and the theoretical z c can be explained by the adhesive used to fix the glass slide within the glass cylinder, increasing the average thickness of the slide.The regression functions in Fig. 4(b) intersect the G z axis at an offset gradient G off z = −0.299mG/cm.In summary, as a result from the measurements shown in Fig. 4, the frequency separation is given by ν 2 − ν 1 = 1.771Hz cm/mG (G z + 0.299 mG/cm).
We will now discuss the cases of mode separations of |ν 2 − ν 1 | < 2 Hz. Figure 5 shows five different scenarios, which correspond to (a) frequency comb, (b) period doubling (p-2), (c) twofold period doubling (p-4), (d) chaos, and (e) mode collapse (synchrony).On the left side in blue the measured spectra are shown together with the corresponding RASER signals (insets).The signals are measured over a time period ranging from 20 s to 85 s.A symmetrical Hamming window is applied to the signal prior to the Fourier transformation to suppress sinc wiggles.The chosen signal lengths and corresponding Hamming windows are a good compromise between the lowest possible influence of magnetic field fluctuations and sufficient spectral resolution.The five panels on the right (red) are the results from numerical simulations, which are based on the 3D model described by Eqs. ( 4)- (6).The corresponding initial conditions and the simulation parameters β, , T 1 , T * 2 , d 0 , and ν 2 − ν 1 are given in the caption of Fig. 5.The values of these parameters are derived from measurements or from experimental parameters, as described at the end of Sec.III B. As a key result, there is excellent agreement between all five measured and the simulated nonlinear phenomena, given the common set of parameters and the same preset values for ν 2 − ν 1 .
In Fig. 5(a) on the left, at ν 2 − ν 1 = 1.57Hz, a first glimpse of a frequency comb is visible, showing two major lines separated by ν = 1.49Hz accompanied by two smaller lines separated by the same ν from their larger neighbors.Most features of the spectrum and the signal, such as shape, line amplitudes, and separation ν, agree well with the simulation shown on the right.
In Fig. 5(b) on the left, at ν 2 − ν 1 = 1.26 Hz the spectrum of the measured frequency comb shows a separation of ν = 1.12 Hz between two consecutive lines.Additional lines can be identified exactly at half the distance ν/2 between two consecutive major lines.This is the first sign for a period doubling (p-2) process.The simulated spectrum on the right reflects nearly all the basic features of the measured .12,0.116, 0.16, 0.12, 0.12} s −1 .Initial conditions: φ 1 (0) = φ 2 (0) = 0, A(0) = 10 10 , d (0) = 0.The zero time in the insets is not identical to the initial time at t = 0.For high resolution, and to prevent sinc artefacts in the simulated spectra, 100 s of the simulated signal duration is folded with a Gauss window.For the measurements, this window is typically 20 to 40 s, which is a good compromise between magnetic field fluctuations and spectral resolution.spectrum, except for a slightly larger separation ν = 1.16 Hz and simulated amplitudes of the p-2 lines that differ from the amplitudes in the measured spectrum.
In Fig. 5(c) on the left, at ν 2 − ν 1 = 0.44 Hz the frequency comb now consists of 10 lines, and the separation is ν = 0.61 Hz > |ν 2 − ν 1 | = 0.44 Hz.The same is true for the simulated spectrum on the right.This is in stark contrast to the predictions made by the 1D and 2D models, which predict ν < |ν 2 − ν 1 | in the nonlinear regime.Moreover, a twofold period doubling (p-4) is visible in both the measured and simulated spectrum.In each frequency interval ν between two consecutive major lines, three smaller lines appear, which divide the interval ν into four subintervals of ν/4.
In Fig. 5(d) on the left, at ν 2 − ν 1 = 0.88 Hz, the spectrum as well as the RASER signal contains chaotic features.A first hint for chaos is that in the measurement window of 85 s, the RASER signal envelope for a certain interval of several seconds at t = t 0 never repeats at a later time t = t 1 .The second hint for chaos is referring to certain features of the spectrum.Both the experimental and simulated spectrum show a continuum of lines distributed over a frequency interval of about 4 Hz.The envelope of this spectrum roughly decays from the center at ν 0 = 50 Hz in a hyperbolic manner, i.e., proportional to 1/(ν − ν 0 ).This decay is consistent with a criterion for chaos shown by Haken [52].Clear evidence for chaos based on the exponential divergence of initially close trajectories will be shown at the end of this section.In Fig. 5(e), at ν 2 − ν 1 = 0.17 Hz, the measured and simulated spectra have collapsed into one line, and the corresponding RASER signals are pure sinusoidal waves.The lines of the two compartments oscillate in synchrony at one common central frequency (here ν 0 = 50 Hz) and with a constant phase shift.
Figure 6 is an overview of the measured and simulated mode separations ν vs ν 2 − ν 1 for the 1D, 2D, and 3D models.Figure 6(a) shows the results of three numerical simulations in order to compare the 1D (green line), 2D (brown line), and 3D model (red line).To ensure a valid comparison, the same rates R 1D = 44.for the 1D, 2D, and 3D model are far from each other.A possible explanation is that with increasing number of dimensions, the corresponding two RASER trajectories are less restricted in their configuration space.This reduces their effective interaction, with the consequence that the regime of mode collapse becomes smaller.The enlarged view in Fig. 6(b) shows the measured and simulated windows for period-n, chaos, and synchrony.All measurements (blue circles) associated with frequency combs, period-2, period-4, chaos, and mode collapse are marked by arrows.Three chaotic regimes indicated by 1, 2, and 3 (red areas) left and right from the regime of mode collapse (green area) are predicted by the simulation.The experimentally observed chaos at ν 2 − ν 1 = 0.88 Hz lies within the simulated region 3 of chaos with 0.83 Hz < |ν 2 − ν 1 | < 1.09 Hz.In the 3D (d, A, ) space, chaos is characterized by a never-closing strange attractor [44].Close to each chaotic window, the associated trajectory of period-n circles n times per period in the (A,d) space before closing.Note that both chaotic windows denoted by 1(0.22 Hz |ν 2 − ν 1 | < 0.32 Hz) are lying very close to the regime of synchrony given by |ν 2 − ν 1 | < 0.22 Hz.We remark that a larger number of additional period-n windows (n = 3, 5, 6, 8, 10) have been found in simulations, but not experimentally verified (see Appendix C).
Finally, we present a more concise test that the measured and simulated RASER signals in Figs. 5 and 6 indeed represent a chaotic regime.One intrinsic characteristic of chaotic behavior is the extreme sensitivity towards changes in initial conditions.Two trajectories are chaotic when they start extremely close together and very soon diverge in an exponential manner.This is demonstrated in Fig. 7.In each panel two RASER trajectories (red and green lines) are simulated, which differ only by their initial small phase difference red (0) = φ red 2 (0) − φ red 1 (0) = 0 and gr (0) = φ gr 2 (0) − φ gr 1 (0) = 10 −13 .All other simulation parameters and initial values are equal for both trajectories (see caption).The two RASER trajectories in Figs.7(a) and 7(b) were chosen as a reference in order to compare with the adjacent chaotic regimes.The first two RASER signals in Fig. 7(a) represent a period-4 trajectory, with ν 2 − ν 1 = 1.107Hz located close above the chaotic window 3(0.85 Hz < ν 2 − ν 1 < 1.09 Hz). Figure 7(b) corresponds to a period-2 trajectory simulated at ν 2 − ν 1 = 0.35 Hz, which is close to the chaotic window 1(0.21 Hz < ν 2 − ν 1 < 0.32 Hz).If the trajectories in Figs.7(a) and 7(b) were not chaotic, the signals and the spectra (insets) should be insensitive to a small change in the initial conditions.Indeed, neither the two trajectories (red and green) nor the corresponding two spectra in Figs.7(a) and 7(b) can be distinguished.
This high sensitivity to small changes in the initial conditions is shown in Fig. 7(c).Two RASER signals with a phase difference of 10 −13 are simulated at ν 2 − ν 1 = 0.88 Hz, which is inside the chaotic window 3. The two signals are identical for the first 40 s, but for t > 40 s they become significantly different.The corresponding two spectra in the inset, which result from a Fourier transformation on the window 60 s < t < 260s, are broad and very different.The diverging two trajectories can also be observed in Fig. 7(d) with ν 2 − ν 1 = 0.27 Hz, which lies close above the point of collapse (synchrony) at 0.22 Hz.Compared to Fig. 7(c), the two chaotic trajectories start to diverge at a later time t > 60 s, and the two different chaotic spectra are broader.

V. CONCLUSION AND OUTLOOK
In conclusion, the physics of the two-compartment RASER has been investigated with respect to nonlinear phenomena such as frequency combs, nonlinear frequency shifts, multiple-period doublings, chaos, and synchrony (mode collapse).Only the two-mode RASER theory evolving in the 3D (d, A, ) space is in full agreement with all experimental results.Therefore, the 3D-RASER theory should be able to predict chaos and synchrony in regimes not accessible for current experiments.This includes the physics of N > 2 interacting RASER modes at exceedingly high pumping rates and/or close to the maximum possible polarization One remarkable result is that under appropriate conditions, the range of mode separation for chaos and synchronous motion lie very close together (see Appendix C).Currently, it is not evident what defines the number of chaotic windows at specific parameter values and how many chaotic windows exist for N > 2 modes.These questions will be the focus of future investigations.
We believe that the present experimental and theoretical results are important for several reasons.First, for RASER MRI [39], which has been discovered recently, a deeper understanding of the image formation process is necessary.In RASER MRI, many RASER active slices interact in the presence of a magnetic field gradient, so synchrony and eventually chaos play an important role for understanding image formation.
Second, spectra showing chaos and multiple period doubling were demonstrated at mode separations ν 2 − ν 1 much smaller than the line width ν L = 1/(π T * 2 ) = 0.8 Hz.This might open a new avenue for high-resolution NMR spectroscopy, where chemical shift differences and J-couplinginduced splittings well below the line width ν L can be measured.The prerequisite is a mathematical model, which assigns the measured spectral features (i.e., the location of chaos and period-n signals) to all involved J couplings and chemical shift differences.
Third, calculations based on the EHQE approach [32,48] show that RASERs pumped with p-H 2 , by DNP, or by SEOP can be miniaturized to scales smaller than 100 µm with a moderate loss in signal-to-noise ratio.This approach enables the implementation of microfabricated RASER devices, which can be integrated into electronic circuits, for example, as on-chip local frequency synthesizers or magnetic field and rotational sensors, or as chaos-based true random generators.nonlinear frequency shifts, and synchrony (mode collapse), similar to the 1D model of two coupled oscillators [Eq.( 1)].We are now interested in where the collapse of Eqs.(B2) and (B3) occurs and compare this result with the collapse in the 1D and 3D model.
The condition for collapse for the 1D model is simply proportional to d col , i.e., (ν 2 − ν 1 ) col = R/π = (β/π )d col where (ν 2 − ν 1 ) col is the frequency difference of where the collapse should occur.For the 2D model, the functional dependence between (ν 2 − ν 1 ) col and d col is more complex.At the collapse, two stationary conditions have to be fulfilled in Eqs.(B2) and (B3), i.e., dA s /dt = 0 and d s /dt = 0. From the first condition dA s /dt = 0 the squared stationary amplitude can be calculated as . After replacing A 2 s in Eq. (B3) and setting d s /dt = 0, a transcendent condition for the stationary phase s at a given (ν 2 − ν 1 ) col can be calculated, where in Eq. (B5) ν L = 1(π T * 2 ) is the line width.An elegant expression for d col can now be derived by inserting the known values for A 2 s and s into Eq.(A4).After some algebraic transformations, the result is In Eq. (B6), d th = 1/(βT * 2 ) is the threshold population inversion.The pair of Eqs.(B5) and (B6) constitutes two equations which for a given (ν 2 − ν 1 ) col and line width ν L implicitly fix the stationary phase s and d col .Note that the range of the stationary phase must be within π/2 s π .For 0 s < π/2, which is the possible range for the 1D model, Eq. (B6) predicts d col d th .So for the 2D model the range 0 s < π/2 can be excluded since there is no RASER activity below the threshold.Close to threshold, i.e., close to d th = 1/(β T * 2 ), the stationary value for the phase difference is s = π/2, (ν 2 − ν 1 ) col = ν L and d col = d th .Equation (B5) states that for any given frequency difference (ν 2 − ν 1 ) col > ν L there always exists a corresponding stationary phase s , and the population inversion for collapse d col is a monotonically increasing function of (ν 2 − ν 1 ) col .Finally, for the 3D model we have not found an analytical expression for the point of collapse.Therefore, we have to rely on numerical simulations.at frequency ν * μ deviates from the free 1 H RASER frequency ν μ by [57,58] For our two-compartment RASER the 1 H line width for T * 2 = 0.4 s is ν L = 1/(π T * 2 ) = 0.8 Hz, and the two modes oscillate at frequencies Given the quality factor Q = 100 of the LC resonator resonating at ν c = 333 KHz, so ν c = 3.33 kHz and a typical resonance offset (ν c − ν μ ) ∼ 50 Hz, Eq. (D1) predicts for two modes lying close together (a few Hz) a frequency shift towards ν c of (ν μ − ν * μ ) ∼ 12 mHz.Specifically for two RASER modes differing by ν 2 − ν 1 = 1 Hz the relative frequency shift due to cavity pulling is 0.24 mHz, which is much smaller compared to our observed frequency shifts and line splittings (at least a few tens of mHz).Therefore cavity pulling effects can be neglected in all our experiments.
Distant dipolar fields (DDFs) and their nonlinear effects play a significant role in highly polarized nuclear spin ensembles with high spin densities.One recent example is the DNP-pumped solid state 1 H RASER [51], where strong DDF effects (chirped RASER pulses) and frequency shifts in the kHz regime have been observed at low temperatures (T ∼ 1.2 K).The number of polarized 1 H spins in this experiment were about d 0 > 5 × 10 21 for a sample volume V s ∼ 0.8 cm 3 .In a solid the 1 H spins are not moving on the timescale of the RASER experiments so DDF effects do substantially contribute to the spin evolution in addition to radiation damping effects.A review of nonlinear NMR spectroscopy in liquids including DDF and radiation damping effects has been published by Desvaux [50].
In our model for two coupled RASER modes [Eqs.( 4)-( 6)], the radiation damping field β•d(t) is about 18 s −1 , which is nearly 100-fold larger than the DDFs (0.27 Hz, see below).For each experiment the ratio between radiation damping and DDF effects can be evaluated by inspection of the corresponding time constants τ rd and τ d , respectively.The radiation damping rate is given by 1/τ rd = (μ 0 /2)γ H ηQM 0 , with positive or negative values for 1/τ rd depending on the sign of the 1 H magnetization M 0 = ( h/2)γ H n s P H (n s is the 1 H number density and P H the polarization).We define the rate of evolution due to an average DDF field as 1/τ d = | δω d | = |μ 0 ξγ H M 0 |.In both expressions for 1/τ rd and 1/τ d the vacuum permeability μ 0 , the 1 H gyromagnetic ratio γ H and magnetization M 0 are common.For radiation damping the product between the filling-and the quality factor ηQ is important, while for DDFs the shape factor ξ of the sample matters.This shape factor ranges from −1/2 ξ 1 and is a function of the sample geometry with ξ = 0 representing a sphere, ξ = −1/2 a long cylinder, and ξ = 1 a thin disk.The absolute of the ratio R between the two time constants τ rd and τ d can now be defined as a function only depending on ξ , η, and Q: For R = 0 (ξ = 0, spherical sample) any DDF effects are absent, for R 1 radiation damping dominates the spin evolution, and for R 1 DDFs dominate the spin evolution.An extreme example for DDF effects in the 10 mHz regime including slow diffusional motion in the gas phase is reported in [59], where the chemical shift of hyperpolarized 129 Xe in liquid toluene was measured in the Earth magnetic field.At higher polarization the Xe line of the resting liquid phase shows a substantial line distortion, broadening, and frequency shift with respect to the Xe gas line.Since the experimental conditions in [59] were Q ∼ 10, η = 0.04 (2 cm 3 sample volume, 50 cm 3 coil volume), and |ξ | ∼ 0.8 (slanted cylinder), we get R = 4, so DDFs dominate over radiation damping effects.Therefore distant dipolar fields inside the liquid compartment distort and broaden the line shape.However, the shape of the Xe gas line is identical to the Xe gas peak at much lower Xe polarization (with small DDF and radiation damping effects), so apart from a very small average dipolar shift there is no line distortion, although the Xe atoms in the gas phase experience different DDF fields from the liquid phase at different locations in the sample.This is explained by the Xe gas diffusion (coefficient D Xe ∼ 0.1-0.2cm 2 /s at 1 bar and RT), which average out the local distant dipolar fields' interaction on the timescale of the measurement (tens of seconds).In the resting liquid phase this is not the case since D Xe in liquid toluene is about four orders of magnitude smaller compared to the gas phase, so there is no averaging out of DDF inhomogeneities.
In contrast to the above case in the present two-mode RASER experiments the parameters Q = 100, η = 0.66, and ξ = −1/2 (approximately a long half cylinder for each compartment) lead to R = 0.015 1, and radiation damping dominates the spin.Furthermore, the p-H 2 bubbles lead to a rather intense convective and turbulent motion of all volume elements in both compartments.The average velocity of these bubbles was estimated to be several cm/s, so the fast random motion corresponds to a diffusion coefficient D 1cm 2 /s.Consequently, if the RASER is in a stationary state, all DDF inhomogeneities which lead to a line broadening or to a complex magnetization evolution average out and only an average dipolar field B dip given by remains.Given ξ = −1/2, n s = 2.4 × 10 26 /m 3 for a pyrazine concentration of c pyr = 0.1 mol/L and a 1 H polarization P H ∼ −3 × 10 −3 the average dipolar field in each compartment is B dip ∼ 6.3 nT.This corresponds to an average frequency shift of (γ H /2π ) B dip = 0.27 Hz.This shift is in agreement with the measurement of initial relaxation oscillations of the pyrazine two-compartment RASER (at G z = 0).The measured frequency difference between the oscillations before and after the first maximum of the initial RASER signal (where magnetization M z changes sign) was measured to be between 0.2 Hz and 0.4 Hz depending on the pumping conditions.
Since both compartments have nearly the same shape, volume V s , polarization P H , and spin density n s , both corresponding RASER frequencies are shifted by the same amount of 0.27 Hz, so averaged distant dipolar fields can be neglected to first order.Only fluctuations of the bubbling and pumping rates between the two compartments, which are on the order of 10%, could be responsible to dipolar differential shift effects in the order of ±27 mHz.This small shift could be relevant if the frequency difference between the two modes is at the boundary between mode collapse and chaotic window 1 (ν 2 -ν 1 = 0.21 Hz; see Fig. 6).In this case the RASER signal could change randomly from a collapsed state (one single oscillation) to a chaotic signal.
In conclusion, in the experimentally observed time windows neither cavity pulling nor DDF effects have an appreciable influence on the evolution of the two 1 H RASER modes.All measured different scenarios shown in Fig. 6, the chaotic windows, the mode collapse, and multiple period doublings, can be solely explained by the coupling of the two modes by radiation damping, and thus the transitions between these different scenarios can be predicted by the model given by Eqs. ( 4)- (6).

FIG. 1 .
FIG.1.Experimental setup of the SABRE-pumped twocompartment1 H RASER, adapted from Suefke et al.[32], depicting the parahydrogen (p-H 2 ) supply (red) consisting of the generator, valves (V), a needle valve (NV), pressure gauge (p), and the injection capillaries, as well as the external high-quality enhanced (EHQE) NMR circuit48 (blue) and the magnetic field B 0 with gradient-induced additional term z G z (green).Bottom left: para-hydrogen molecule (p-H 2 ) together with SABRE catalyst IMes-Ir used to create negative polarization P H ≈ −3 × 10 −3 of pyrazine (pz).Bottom right: Example for a 1 H chaotic signal of pyrazine measured over a time period of 20 s.

Figure 2 (
Figure 2(b) shows a plot of the numerically evaluated frequency difference ν as a function of the frequency difference ν 2 − ν 1 .At a fixed value for the coupling strength R = 12.6 s −1 , the transition into the synchronized motion (the collapse into one line with ν = 0) occurs at |ν 2 − ν 1 | = R/π = 4.01 s −1 .Far from this region, i.e., for|ν 2 − ν 1 | R, ν is linear in the argument ν 2 − ν 1 , more specifically | ν/(ν 2 − ν 1 )| = 1.For values of ν 2 − ν 1 between the linear and the collapse regime, ν vs ν 2 − ν 1 is a nonlinear function with a corresponding differential slope in absolute units larger than one.In summary, for the two nonlinear coupled oscillators moving in 1D space, either a frequency comb with an associated nonlinear frequency shift or a synchronized (collapsed) regime is observable.Multiple frequency doubling and chaos could not be observed in the 1D space.This is expected by the Poincaré-Bendixson theorem[44], which predicts the existence of chaos to require at least three dimensions.

FIG. 4 .Figure 4
FIG. 4. Measured separation between two RASER modes ν vs the gradient G z .(a) Top: 1 H RASER signal of SABREpumped pyrazine measured at B 0 = 7.8 mT (333.3 kHz 1 H Larmor frequency) and for G z = 1.152 mG/cm.Bottom: The corresponding Fourier transformed spectrum shows two lines separated by ν = 2.52 Hz.The RASER signal has been folded with a Hamming window in order to avoid sinc artifacts.The off-resonance center frequency is ν 0 = 50 Hz.(b) Measured mode separation ν (circles) vs gradient G z .All data points with |G z | > 0.77 mG/cm fit well to two linear functions (solid lines) with a slope given by |(2π ) −1 γ H 2z c | = 1.771Hz cm/mG.