We present a comprehensive analysis of the periodic flares observed in the 6.7 GHz methanol transition in G22.3560.066 utilizing the MaxwellBloch equations MBEs as a framework to model these phenomena. By solving the onedimensional MBEs we describe the behavior of both the quasisteadystate maser and transient superradiance regimes. Our findings indicate that the observed periodic flares with varying timescales across different velocities are consistent with the characteristics of Dickes superradiance triggered by a common radiative pump in regions of varying inverted column densities. This work provides new insights into the physical processes governing variability in maserhosting regions and underscores the significance of superradiance as a powerful radiation mechanism in astrophysical environments.
['Astrophysics - High Energy Astrophysical Phenomena']
Superradiance and Periodic 6.7 GHz Methanol Flaring in G22.3560.066
{'No Header': ', 1-5 (0000)', 'T. Rashidi, 1 V. Anari, 1 A. Bartkiewicz, 2 P. Wolak, 2 M. Szymczak, 2 F. Rajabi, 1 ★': '1 Department of Physics and Astronomy, McMaster University, 1280 Main Street West Hamilton, Ontario, L8S 4M1, Canada  \n2 Institute of Astronomy, Faculty of Physics, Astronomy and Informatics, Nicolaus Copernicus University, Grudziadzka 5, 87-100 Toruń, Poland', 'ABSTRACT': "We present a comprehensive analysis of the periodic flares observed in the 6.7 GHz methanol transition in G22.356+0.066, utilizing the Maxwell-Bloch equations (MBEs) as a framework to model these phenomena. By solving the one-dimensional MBEs, we describe the behavior of both the quasi-steady-state maser and transient superradiance regimes. Our findings indicate that the observed periodic flares, with varying timescales across different velocities, are consistent with the characteristics of Dicke's superradiance, triggered by a common radiative pump in regions of varying inverted column densities. This work provides new insights into the physical processes governing variability in maser-hosting regions and underscores the significance of superradiance as a powerful radiation mechanism in astrophysical environments.  \nKey words: masers - radiation mechanisms: non-thermal - ISM: individual objects: G22.356+0.066", '1 INTRODUCTION': "G22.356+0.066 (also known as G22.357+0.066 and IRAS182900924) is a high-mass young stellar object (HMYSO) located at a distance of 4 . 3 + 3 . 8 -1 . 3 kpc (Reid et al. 2019). It is part of the Giant Molecular Complex G23.3-0.3, which comprises several HII regions and supernova remnants, serving as star-forming sites (Zhang et al. 2023).  \nHMYSOs are known to host variable methanol masers (Goedhart et al. 2004; van der Walt et al. 2009; Szymczak et al. 2012). In 2011, Szymczak et al. (2011) reported the discovery of periodic flares of the 6.7 GHz methanol in G22.356+0.066 using the Toruń 32-m radio telescope. The flares appeared with a period of approximately 178 days (178.2 ± 1.9), with a rise-to-decay ratio of 0.34 (Szymczak et al. 2015).  \nThe maser activity in G22.356+0.066 is not limited to 6.7 GHz methanol and is also reported in some velocity components of the 12.2 GHz methanol (Breen et al. 2016), 22 GHz water (Breen & Ellingsen 2011), and 1665 MHz OH (Beuther et al. 2019) transitions. However, the flux density in these lines is either faint or the available data are too sparse, with no reports of periodicity.  \nA number of scenarios have been proposed to explain periodic flaring activities in maser-hosting regions, most of which are based on the maser regime of radiation. Some of these scenarios are based on the disturbance of the masing region by shock waves or clumps. Others are based on variations in the pump photons and/or  \nthe background radiation field. However, when the flux densities return to nearly the same quiescent level after flaring, this indicates that the masing region was not significantly affected by whatever mechanism caused the flare. Therefore, the disturbance scenarios involving shock waves or clumps are dismissed in such incidents (van der Walt et al. 2009) and only mechanisms involving radiative coupling with the masing region are considered viable. Examples of such mechanisms include the colliding wind binary (van der Walt et al. 2009; van der Walt 2011), cyclic accretion onto a young binary system (Araya et al. 2010), and stellar pulsation (Inayoshi et al. 2013).  \nMore recently, a new regime of radiation based on Dicke's superradiance has been introduced in maser-hosting regions to model periodic (Rajabi et al. 2023; Houde et al. 2024) and non-periodic (Rajabi & Houde 2016, 2017; Rajabi et al. 2019) flaring events. In essence, superradiance is a transient, highly efficient radiation mechanism that manifests itself through sharp rises in flux densities, with burst durations scaling with some fundamental properties of the transition, such as the wavelength of the transition, and physical parameters characterizing the source, such as the inverted column density. In Rajabi et al. (2023), superradiance provides a natural explanation for periodic flaring in G9.62+0.20E at four different transitions (6.7 GHz and 12.2 GHz methanol lines, as well as OH 1665 MHz and 1667 MHz) exhibiting different timescales in response to a common radiative pump. In Szymczak et al. (2011), it is noted that the profiles of periodic methanol flares in G9.62+0.20E and G22.356+0.066 are very similar. This suggests that the same underlying physical mechanism could be responsible for the periodicity  \nin both sources. It is, therefore, a natural step to extend the superradiance analysis performed on G9.62+0.20E to G22.356+0.066. However, in G22.356+0.066, the periodicity is only reported in the 6.7 GHz methanol transition line, which will therefore be the focus of this paper.  \nThis paper is structured as follows: In Sec. 2, we briefly discuss the two distinct regimes of radiation, namely the maser action and superradiance. Sec. 3 provides a detailed discussion of the modelling and analysis. In Sec. 4, we evaluate the validity of the superradiance model presented and conclude the paper. The MBEs used in our analysis are detailed in Appendix A. Finally, the autocorrelation analysis developed to quantify and identify variations in the timescales of flares across different velocity components is summarized in Appendix B.", '2 THE MASER AND SUPERRADIANCE REGIMES': "As discussed extensively in Rajabi & Houde (2020) and Rajabi et al. (2023), the Maxwell-Bloch equations (MBEs) (see Appendix A) provide a comprehensive description of a compound system comprising an extended group of 𝑁 two-level molecules (or atoms) and a radiation field. This set of equations tracks the underlying difference in population between the energy levels, the induced polarization in the system, and the evolution of the radiation field (MacGillivray & Feld 1976; Gross & Haroche 1982; Benedict et al. 1996). A detailed study of the MBEs for a gas with an initial population inversion reveals two complementary regimes: the quasi-steady-state maser regime, where the stimulated emission process dominates, and the fast transient superradiance regime (Feld & MacGillivray 1980; Rajabi & Houde 2020; Rajabi et al. 2023; Houde et al. 2024).  \nIn the quasi-steady-state maser regime, temporal variations in the population inversion, polarization, and radiation field occur on an evolution timescale ( 𝑇 e) longer than the noncoherent decaying and dephasing timescales ( 𝑇 1 and 𝑇 2 ) that characterize the system's response. Consequently, a maser system in this regime can quickly adjust to the excitation signal. For example, a saturated maser's response to an inversion pump closely tracks the profile of the pump itself (Rajabi et al. 2023).  \nIn the fast transient superradiance regime, the system responds differently to variations in the pump signal whenever 𝑇 e ≈ 𝑇 R < 𝑇 1 , 𝑇 2 , with 𝑇 R being the characteristic superradiance timescale. The transient intensity triggered by the excitation is thus a characteristic response of the system and is largely independent of the shape of the input signal, whether it is a pump or a trigger (Rajabi et al. 2023). Another distinctive aspect of the superradiance is that its peak intensity scales with the square of the number of molecules involved ( 𝑁 2 ), in contrast to the linear scaling ( 𝑁 ) observed in a saturated maser. The characteristics of the transient response, such as its duration, profile, and strength, are influenced by the physical conditions of the environment, including the relative significance of 𝑇 1 and 𝑇 2 compared to 𝑇 R , as well as the inverted column density. Additionally, the parameters that define the radiative transition, such as the wavelength and the Einstein spontaneous emission coefficient, play a crucial role. More precisely, the characteristic superradiance timescale is defined as:  \n𝑇 R = 8 𝜋 3 𝑛𝐿𝜆 2 𝜏 sp , (1)  \nwhere 𝜏 sp = Γ -1 is the spontaneous emission timescale of the transition, 𝜆 is the wavelength of radiation, 𝑛 is the population  \ninversion density, and 𝐿 is the length of the system ( 𝑛𝐿 is the inverted column density) 1 . It should be noted that while 𝑇 R can be interpreted simply in single superradiance burst events, it is more complex to define in the case of periodic or multi-flaring events. We thus define an initial 𝑇 R , 0 before the onset of the first modelled flare, using the initial inverted column density to establish a relevant timescale for interpreting the MBEs solutions (see Sec. 5 of Rajabi et al. 2023).  \nAs we will show in Sec. 3.2, the observational data of the periodic flaring in G22.356+0.066 support the superradiance model. Weobserve that the speed at which a superradiance system responds to an excitation varies among different velocity components. This is quantified through an autocorrelation analysis (see Appendix B), which reveals the presence of different timescales at various velocities. This is a phenomenon that can be readily explained within the context of superradiance (Rajabi & Houde 2017; Rajabi et al. 2023).  \nMore precisely, this implies that, when a medium is pumped, different velocity components of the 6.7 GHz methanol transition might couple differently to the same pumping source resulting in an inversion density 𝑛 that can vary throughout the medium.", '3.1 Data': 'The data sets used in our analysis were obtained with the Toruń 32-m telescope between June 2009 and April 2014 at irregular intervals; see Szymczak et al. (2015) for more details. The observations revealed cyclic variations in flux density for features within the velocity range of 78.95 km s -1 to 83.50 km s -1 , which are the focus of this investigation. The periodic features in G22.356+0.066 consist of a non-zero quiescent level, upon which the periodic flares are superimposed. In Appendix B, we selected only those features with a nearly constant quiescent level for modelling.', '3.2 Model and results': 'To model the flares in G22.356+0.066, similar to Rajabi et al. (2019) and Rajabi et al. (2023), we solve the MBEs using the fourth-order Runge-Kutta method for a one-dimensional sample with an inverted column density 𝑛𝐿 (Mathews 2017). To simulate periodic flaring, we employ a periodic pump:  \nΛ 𝑁 ( 𝑧, 𝜏 ) = Λ 0 + ∞ ∑︁ 𝑚 = 0 Λ 1 ,𝑚 cosh 2 GLYPH<2> ( 𝜏 -𝜏 0 -𝑚𝜏 1 ) / 𝑇 𝑝 GLYPH<3> . (2)  \nThis pump propagates along the symmetry axis ( 𝑧 -axis) of the sample. It consists of a constant term Λ 0 , modelling quiescent flux levels, and a sequence of pump pulses with amplitudes Λ 1 ,𝑚 , which can vary from pulse to pulse. The period of the pulses, 𝜏 1 , is set to 178.5 days, closely matching the reported period of the flares at the source. The duration of the pump pulses is determined by 𝑇 𝑝 , and a time delay 𝜏 0 is chosen to fit superradiance to the data.  \nThe results of our simulations for four selected velocities are presented in Figures 1-4. The black dots in the top panels represent the observational data, while the solid blue curves correspond to the flux density obtained by solving the MBEs. The bottom panels  \nFigure 1. Top: Periodic flares of the 6.7 GHz methanol transition line in G22.356+0.066 for 𝑣 lsr = 79 . 30 km s -1 (black dots) and the flux density yielded from MBEs (solid blue curve) are plotted as a function of time in Julian Day-2450000. Bottom: The selected periodic pump with Λ 0 = 3 . 97 × 10 -19 cm -3 s -1 and Λ 𝑖 = 4 . 45 × 10 -18 cm -3 s -1 , along with pump amplitudes 𝑎 𝑚 of 0.84, 0.73, 1.10, 0.91, and 0.76 (green), respectively, and the population inversion density ( 𝑛 ) response to the pump (vermilion). The model yields an initial inverted column density of 𝑛 0 𝐿 = 5 . 28 × 10 3 cm -2 , corresponding to 𝑇 R , 0 = 5 . 03 × 10 4 s.  \n<!-- image -->  \ndisplay the population inversion (in vermilion) and the periodic pump (in green), with the horizontal axis showing time in Julian date.  \nIn these figures, after each pump pulse, the population inversion increases, followed by a decay, with the peaks of the population inversion occurring with a delay relative to the peaks of the inversion pump. It is also seen that the peaks of the inversion pump are narrower than those of the population inversion.  \nTo fit the superradiance models to the observed data, we varied Λ 0 and Λ 1 ,𝑚 in equation (2). We kept the parameters 𝑇 1 = 1 . 77 × 10 7 s ( ≃ 205 days), 𝑇 2 = 1 . 02 × 10 6 s ( ≃ 11 . 8 days), 𝜏 1 = 1 . 52 × 10 7 s ( ≃ 178 . 5 days), and 𝑇 𝑝 = 3 . 89 × 10 5 s ( ≃ 4 . 5 days) constant across all velocities. To determine Λ 1 ,𝑚 , we used a constant Λ 𝑖 for each velocity and introduced a varying amplitude 0 < 𝑎 𝑚 ≤ 1 for each flare, defining Λ 1 ,𝑚 = 𝑎 𝑚 Λ 𝑖 .  \nDue to the non-linear nature of the MBEs and the fact that the response of the superradiant system does not precisely follow the pump (see Sec. 2), the choice of parameters Λ 0 , Λ 𝑖 , and 𝑎 𝑚 is not unique. Specific values for these parameters are provided in the figure captions. For each velocity component, the initial inversion density 𝑛 0 is calculated as 𝑛 0 = Λ 0 𝑇 1 , which in all cases exceeds the threshold needed for superradiance and the ensuing saturated maser regime (Rajabi & Houde 2020; Rajabi et al. 2023). Thus each model yields an initial inverted column density 𝑛 0 𝐿 , corresponding to an initial characteristic superradiance timescale 𝑇 R , 0 , as defined by equation (1). This timescale is determined before the onset of the first superradiance burst, when only very weak levels of radiation are present in the simulations. It is important to note that 𝑇 R , 0 ≪ 𝑇 1 , 𝑇 2 in all cases, suggesting the potential for transient superradiance to follow.  \nAs discussed in Appendix B, we used an autocorrelation analysis to compare the flare duration for different velocities. The width of the central peak in the autocorrelation of a periodic function is a measure of the average widths of the flares. These widths,  \nN  \nFigure 2. Top: Same as Figure 1, but for 𝑣 lsr = 80 . 35 km s -1 . Bottom: Same as Figure 1, with Λ 0 = 4 . 53 × 10 -19 cm -3 s -1 and Λ 𝑖 = 2 . 93 × 10 -18 cm -3 s -1 , and pump amplitudes 𝑎 𝑚 of 0.84, 0.81, 0.97, 1.00, and 0.79, respectively. The model yields an initial inverted column density of 𝑛 0 𝐿 = 6 . 02 × 10 3 cm -2 , corresponding to 𝑇 R , 0 = 4 . 41 × 10 4 s.  \n<!-- image -->  \nFigure 3. Top: Same as Figure 1, but for 𝑣 lsr = 81 . 32 km s -1 . Bottom: Same as Figure 1, with Λ 0 = 3 . 51 × 10 -19 cm -3 s -1 and Λ 𝑖 = 3 . 62 × 10 -18 cm -3 s -1 , and pump amplitudes 𝑎 𝑚 of 0.97, 0.92, 0.96, 0.95, and 0.82, respectively. The model yields 𝑛 0 𝐿 = 4 . 66 × 10 3 cm -2 , corresponding to 𝑇 R , 0 = 5 . 69 × 10 4 s.  \n<!-- image -->  \nFWHM A , are compared in Figure 5, where the vertical axis represents the FWHM A in days and the horizontal axis shows velocity in km s -1 .', '4 DISCUSSION': 'Our analysis demonstrates that our model based on the MBEs effectively captures and reproduces the duration of flare profiles of the 6.7 GHz methanol maser at four distinct velocity components in G22.356+0.066. A key feature of the model is its ability to use a common pump across different velocity components, with only the pump amplitudes being modified. Specifically, while the pump-  \nFigure 4. Top: Same as Figure 1, but for 𝑣 lsr = 81 . 49 km s -1 . Bottom: Same as Figure 1, with Λ 0 = 3 . 01 × 10 -19 cm -3 s -1 and Λ 𝑖 = 4 . 62 × 10 -18 cm -3 s -1 , and pump amplitudes 𝑎 𝑚 of 1.00, 0.94, 0.98, 0.98, and 0.80, respectively. The model yields 𝑛 0 𝐿 = 4 . 00 × 10 3 cm -2 , corresponding to 𝑇 R , 0 = 6 . 63 × 10 4 s.  \n<!-- image -->  \nFigure 5. Full width at half maximum (FWHM) of the central peak of the autocorrelation curve for four flaring features in G22.356+0.066. The vertical axis shows the FWHMA in days, and the horizontal axis shows the velocity of the feature in km s -1 . The velocity corresponding to each FWHMA is explicitly labeled on the graph.  \n<!-- image -->  \nrelated timescales 𝜏 1 and 𝑇 𝑝 were held constant, the pump parameters Λ 0 and Λ 1 ,𝑚 were adjusted to fit the model to the observed data. These adjustments in pump amplitude successfully reproduced the varying flare durations, which show a notable increase from 𝑣 lsr = 79 . 30 km s -1 to 𝑣 lsr = 81 . 49 km s -1 , as quantified by autocorrelation analysis and illustrated in Figure 5. The variations in pump amplitude can be interpreted as differences in the coupling of the pump to the systems responsible for the flares, potentially due to factors such as orientation, sample geometry, and other local conditions (Rajabi et al. 2023).  \nIt is commonly understood that 6.7 GHz methanol masers are driven by infrared (IR) radiation through radiative pumping (Sobolev et al. 1997). A nearby variable IR source within the flar-  \ning region of G22.356+0.066 was cataloged by Lucas et al. (2008). Szymczak et al. (2015) mapped the 6.7 GHz methanol flares alongside this IR source, highlighting its location a few hundred AU east of the periodic flaring region. Although IR light curves are not available, the superradiance flares are not sensitive to the exact profile of the pumps used (e.g., symmetric or asymmetric), and similar results can be obtained with different pump profiles while maintaining a constant flare duration (Houde et al. 2024). The key point is that, within such a small physical proximity of the flares to the pump, the variable flare durations cannot be explained by traditional maser theory, where maser flares are expected to closely mimic changes in the pump for saturated masers. A common pump cannot account for the presence of substantial differences (a variation of more than 50%) in flare timescale. However, superradiance highlights a transient regime of radiation where the flare duration, even in response to the same pump duration, can vary depending on the physical conditions realized at the source.  \nThe maps of 6.7 GHz methanol flares in G22.356+0.066 presented by Bartkiewicz et al. (2009) and Szymczak et al. (2015) indicate that periodic flaring sources at different velocities are concentrated within a compact region of about 80 AU. Consequently, to fit the MBEs model to the light curves of various velocity channels, the dephasing and relaxation timescales, 𝑇 1 and 𝑇 2 , were kept constant, implying a similar environment across the source.  \nRemarkably, the timescales 𝑇 1 ≃ 205 days and 𝑇 2 ≃ 11 . 8 days, derived here closely align with those reported by Rajabi et al. (2023) for modelling the periodic flares in G9.62+0.20E, further suggesting comparable environmental conditions between the two sources. Finally, G22.356+0.066 is the second example of periodic flaring, alongside G9.62+0.20E, where we successfully model flares of different durations originating from a close physical vicinity, using the MBE framework and the superradiance mechanism.', 'APPENDIX A: MAXWELL-BLOCH EQUATIONS': "The models are derived from the MBEs (Gross & Haroche 1982; Benedict et al. 1996; Rajabi et al. 2023):  \n𝜕𝑛 ' 𝜕𝜏 = 𝑖 ℏ GLYPH<16> 𝑃 + 0 𝐸 + 0 -𝑃 -0 𝐸 -0 GLYPH<17> -𝑛 ' 𝑇 1 + Λ 𝑁 (A1)  \n𝜕𝑃 + 0 𝜕𝜏 = 2 𝑖𝑑 2 ℏ 𝐸 -0 𝑛 ' -𝑃 + 0 𝑇 2 + Λ 𝑃 (A2)  \n𝜕𝐸 + 0 𝜕𝑧 = 𝑖𝜔 0 2 𝜖 0 𝑐 𝑃 -0 , (A3)  \nwhere the slow-varying-envelope and rotating wave approximations are applied. Here, 𝜏 ≡ 𝑡 -𝑧 / 𝑐 represents the retarded time. The variable 𝑛 ' denotes half of the inverted population density ( 𝑛 = 2 𝑛 ' ) , while 𝑃 + 0 and 𝐸 + 0 are the amplitudes of the polarization and electric field, respectively.  \nTo account for the inversion pump, the term Λ 𝑁 , as defined in equation (2), is included in equation (A1). Internal fluctuations in 𝑛 ' and 𝑃 + 0 are introduced to initiate the temporal evolution of the system. To represent these internal fluctuations in polarization, a phenomenological polarization pump term Λ 𝑃 is added to equation (A2). Non-coherent relaxation and polarization dephasing are modelled by the decay terms -𝑛 ' / 𝑇 1 and -𝑃 + 0 / 𝑇 2 , respectively. The polarization and electric field vectors are:  \nFigure B1. Removal of the non-zero quiescent level from the SR fit for the 𝑣 = 80 . 35 km s -1 velocity channel in G22.356+0.066. The vertical axis shows the flux density in Jy, while the horizontal axis is for time in Julian Day -2450000. The solid blue curve illustrates the fit produced by the MBEs. The dashed orange line shows a third-order polynomial fit to the quiescent level. The solid vermilion curve displays the adjusted MBEs fit, which is obtained by subtracting the third-order polynomial fit from the original MBEs fit.  \n<!-- image -->  \n𝑷 ± ( 𝜏, 𝑧 ) = 𝑃 ± 0 ( 𝜏, 𝑧 ) 𝑒 ± 𝑖 𝜔 0 𝜏 𝝐 𝑑 (A4)  \n𝑬 ± ( 𝜏, 𝑧 ) = 𝐸 ± 0 ( 𝜏, 𝑧 ) 𝑒 ∓ 𝑖 𝜔 0 𝜏 𝝐 𝑑 , (A5)  \nwhere 𝝐 𝑑 = 𝒅 / 𝑑 is the unit vector for the dipole of the molecular transition, 𝜔 0 = 𝑐𝑘 is the transition angular frequency, and 𝑑 = | 𝒅 | is the transition electric dipole moment. The superscripts + and -denote transitions between the lower and upper levels and vice versa, respectively.", 'APPENDIX B: AUTOCORRELATION': "To quantitatively assess the differences in flare durations across features at different velocities, we perform an autocorrelation analysis of the flares. The width of the central peak in the autocorrelation function provides a measure of the mean flare duration. By comparing the full widths at half maximum (FWHM A ) of the autocorrelation functions for different velocities, we can demonstrate that the flaring features exhibit varying durations.  \nTo obtain smooth autocorrelation curves, we use fits to the data yielded by MBEs (solid blue curves in Figs. 1 to 4), as the data are limited and unevenly spaced. We first subtract a third-order polynomial function, which fits the quiescent maser level between flares, from our models, as shown in Fig. B1. We then compute the normalized autocorrelation using Python's scipy.signal.correlate function in full mode. The resulting autocorrelation is depicted in Fig. B2, where the vertical axis shows the normalized autocorrelation, and the horizontal axis indicates lags in days. The vertical distance between consecutive peaks is 178.5 days, corresponding to the assumed flare period in our model.  \nThe FWHM of the central peak for each autocorrelation curve is measured and reported in Fig. 5.", 'ACKNOWLEDGEMENTS': "F.R.'s research is supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) Discovery Grant RGPIN2024-06346. AB, PW, MS acknowledge support from the National Science Centre, Poland through grant 2021/43/B/ST9/02008. We  \nFigure B2. Autocorrelation curves for flaring features in G22.356+0.066 at velocities 79.30, 80.35, 81.32, and 81.49 km s -1 . The vertical axis shows the normalized autocorrelation, and the horizontal axis represents lags in days. The vertical dashed line marks a lag of 178.5 day, the period of the flaring events.  \n<!-- image -->  \nalso appreciate the valuable feedback provided by Prof. Martin Houde.", 'REFERENCES': 'Araya E. D., Hofner P., Goss W. M., Kurtz S., Richards A. M. S., Linz H., Olmi L., Sewiło M., 2010, ApJ, 717, L133  \nBartkiewicz A., Szymczak M., van Langevelde H. J., Richards A. M. S., Pihlström Y. M., 2009, A&A, 502, 155  \nBenedict M. G., et al., 1996, Super-radiance: Multiatomic Coherent Emission. IOP Publishing Ltd  \nBeuther H., et al., 2019, A&A, 628, A90  \nBreen S. L., Ellingsen S. P., 2011, MNRAS, 416, 178  \nBreen S. L., Ellingsen S. P., Caswell J. L., Green J. A., Voronkov M. A.,  \nAvison A., Fuller G. A., Quinn L. J., 2016, MNRAS, 459, 4066 Feld M. S., MacGillivray J. C., 1980, Topics in Current Physics, 21, 7 Goedhart S., Gaylard M. J., van der Walt D. J., 2004, MNRAS, 355, 553 Gross M., Haroche S., 1982, Phys. Rep., 93, 301  \nHoude M., Rajabi F., MacLeod G. C., Goedhart S., Tanabe Y., van den Heever S. P., Wyenberg C. M., Yonekura Y., 2024, in Hirota T., Imai H., Menten K., Pihlström Y., eds, IAU Symposium Vol. 380, Cosmic Masers: Proper Motion Toward the Next-Generation Large Projects. pp 399-413 ( arXiv:2307.05274 ), doi:10.1017/S1743921323002764  \nInayoshi K., Sugiyama K., Hosokawa T., Motogi K., Tanaka K. E. I., 2013, ApJ, 769, L20  \nLucas P. W., et al., 2008, MNRAS, 391, 136  \nMacGillivray J. C., Feld M. S., 1976, Phys. Rev. A, 14, 1169  \nMathews A., 2017, PhD thesis, Honours thesis, The University of Western Ontario  \nRajabi F., Houde M., 2016, ApJ, 826, 216  \nRajabi F., Houde M., 2017, Science Advances, 3, e1601858  \nRajabi F., Houde M., 2020, MNRAS, 494, 5194  \nRajabi F., Houde M., Bartkiewicz A., Olech M., Szymczak M., Wolak P., 2019, MNRAS, 484, 1590  \nRajabi F., Houde M., MacLeod G. C., Goedhart S., Tanabe Y., van den Heever S. P., Wyenberg C. M., Yonekura Y., 2023, MNRAS, 526, 443 Reid M. J., et al., 2019, ApJ, 885, 131  \nSobolev A. M., Cragg D. M., Godfrey P. D., 1997, MNRAS, 288, L39 Szymczak M., Wolak P., Bartkiewicz A., van Langevelde H. J., 2011, A&A,  \n531, L3  \nSzymczak M., Wolak P., Bartkiewicz A., Borkowski K. M., 2012, AN, 333, 634  \nSzymczak M., Wolak P., Bartkiewicz A., 2015, MNRAS, 448, 2284  \nZhang S., et al., 2023, MNRAS, 520, 322  \nvan der Walt D. J., 2011, AJ, 141, 152  \nvan der Walt D. J., Goedhart S., Gaylard M. J., 2009, MNRAS, 398, 961'}