Accurately localizing fast radio bursts (FRBs) is essential for understanding their birth environments and for their use as cosmological probes. Recent advances in radio interferometry, particularly with MeerKAT, have enabled the localization of individual bursts with arcsecond precision. In this work, we present the localization of 15 apparently non-repeating FRBs detected with MeerKAT. Two of the FRBs, discovered in 2022, were localized in 8 s images from the projects that MeerTRAP was commensal to, while eight were localized using the transient buffer (TB) pipeline, and another one through SeeKAT, all with arcsecond precision. Four additional FRBs lacked TB triggers and sufficient signal, limiting their localization only to arcminute precision. For eight of the FRBs in our sample, we identify host galaxies with greater than 90 per cent confidence, and one with 80 per cent confidence, while two FRBs have ambiguous associations. We measured spectroscopic redshifts for six host galaxies, ranging from 0.33 to 0.85, demonstrating MeerKAT's sensitivity to high-redshift FRBs. We modelled the spectral energy distributions of host galaxies with sufficient photometric coverage to derive their stellar population and star formation properties. This work represents one of the largest uniform samples of well-localized distant FRBs to date, laying the groundwork for using MeerKAT FRBs as cosmological probes and understand how FRB hosts evolve at high redshift.
We investigate the presence and spatial characteristics of the jet base emission in M87* at 230 GHz, enabled by the significantly enhanced (u,v) coverage in the 2021 Event Horizon Telescope (EHT) observations. The integration of the 12−m Kitt Peak Telescope (USA) and NOEMA (France) stations into the array introduces two critical intermediate-length baselines to SMT (USA) and IRAM 30−m (Spain), providing sensitivity to emission structures at spatial scales of ∼250 μas and ∼2500 μas (∼ 0.02 pc and ∼ 0.02 pc). Without these new baselines, previous EHT observations of the source in 2017 and 2018 lacked the capability to constrain emission on large scales, where a “missing flux” of order ∼1 Jy is expected to reside. To probe these scales, we analyzed closure phases–robust against station-based gain calibration errors–and model the jet base emission using a simple Gaussian component offset from the compact ring emission at spatial separations > 100 μas. Our analysis revealed a Gaussian feature centered at (ΔRA ≈ 320 μas, ΔDec. ≈ 60 μ as), projected separation of ≈ 5500 AU, with an estimated flux density of only ∼60 mJy, implying that most of the missing flux identified in previous EHT studies had to originate from different, larger scales. Brighter emission at the relevant spatial scales is firmly ruled out, and the data do not favor more complex models. This component aligns with the inferred position of the large-scale jet and is therefore physically consistent with the emission of the jet base. While our findings point to detectable jet base emission at 230 GHz, the limited coverage provided by only two intermediate baselines limits our ability to robustly reconstruct its morphology. Consequently, we treated the recovered Gaussian as an upper limit on the jet base flux density. Future EHT observations with expanded intermediate baseline coverage will be essential to constrain the structure and nature of this component with higher precision.
Sagittarius A* ( the supermassive black hole at the center of the Milky Way, provides a unique laboratory to study accretion dynamics and plasma processes near the event horizon. We investigated the variability and polarization properties of using ALMA observations during the 2018 Event Horizon Telescope campaign. We analyzed high-cadence full-polarization light curves from ALMA at millimeter wavelengths, performed time-series analysis, and investigated the temporal behavior during an X-ray flare observed by on 2018 April 24. The variability characteristics are compared with expectations from standard accretion flow models. Chandra We find low variability in total intensity (σ/μ < 10%), but significantly higher variability in linear and circular polarization (∼ 30% and ∼ 50%, respectively). A time-series analysis reveals red-noise variability, with power spectral densities between -2 and -3 across all Stokes parameters. Polarized intensity shows stable intra-day timescales, while total intensity exhibits more variable timescales, suggesting distinct emission regions, with polarization likely arising from a coherent structure. On April 24, a statistically significant inter-band delay in polarized intensity coincides with a near-simultaneous X-ray and millimeter peak that deviates from the typical delayed flare scenario. This event also features enhanced millimeter variability and coherent polarization loop evolution. The observed simultaneity challenges standard models of transient synchrotron emission with cooling delays, favoring instead a scenario of continuous energy injection in an optically thin region. Our results offer new constraints on the physical mechanisms driving variability in and provide key observational input for refining theoretical models of accretion and plasma behavior in the vicinity of supermassive black holes.
Context. One of the imminent science goals of pulsar timing arrays (PTAs) is the detection of a continuous gravitational wave (CGW) emitted by an individual supermassive black hole binary (SMBHB). SMBHBs that cause CGWs with GW frequencies fGW > 10 nHz (high-frequency end of the nano-Hertz GW spectrum) have undergone significant orbital evolution and hence a change in fGW over time. In PTA datasets with a sufficiently long observational time span, this means that the contributions of the Earth and the pulsar terms to the CGW signal signature can eventually become resolvable. Since the pulsar term is accumulated incoherently and thus often treated as an additional source of noise, this separation can prove to be beneficial for the detection of the CGW signal in the PTA dataset. Aims. We investigate to which extent resolvable Earth and pulsar terms affect currently used techniques for CGW searches with PTA datasets, which treat the pulsar term as an additional source noise. We focus on the dependence of the pulsar term frequencies on the pulsar distance. We aim to answer the question of whether adding more distant pulsars to a PTA dataset can mitigate biases and improve the detection of CGWs. Methods. We used simulated PTA datasets based on the EPTA DR2 and IPTA DR2 pulsars in order to study the performance of the Earth-term-only Bayesian parameter estimation of the circular SMBHB model parameters and the frequentist narrow-band optimal statistic in the light of resolved pulsar terms due to larger pulsar distances. Results. We show that under ideal conditions, more distant pulsars can facilitate the CGW search with PTA datasets. The Bayesian parameter estimation is yielding better parameter constraints, and the frequentist search becomes more stable. Based on realistic dataset simulations, however, we found that other configuration parameters of a PTA, such as the anisotropic distribution of pulsars and the effective number of pulsars in a PTA, can play a crucial role in the importance of this effect.
We perform searches for gravitational-wave memory in the data of two major Pulsar Timing Array (PTA) experiments located in Europe and Australia. Supermassive black hole binaries (SMBHBs) are the primary sources of gravitational waves in PTA experiments. We develop and carry out the first search for late inspirals and mergers of these sources based on full numerical relativity waveforms with null (nonlinear) gravitational-wave memory. Additionally, we search for generic bursts of null gravitational-wave memory, exploring the possibilities of reducing the computational costs of these searches through kernel density estimation and normalizing flow approximations of the posteriors. We rule out mergers of SMBHBs with a chirp mass of 10 ^10 M _⊙ up to 700 Mpc over 18 yr of observation at 95% credibility. We rule out the observation of generic displacement memory bursts with strain amplitudes >10 ^−14 in brief periods of observation time but across the sky, or over the whole observation time but for certain preferred sky positions, at 95% credibility.
The discovery of persistent radio sources (PRSs) associated with repeating fast radio bursts (FRBs) has shed light on the immediate environments and possible progenitors of these FRBs. The confirmed PRSs may support the theory that FRB progenitors are compact central engines, whilst the non-detections suggest diversity of FRB's local environment. We perform a subarcsecond-resolution MeerKAT search at 1.28 GHz on 25 well-localized FRB positions provided by ASKAP and MeerTRAP. We detect 14 radio sources and provide flux upper limits for 12 non-detections (both these numbers include a source that was detected during two epochs of observation, and not detected during one epoch, adding up to 26). One radio source shows variability as seen in flux variations over three epochs of observation. Archival optical data reveal excesses in the direction of 13 detected radio sources. Similarly for four sources in the X-ray band, with one possibly being a high-energy signature of a radio galaxy core. Since we cannot definitively classify our detected radio sources as PRSs, future high-resolution observations with e-MERLIN will be required to resolve the radio emission and pronounce on the presence of compact PRSs associated with the 14 detected sources presented here.
Event Horizon Telescope (EHT) images of the supermassive black hole M87* depict an asymmetric ring of emission. General relativistic magnetohydrodynamic (GRMHD) models of M87* and its accretion disk predict that the amplitude and location of the ring's peak brightness asymmetry should fluctuate due to turbulence in the source plasma. We compare the observed distribution of brightness asymmetry amplitudes to the simulated distribution in GRMHD models, across varying black hole spin a(*). We show that, for strongly magnetized (MAD) models, three epochs of EHT data marginally disfavor divided by a(*)divided by less than or similar to 0.2. This is consistent with the Blandford-Znajek model for M87's jet, which predicts that M87* should have nonzero spin. We show quantitatively how future observations could improve spin constraints and discuss how improved spin constraints could distinguish between differing jet-launching mechanisms and black hole growth scenarios.
Millisecond pulsars (MSPs) are powerful probes of globular clusters (GCs), tracing stellar evolution, cluster dynamics, and the local gravitational potential. We investigate the MSP population in GC Omega Centauri. We perform Fourier-domain acceleration and jerk searches on MeerKAT observations, and carry out pulsar timing using MeerKAT and Parkes Murriyang data spanning 2021-2025. We fold Fermi LAT and NICER photons using updated radio ephemerides to search for high-energy pulsations. We discover a new isolated MSP, PSR J1326-4728S (hereafter S), with a spin period of 4.538 ms and a dispersion measure of 96.24 cm^3pc. We update the orbital parameters of all known binary systems, with those of I, N, and Q differing significantly from previous estimates, and obtain new timing solutions for G, H, and K. Pulsars B, G, H, K, and L exhibit black widow-like properties, I, N and Q are found in wider binaries, with N and Q having >0.2 M_⊙ companions, and N showing a significant orbital eccentricity (e=0.093). Significant spin period derivatives are measured for eight pulsars and interpreted as arising from the cluster gravitational potential. No pulsed high-energy emission is detected from individual pulsars. The inferred line-of-sight accelerations are consistent with a King-model gravitational potential. While our measurements are insensitive to an intermediate-mass black hole with mass 10^3-10^4 M_⊙, they place an upper limit of <10^5 M_⊙ at 90
We present the first Event Horizon Telescope 1.3 mm observations of the supermassive binary black hole candidate OJ 287. The observations achieved an unprecedented angular resolution of 18 μas and reveal significant structural and polarization variability over just five days, marking the shortest timescale on which such changes have been directly imaged in this source. The inner jet exhibits a twisted ridgeline structure, with features displaying apparent superluminal motions up to about 22 c. The linear polarization maps reveal three main polarized features whose electric-vector position angles (EVPAs) change substantially over the time span of our observations, including a component with a radial polarization consistent with being produced by a recollimation shock. Most notably, we directly resolved two innermost jet components whose EVPAs rotate in opposite directions. The faster component, moving at 2.4 ± 0.9 μas/day (17.4 ± 6.5 c), exhibits counterclockwise EVPA swings of roughly 3.7° per day, while the slower component, with a proper motion of 1.4 ± 0.3 μas/day (10.2 ± 2.2 c), rotates clockwise at approximately 2.5° per day. Previous studies inferred helical magnetic fields in AGN jets from time-resolved or integrated polarization variability but lacked the angular resolution to directly image this effect. Our results provide spatially resolved evidence that a helical magnetic field threads the jet’s collimation and acceleration zone, ruling out models based on the superposition of unresolved components. Our analysis suggests that propagating shocks interact with a Kelvin–Helmholtz plasma instability, illuminating different phases of the helical magnetic field and producing the observed polarization spatial and temporal variability. Moreover, our model naturally accounts for the more rapid polarization rotation observed in the faster moving component. Our model predicts even more rapid swings in polarization, which could be tested with future observations featuring a more densely sampled time coverage.
NICER has enabled mass-radius inferences for pulsars using pulse profile modeling (PPM), providing constraints on the equation of state (EOS) of cold, dense matter. To date, PPM and EOS inference have been carried out as two separate steps, with the former using EOS-agnostic priors. This approach has several drawbacks. Ideally, one would perform a fully hierarchical Bayesian inference where the pulse profile and EOS model parameters are jointly fit, but implementing such a framework is complex and computationally demanding. Here, we present an intermediate solution introducing an EOS-informed prior on mass radius into the existing PPM pipeline using normalizing flows. By focusing on the parameter space consistent with certain EOSs, this approach both tightens constraints on neutron star parameters while reducing computational costs and requiring minimal additional implementation effort. We test this approach on two pulsars, PSR J0740 & thorn; 6620 and PSR J0437-4715, and with two EOS model families: a model based on the speed of sound inside the neutron star interior (CS) and a piecewise-polytropic (PP) model. Both EOS models implement constraints from chiral effective field theory calculations of dense matter. For both pulsar datasets, the inferred radius credible intervals are narrower than in the EOS-agnostic case, with CS favoring smaller radii and PP favoring larger radii. For PSR J0437-4715, the EOS-informed priors reveal a new, more extreme geometric mode that is statistically favored but physically questionable. Including the PPM posteriors in the subsequent EOS inference further tightens the mass-radius posteriors through the chiral effective field theory constraints. However, there is also a sensitivity to the high-density extensions, where the PP (CS) model produces a shift toward larger (smaller) radii and corresponding stiffening (softening) of the pressure-energy density relation.
Context. Pulsar timing arrays (PTAs) have recently reported compelling evidence for the presence of a gravitational-wave background signal. Mapping the gravitational-wave background is key to understanding how it is formed, since anisotropy is a tracer for, for example, a supermassive black hole binary origin. Aims. In this work we refine the frequentist regularised gravitational-wave mapping analysis developed in our previous work (as part of the MeerKAT PTA 4.5-year data release). We derive a point-spread function describing the angular resolution of a PTA. Methods. We investigate how the point spread function changes for different PTA constellations and determine the best possible angular resolution achievable within our framework. Using simulated data, we demonstrate that previous methods do not capture the actual resolution - especially in regions of the sky with a high density of pulsars. Results. We propose an improved scheme that accounts for a variable local resolution and test it using realistic simulations of the latest MeerKAT dataset. We demonstrate that we are able to identify a continuous gravitational wave signal in a region with good pulsar sky coverage with approximately a factor of two increase in significance compared to our previous method.
Pulsars are typically characterized by their stable, highly magnetized, and fast-rotating nature, which underpins their persistent radio emissions. However, the discovery of prolonged radio-quiet (“off”) states in intermittent pulsars, such as PSR B1931+24, has been challenging the fundamental theory of pulsar magnetospheric emission. Despite long-term monitoring with several telescopes over 20 yr, the observations of PSR B1931+24 had not detected any significant emission during its “off” states. Recently, sensitive observations via the Five-hundred-meter Aperture Spherical radio Telescope revealed the mysterious weak emissions containing occasional faint bursting pulses during the “off” states of PSR B1931+24, signifying the theory of pulsar magnetospheric emission that the pulsar cannot remain off for long unless it truly dies. Along with a substantial decrease in flux density, a contraction and change in the integrated pulse profile are measured in the “off” state compared to the radio-loud (“on”) state, indicating alterations in the plasma supply and magnetospheric structure. Additionally, previously unobserved asynchronous, nonuniform emission patterns are found in both states, implying a spatially inhomogeneous pair-cascade associated with the coherent radio emission of the pulsar. Furthermore, statistical comparisons between the faint off-state bursts and the on-state emission show that the faint bursts cannot be explained simply as a scaled-down version of the on-state emission. Instead, they constitute a distinct emission mode, plausibly linked to changes in plasma supply and the magnetospheric configuration of PSR B1931+24. These findings greatly advance our understanding of pulsar magnetospheric dynamics and emission mechanisms.
GPM 1839-10 is an intriguing long-period radio transient (LPT), distinguished by its activity spanning at least three decades and its highly unusual emission characteristics. These features include orthogonal polarization mode (OPM) switches, down-drifting sub-structures, and distinct linear-to-circular polarization conversion behaviors. In this work, we present follow-up observations utilizing the FAST telescope at L-band, yielding a total of seven detected radio pulses. We find a consistent association between OPM switches and a decrease in polarized intensity. This feature strongly supports the hypothesis that the OPM switches are generated by the incoherent summation of OPMs. Our measured Rotation Measures (RMs) are consistent with previous observations, indicating that the magneto-ionic environment is stable. If the source is in a binary system, such stability suggests it may host a weakly magnetized companion. Crucially, we firstly observe clear evidence of a cyclotron absorption feature in one radio pulse, a signature rarely observed in radio sources. This feature allows us to infer that the magnetic field strength at the absorption site has a lower limit of tens of Gauss, which is necessary for the phenomenon to occur. This characteristic can be explained in a scenario where GPM 1839-10 possesses a weakly magnetized companion star.
This paper presents a comprehensive wideband study of FRB 20240619D focusing on its hyperactivity, rotation measure evolution, and the search for an associated persistent radio source. Using data from the MeerKAT, Murriyang, and Lovell telescopes, we analysed the spectral, temporal, and polarimetric properties of 1539 bursts. Our observations reveal a remarkably high burst rate of 161 bursts per hour in early August above a fluence value of 1.6 Jy ms as well as significant secular variations in rotation measure and diverse polarization characteristics, including high linear polarization fractions and occasional circular polarization. The burst activity also showed frequency dependence with approximately 61 per cent of the total number of bursts detected between 1300 and 1800 MHz. The burst activity of FRB 20240619D ceased abruptly after a period of intense activity lasting approximately 80 d, suggesting an episodic behaviour. Follow-up observations with MeerKAT and Australia Telescope Compact Array did not reveal an associated compact persistent radio source. Altogether, our results highlight the importance of continued long-term monitoring and multiwavelength observations in understanding the emission mechanisms and diversity of progenitor populations of fast radio bursts.
A striking aspect of the radio profiles of many millisecond pulsars (MSPs) is that they consist of components separated from each other by regions lacking in emission. We devise a technique for determining 'disjoint' from 'contiguous' components and show that 39 per cent of MSPs have disjoint components as opposed to only 3 per cent of the slow pulsar population. We surmize that the pulsars with these disjoint components show evidence for both emission above the polar cap and from the current sheet beyond the light cylinder (LC), colocated with gamma-ray emission. For some of the radio MSPs only the LC emission is being observed. It is our contention that almost all of the current population of gamma-ray MSPs show evidence for colocated radio emission. A simple geometric explanation allows the presence (or not) of LC emission and the relationship (or not) between the gamma-ray and radio profiles to be determined. The LC components have frequently very high polarization and typically flat position-angle traverses thus helping to explain the difficulties in determining the geometry of MSPs. In cases where the geometry can be determined the values broadly align with expectations. In this picture, the number of potentially detectable radio MSPs is higher than previously thought, although the actual detectability of LC components depends on their luminosity function. A mechanism is required to produce coherent radio emission far from the stellar surface. These ideas have implications for our understanding of the populations of radio-loud and radio-quiet rotation-powered millisecond pulsars, and may have implications for the long-term timing stability of some of these sources.
In very-long baseline interferometric arrays, nearly co-located stations probe the largest scales and typically cannot resolve the observed source. In the absence of a large-scale structure, closure phases constructed with these stations are zero and, since they are independent of station-based errors, they can be used to probe data issues. Here, we show how these trivial closure phases become nonzero with a brightness distribution on smaller scales than their short baseline would suggest. When applied to sources that are made up of a bright compact and large-scale diffuse component, the trivial closure phases directly measure the centroid relative to the compact source and higher-order image moments. We present a technique to measure these image moments with minimal model assumptions and validate it on synthetic Event Horizon Telescope (EHT) data. We then apply this technique to 2017 and 2018 EHT observations of M87* and find a weak preference for extended emission in the direction of the large-scale jet. We also apply it to 2021 EHT data and measure the source centroid about 1 mas northwest of the compact ring, which is consistent with the jet observed at lower frequencies.
The Fermi Large Area Telescope (LAT) provides advantages for radio pulsar searches by enabling efficient target selection. We can confidently point radio telescopes to the positions of Fermi unidentified gamma-ray sources that have a high probability of hosting a pulsar. As part of Transients and Pulsars with MeerKAT (TRAPUM), we conducted a survey of Fermi-LAT sources using the Ultra High Frequency (UHF; 544-1088 MHz) receiver of the MeerKAT radio telescope. We observed 79 sources that were identified as pulsar-like candidates using a random forest technique from the Fermi-LAT Fourth Source Catalogue. We observed each target for 10 min at two separate epochs. As a result, we discovered nine new millisecond pulsars (MSPs) and six slow pulsars. Based on the radio discoveries, we also searched for gamma-ray pulsations, confirming that seven of the newly discovered MSPs are associated with Fermi-LAT sources, and performed joint radio and gamma-ray pulsar timing. Companion mass estimates and evidence of radio eclipses indicate that among the nine MSPs there are three black widows and three redbacks. Lastly, we compared the discovered pulsars in the MeerKAT UHF survey against the previous Fermi sources TRAPUM survey at L band, concluding the superiority of UHF observations in sensitivity to fainter pulsars and in detection rate than L band for finding new gamma-ray MSPs.
Context. The 2017 observing campaign of the Event Horizon Telescope (EHT) delivered the first very long baseline interferometry (VLBI) images at the observing frequency of 230 GHz, leading to a number of unique studies on black holes and relativistic jets from active galactic nuclei (AGN). In total, eighteen sources were observed, including the main science targets, Sgr A* and M 87, and various calibrators. Sixteen sources were AGN. Aims. We investigated the morphology of the sixteen AGN in the EHT 2017 data set, focusing on the properties of the VLBI cores: size, flux density, and brightness temperature. We studied their dependence on the observing frequency in order to compare it with the Blandford-Konigl (BK) jet model. In particular, we aimed to study the signatures of jet acceleration and magnetic energy conversion. Methods. We modeled the source structure of seven AGN in the EHT 2017 data set using linearly polarized circular Gaussian components (1749+096, 1055+018, BL Lac, J0132-1654, J0006-0623, CTA 102, and 3C 454.3) and collected results for the other nine AGN from dedicated EHT publications, complemented by lower frequency data in the 2-86 GHz range. Combining these data into a multifrequency EHT+ data set, we studied the dependences of the VLBI core component flux density, size, and brightness temperature on the frequency measured in the AGN host frame (and hence on the distance from the central black hole), characterizing them with power law fits. We compared the observations with the BK jet model and estimated the magnetic field strength dependence on the distance from the central black hole. Results. Our observations spanning event horizon to parsec scales indicate a deviation from the standard BK model, particularly in the decrease of the brightness temperature with the observing frequency. Only some of the discrepancies may be alleviated by tweaking the model parameters or the jet collimation profile. Either bulk acceleration of the jet material, energy transfer from the magnetic field to the particles, or both are required to explain the observations. For our sample, we estimate a general radial dependence of the Doppler factor delta proportional to r(<= 0.5). This interpretation is consistent with a magnetically accelerated sub-parsec jet. We also estimate a steep decrease of the magnetic field strength with radius B proportional to r(-3), hinting at jet acceleration or efficient magnetic energy dissipation.
We report the discovery of a new black-widow millisecond pulsar, PSR J1544-2555, associated with the Fermi-Large Area Telescope (LAT) source 4FGL J1544.2-2554. Optical, radio, and gamma-ray observations confirmed its nature as a compact spider binary system. Optical photometry from ULTRACAM revealed a similar to 2.7-h orbital period, guiding MeerKAT observations that detected similar to 2.4-ms radio pulsations. Subsequent timing campaigns using the Murriyang Parkes Telescope, the Effelsberg 100-m Radio Telescope, and the Nan & ccedil;ay Radio Telescope allowed us to obtain a preliminary timing solution, which enabled us to find gamma-ray pulsations. The final timing solution, spanning 16 yr of Fermi-LAT gamma-ray data, also displays orbital period variations typical of spider pulsars. X-ray observations from eROSITA indicate non-thermal emission, but the relatively low count rate prohibits the search for X-ray pulsations. Optical light curve modelling using Icarus suggests the asymmetry is best explained by a spot model, where uneven heating creates localized temperature variations on the companion. While the optical spectra we obtained are compatible with the physical properties we infer for the companion star, they were not of sufficient signal-to-noise to allow for radial velocity measurements, thus limiting constraints on the neutron star's mass. The observed bluer colour near the light curve minimum suggests possible non-thermal emission from intrabinary shocks, supported by the presence of an X-ray source. This discovery exemplifies the proven capability of the Fermi-LAT catalogue in identifying millisecond pulsar candidates and highlights the role of optical surveys in detecting variable sources suitable for radio follow-up.
We report three epochs of polarized images of M87* at 230 GHz using data from the Event Horizon Telescope (EHT) taken in 2017, 2018, and 2021. The baseline coverage of the 2021 observations is significantly improved through the addition of two new EHT stations: the 12 m Kitt Peak Telescope and the Northern Extended Millimetre Array (NOEMA). All observations result in images dominated by a bright, asymmetric ring with a persistent diameter of 43.9 +/- 0.6 mu as, consistent with expectations for lensed synchrotron emission encircling the apparent shadow of a supermassive black hole. We find that the total intensity and linear polarization of M87* vary significantly across the three epochs. Specifically, the azimuthal brightness distribution of the total intensity images varies from year to year, as expected for a stochastic accretion flow. However, despite a gamma-ray flare erupting in M87 quasi-contemporaneously to the 2018 observations, the 2018 and 2021 images look remarkably similar. The resolved linear polarization fractions in 2018 and 2021 peak at similar to 5%, compared to similar to 15% in 2017. The spiral polarization pattern on the ring also varies from year to year, including a change in the electric vector position angle helicity in 2021 that could reflect changes in the magnetized accretion flow or an external Faraday screen. The improved 2021 coverage also provides the first EHT constraints on jet emission outside the ring, on scales of less than or similar to 1 mas. Overall, these observations provide strong proof of the reliability of the EHT images and probe the dynamic properties of the horizon-scale accretion flow surrounding M87*.